rblme function

Fits robust Bayesian Linear Mixed-effects Models (LMM) to data via robust REML estimating functions.

This function fits robust Bayesian LMMs to data via robust REML estimating functions. The latters are those proposed by Richardson & Welsh (1995), which are robustified versions of restricted maximum likelihood (REML) estimating equations. Posterior sampling is done with an ABC-MCMC algorithm, where the data are summarised through a rescaled version of the aforementioned estimating functions; see Ruli et al. (2017) for the properties and details of the method. The current package version (0.1.2) supports only models with a single random effects. Extensions to more general settings will be provided in the future versions of the package.

rblme(nabc, h.obj, chain.control = list(trace.init = NULL, thin.by = NULL),
        n.cores = 1)

Arguments

• nabc: number of posterior samples.
• h.obj: list of objects as returned by the tune.h function. Hence tune.h must be called first.
• chain.control: parameters that control the tracing and the thinning of the chain(s).
• n.cores: number of cores for parallel computation on non Windows machines. For n.cores>2, n.cores chains are run each on a different core with using the same parameters but with a different random seed.

Returns

list or list of lists with elements abc and effi. In case of n.cores=1, effi is the actual acceptance rate of the ABC-MCMC algorithm whereas in abc are stored the posterior samples. The latters are stored as a $(q + c) \times$nabc matrix, where $q$ is the number of fixed effects, i.e. the number of columns in the design matrix and $c = 2$ is the number of variance components. Hence, the first $q$ rows of the matrix abc give the posterior samples for the fixed effects and the last two rows give the posterior samples for the log-variances of the fixed effects and the residual term, respectively. If n.cores > 1, i.e. if simulations are performed in parallel, then a list of lists is returned, where each element of the list is a list with elements abc and effi, where abc and effi are as those aforementioned.

Examples

## The following example is meant for function documentation.
## For realistic use probably you'll need to take a larger sample and choose a
## "better" bandwidth h.

data(ergoStool)

require(lme4)
fm1 <- lmer(effort~Type + (1| Subject), data = ergoStool)

## tune h to get 0.8% acceptance
hopt <- tune.h(effort~Type + (1|Subject), data = ergoStool, n.samp = 1e+4,
               acc.rate = 0.01, n.sim.HJ = 100, grid.h = seq(0.3, 0.7, len = 3),
               prior = list(beta.sd = 10, s2.scale = 5), n.cores = 1)

## draw posterior samples with hopt.
abc.tmp <- rblme(nabc = 1e+4, h.obj = hopt,
                 n.cores = 1)

# process ABC samples
abc.sim <- t(abc.tmp$abc)
abc.sim[,c(5,6)] <- exp(abc.sim[,c(5,6)])

# ABC posterior
colMeans(abc.sim)

# REML estimates
summary(fm1)

References

Ruli E., Sartori N. & Ventura L. (2017) Robust approximate Bayesian inference with an application to linear mixed models. https://arxiv.org/abs/1706.01752

Richardson A. M. & Welsh A. H. (1995) Robust restricted maximum likelihood in mixed linear models. Biometrics 51 , 1429-1439.

See Also

tune.h, ergoStool.