[R] mcmcsamp - How does it work?

[Douglas Bates]

On 10/20/06, Cleber N. Borges <cborges at iqm.unicamp.br> wrote:

> Hello,

> I am a chemical student and I make use of 'lme/lmer function'
> to handle experiments in split-plot structures.
> I know about the mcmcsamp and I think that it's very promising

> I would like knowing "the concept behind" of the mcmcsamp function.
> I do not want the C code of the MCMCSAMP function.
> I would like to get the "pseudo-algorithm" to understanding
> that it does.

The 'mcmc' in the name "mcmcsamp" stands for Markov chain Monte Carlo,
a technique to derive a sample from the posterior distribution of the
parameters in a (Bayesian) statistical model.  You may want to look at
the Wikipedia entries on "Markov chain Monte Carlo" and also on "Gibbs
sampling" for more information on this.

In Gibbs sampling we partition the parameter vector into disjoint
subsets and cycle through sampling from one subset conditional on the
current values of all the others.  Often the subsets are of size 1
which makes the sampling easier but causes problems with serial
correlation in the resulting chain.  For a linear mixed effects model
we can divide the parameters into three subsets

 1) the variance, \sigma^2, of the per-observation noise term
 2) the parameters that determine the variance-covariance matrix of
the random effects
 3) the random effects and the fixed effects (in the Bayesian
formulation on the model the random effects are regarded as
parameters).

Conditional on the other two subsets and on the data, we can sample
directly from the posterior distribution of the remaining subset.  For
the first subset we sample from a chi-squared distribution conditional
on the current residuals.  The prior for the variances and covariances
of the random effects is chosen so that for the second subset we
sample from a Wishart distribution.  Finally, conditional on the first
two subsets and on the data the sampling for the third subset is from
a multivariate normal distribution.

Starting from the ML or REML estimates of the parameters in the model
we cycle through these steps many times to generate a sample from the
posterior distribution of the parameters.  Code in some of the threads
mentioned below is then used to make inferences based on this sample.

> I accompanied the threads [1] and [2] and 'googled' several terms but
> I do not obtained success about this matter.
>
> [ 1 ] fixed effects following lmer and mcmcsamp - which to present? -
> 8/8/2006
> [ 2 ] mcmcsamp-CI instead of P-values - references? - 3/10/2006