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

Description

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.

Usage

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.

Value

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.

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.

Examples

Run this code

# NOT RUN {
## 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)

# }