Matrix (version 0.995-15)

mcmcsamp: Generate an MCMC sample

Description

This generic function generates a sample from the posterior distribution of the parameters of a fitted model using Markov Chain Monte Carlo methods.

Usage

mcmcsamp(object, n, verbose, ...)

Arguments

object

An object of a suitable class - usually an lmer or glmer object.

n

integer - number of samples to generate. Defaults to 1.

verbose

logical - if TRUE verbose output is printed. Defaults to FALSE.

...

Some methods for this generic function may take additional, optional arguments. The method for lmer objects takes the optional argument saveb which, if TRUE, causes the values

Value

• An object of (S3) class "mcmc" suitable for use with the functions in the "coda" package.

Examples

require("lattice", quietly = TRUE, character = TRUE)
data(sleepstudy)
(fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
set.seed(101)
samp1 <- mcmcsamp(fm1, n = 1000)
frm <-
    data.frame(vals = c(samp1), iter = rep(1:nrow(samp1), ncol(samp1)),
    par = factor(rep(1:ncol(samp1), each = nrow(samp1)),labels = colnames(samp1)))
densityplot(~ vals | par, frm, plot = FALSE,
            scales = list(relation = 'free', x = list(axs='i')))
xyplot(vals ~ iter | par, frm, layout = c(1, ncol(samp1)),
       scales = list(x = list(axs = "i"), y = list(relation = "free")),
       main = "Trace plot", xlab = "Iteration number", ylab = "",
       type = "l")
qqmath(~ vals | par, frm, type = 'l',
       scales = list(y = list(relation = 'free')))
if (require("coda", quietly = TRUE, character = TRUE)) {
   print(summary(samp1))
   print(autocorr.diag(samp1))
}
(eDF <- mean(samp1[,"deviance"]) - deviance(fm1, REML=FALSE)) # potentially useful approximate D.F.
stopifnot(abs(eDF - 7) < 1)