MCMC Estimation for Mixed Effects Model

Fits a mixed effects model via MCMC. The outcome can be normally distributed or ordinal (Goldstein, 2011; Goldstein, Carpenter, Kenward & Levin, 2009).

ml_mcmc( formula, data, iter=3000, burnin=500, print_iter=100, outcome="normal",
     nu0=NULL, s0=1, psi_nu0_list=NULL, psi_S0_list=NULL, inits_lme4=FALSE,
     thresh_fac=5.8, ridge=1e-5)

# S3 method for ml_mcmc
summary(object, digits=4, file=NULL, ...)

# S3 method for ml_mcmc
plot(x, ask=TRUE, ...)

# S3 method for ml_mcmc
coef(object, ...)

# S3 method for ml_mcmc
vcov(object, ...)

ml_mcmc_fit(y, X, Z_list, beta, Psi_list, sigma2, alpha, u_list, idcluster_list,
    onlyintercept_list, ncluster_list, sigma2_nu0, sigma2_sigma2_0, psi_nu0_list,
    psi_S0_list, est_sigma2, est_probit, parameter_index, est_parameter, npar, iter,
    save_iter, verbose=TRUE, print_iter=500, parnames0=NULL, K=9999, est_thresh=FALSE,
    thresh_fac=5.8, ridge=1e-5, parm_summary=TRUE )

## exported Rcpp functions
miceadds_rcpp_ml_mcmc_sample_beta(xtx_inv, X, Z_list, y, u_list, idcluster_list, sigma2,
     onlyintercept_list, NR, ridge)
miceadds_rcpp_ml_mcmc_sample_u(X, beta, Z_list, y, ztz_list, idcluster_list,
    ncluster_list, sigma2, Psi_list, onlyintercept_list, NR, u0_list, ridge)
miceadds_rcpp_ml_mcmc_sample_psi(u_list, nu0_list, S0_list, NR, ridge)
miceadds_rcpp_ml_mcmc_sample_sigma2(y, X, beta, Z_list, u_list, idcluster_list,
     onlyintercept_list, nu0, sigma2_0, NR, ridge)
miceadds_rcpp_ml_mcmc_sample_latent_probit(X, beta, Z_list, u_list, idcluster_list, NR,
     y_int, alpha, minval, maxval)
miceadds_rcpp_ml_mcmc_sample_thresholds(X, beta, Z_list, u_list, idcluster_list, NR, K,
     alpha, sd_proposal, y_int)
miceadds_rcpp_ml_mcmc_predict_fixed_random(X, beta, Z_list, u_list, idcluster_list, NR)
miceadds_rcpp_ml_mcmc_predict_random_list(Z_list, u_list, idcluster_list, NR, N)
miceadds_rcpp_ml_mcmc_predict_random(Z, u, idcluster)
miceadds_rcpp_ml_mcmc_predict_fixed(X, beta)
miceadds_rcpp_ml_mcmc_subtract_fixed(y, X, beta)
miceadds_rcpp_ml_mcmc_subtract_random(y, Z, u, idcluster, onlyintercept)
miceadds_rcpp_ml_mcmc_compute_ztz(Z, idcluster, ncluster)
miceadds_rcpp_ml_mcmc_compute_xtx(X)
miceadds_rcpp_ml_mcmc_probit_category_prob(y_int, alpha, mu1, use_log)
miceadds_rcpp_pnorm(x, mu, sigma)
miceadds_rcpp_qnorm(x, mu, sigma)
miceadds_rcpp_rtnorm(mu, sigma, lower, upper)

Arguments

formula	An R formula in lme4-like specification
data	Data frame
iter	Number of iterations
burnin	Number of burnin iterations
print_iter	Integer indicating that every `print_iter`th iteration progress should be displayed
outcome	Outcome distribution: `"normal"` or `"probit"`
nu0	Prior sample size
s0	Prior guess for variance
inits_lme4	Logical indicating whether initial values should be obtained from fitting the model in the lme4 package
thresh_fac	Factor for proposal variance for estimating thresholds which is determined as `thresh_fac`\(/N\) (\(5.8/N\) as default).
ridge	Ridge parameter for covariance matrices in sampling
object	Object of class `ml_mcmc`
digits	Number of digits after decimal used for printing
file	Optional file name
...	Further arguments to be passed
x	Object of class `ml_mcmc`
ask	Logical indicating whether display of the next plot should be requested via clicking
y	Outcome vector
X	Design matrix fixed effects
Z_list	Design matrices random effects
beta	Initial vector fixed coefficients
Psi_list	Initial covariance matrices random effects
sigma2	Initial residual covariance matrix
alpha	Vector of thresholds
u_list	List with initial values for random effects
idcluster_list	List with cluster identifiers for random effects
onlyintercept_list	List of logicals indicating whether only random intercepts are used for a corresponding random effect
ncluster_list	List containing number of clusters for each random effect
sigma2_nu0	Prior sample size residual variance
sigma2_sigma2_0	Prior guess residual variance
psi_nu0_list	List of prior sample sizes for covariance matrices of random effects
psi_S0_list	List of prior guesses for covariance matrices of random effects
est_sigma2	Logical indicating whether residual variance should be estimated
est_probit	Logical indicating whether probit model for ordinal outcomes should be estimated
parameter_index	List containing integers for saving parameters
est_parameter	List of logicals indicating which parameter type should be estimated
npar	Number of parameters
save_iter	Vector indicating which iterations should be used
verbose	Logical indicating whether progress should be displayed
parnames0	Optional vector of parameter names
K	Number of categories
est_thresh	Logical indicating whether thresholds should be estimated
parm_summary	Logical indicating whether a parameter summary table should be computed
xtx_inv	Matrix
NR	Integer
ztz_list	List containing design matrices for random effects
u0_list	List containing random effects
nu0_list	List with prior sample sizes
S0_list	List with prior guesses
sigma2_0	Numeric
y_int	Integer vector
minval	Numeric
maxval	Numeric
sd_proposal	Numeric vector
N	Integer
Z	Matrix
u	Matrix containing random effects
idcluster	Integer vector
onlyintercept	Logical
ncluster	Integer
mu1	Vector
use_log	Logical
mu	Vector
sigma	Numeric
lower	Vector
upper	Vector

Details

Fits a linear mixed effects model \(y=\bm{X}\bm{beta}+ \bm{Z}\bm{u}+e\) with MCMC sampling. In case of ordinal data, the ordinal variable \(y\) is replaced by an underlying latent normally distributed variable \(y^\ast\) and the residual variance is fixed to 1.

Value

List with following entries (selection)

sampled_values

Sampled values

par_summary

Parameter summary

References

Goldstein, H. (2011). Multilevel statistical models. New York: Wiley. doi: 10.1002/9780470973394

Goldstein, H., Carpenter, J., Kenward, M., & Levin, K. (2009). Multilevel models with multivariate mixed response types. Statistical Modelling, 9(3), 173-197. doi: 10.1177/1471082X0800900301

Examples