This function simulates from the conditional distribution of the random effects of binomial and Poisson models.
Laplace.sampling.lr(
mu,
sigma2,
K,
y,
units.m,
control.mcmc,
messages = TRUE,
plot.correlogram = TRUE,
poisson.llik = FALSE
)
mean vector of the linear predictor.
variance of the random effect.
random effect design matrix, or kernel matrix for the low-rank approximation.
vector of binomial/Poisson observations.
vector of binomial denominators, or offset if the Poisson model is used.
output from control.mcmc.MCML
.
logical; if messages=TRUE
then status messages are printed on the screen (or output device) while the function is running. Default is messages=TRUE
.
logical; if plot.correlogram=TRUE
the autocorrelation plot of the conditional simulations is displayed.
logical; if poisson.llik=TRUE
a Poisson model is used or, if poisson.llik=FALSE
, a binomial model is used.
A list with the following components
samples
: a matrix, each row of which corresponds to a sample from the predictive distribution.
h
: vector of the values of the tuning parameter at each iteration of the Langevin-Hastings MCMC algorithm.
Binomial model. Conditionally on \(Z\), the data y
follow a binomial distribution with probability \(p\) and binomial denominators units.m
. Let \(K\) denote the random effects design matrix; a logistic link function is used, thus the linear predictor assumes the form $$\log(p/(1-p))=\mu + KZ$$ where \(\mu\) is the mean vector component defined through mu
.
Poisson model. Conditionally on \(Z\), the data y
follow a Poisson distribution with mean \(m\lambda\), where \(m\) is an offset set through the argument units.m
. Let \(K\) denote the random effects design matrix; a log link function is used, thus the linear predictor assumes the form $$\log(\lambda)=\mu + KZ$$ where \(\mu\) is the mean vector component defined through mu
.
The random effect \(Z\) has iid components distributed as zero-mean Gaussian variables with variance sigma2
.
Laplace sampling. This function generates samples from the distribution of \(Z\) given the data y
. Specifically, a Langevin-Hastings algorithm is used to update \(\tilde{Z} = \tilde{\Sigma}^{-1/2}(Z-\tilde{z})\) where \(\tilde{\Sigma}\) and \(\tilde{z}\) are the inverse of the negative Hessian and the mode of the distribution of \(Z\) given y
, respectively. At each iteration a new value \(\tilde{z}_{prop}\) for \(\tilde{Z}\) is proposed from a multivariate Gaussian distribution with mean $$\tilde{z}_{curr}+(h/2)\nabla \log f(\tilde{Z} | y),$$
where \(\tilde{z}_{curr}\) is the current value for \(\tilde{Z}\), \(h\) is a tuning parameter and \(\nabla \log f(\tilde{Z} | y)\) is the the gradient of the log-density of the distribution of \(\tilde{Z}\) given y
. The tuning parameter \(h\) is updated according to the following adaptive scheme: the value of \(h\) at the \(i\)-th iteration, say \(h_{i}\), is given by $$h_{i} = h_{i-1}+c_{1}i^{-c_{2}}(\alpha_{i}-0.547),$$
where \(c_{1} > 0\) and \(0 < c_{2} < 1\) are pre-defined constants, and \(\alpha_{i}\) is the acceptance rate at the \(i\)-th iteration (\(0.547\) is the optimal acceptance rate for a multivariate standard Gaussian distribution).
The starting value for \(h\), and the values for \(c_{1}\) and \(c_{2}\) can be set through the function control.mcmc.MCML
.