mcmcsamp               package:Matrix               R Documentation

_G_e_n_e_r_a_t_e _a_n _M_C_M_C _s_a_m_p_l_e

_D_e_s_c_r_i_p_t_i_o_n:

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

_U_s_a_g_e:

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

_A_r_g_u_m_e_n_t_s:

  object: An object of a suitable class - usually an 'lmer' 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
          of the random effects in each sample to be saved.  Note that
          this can result in very large objects being saved.  Use with
          caution. A second optional argument is 'trans' which, if
          'TRUE' (the default), returns a sample of transformed
          parameters.  All variances are expressed on the logarithm
          scale and any covariances are converted to Fisher's "z"
          transformation of the corresponding correlation.

_V_a_l_u_e:

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

_M_e_t_h_o_d_s:

     _o_b_j_e_c_t = "_l_m_e_r" generate MCMC samples from the posterior
          distribution of the parameters of a linear mixed model or a
          generalized linear mixed model.  The prior on the fixed
          effects parameters is taken to be locally uniform.  The prior
          on the variance-covariance matrices of the random effects is
          taken to be the locally non-informative prior described in
          Box and Tiao (1973). Conditional on the current values of the
          random effects these are sampled from a Wishart distribution.

_E_x_a_m_p_l_e_s:

     require("lattice", quietly = TRUE, character = TRUE)
     (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)) # potentially useful approximate D.F.

