[R-sig-ME] standard errors for predictions from lmer and/or glmm.admb

[Ben Bolker]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

On 06/02/2011 02:45 PM, Rafael Mares wrote:
> Dear all
> 
> I'm trying to obtain the standard errors for predictions from two
> different mixed effects models:
> 
> - a glmm using lmer (lme4 ver. 0.999375-39) with binomial errors
> (binary response variable), both continuous and categorical
> explanatory variables and two crossed random effects
> - a glmm using glmm.admb (glmmADMB ver. 0.5-2) with negative binomial
> errors (count response variable), both continuous and categorical
> explanatory variables and one random effect
> 
> I'm not comparing the two models, but I'd like to find a way to obtain
> the standard errors for the predictions (using new data sets) that
> would be applicable to both the lmer and the glmm.admb models... not
> necessarily applicable to both, but I think it would help me
> understand what I'm doing. I've read some of the older posts on this
> mailing list regarding prediction errors and intervals, but I'm either
> not understanding them correctly or the methods aren't applicable to
> neither lme4 nor glmmADMB.
> 
> So far, I am able to obtain predicted values using new data sets, for
> both models, but I don't know how to get the standard errors. For the
> predictions, I'm using:
> 
> for the lmer model:
> 	mm = model.matrix(terms(my.binomial.model),new.data1)
> 	predictions = mm %*% fixef(my.binomial.model)	# not transformed to
> the original scale
> 	
> for the glmm.admb model:
> 	mm = model.matrix(my.negbinom.model$fixed,new.data2)
> 	predictions = mm%*%my.negbinom.model$b	# not transformed to the original scale
> 	
> Thank you in advance for any help. I've copied the model calls and
> parts of the outputs below, in case that helps.
> 
> All the best,
> 
> Raff
> 

  See the code on http://glmm.wikidot.com/faq.

  The basic idea is that if fm1 is a fitted model and mm is a model
matrix based on the new data, then

 mm %*% fixef(fm1)  [for nlme/lme4-alikes]  or
 mm %*% coef(fm1)   [for some other models]

 gives the point predictions.

pvar1 <- diag(mm %*% tcrossprod(vcov(fm1),mm))

  Gives the prediction variances based only on the uncertainty in the
fixed effect predictors.

 pvar1+VarCorr(fm1)$Subject[1]

  adds the among-subject variance to this variance (will work only for
lme4 -- the bleeding-edge version of glmm.admb has a VarCorr extractor
function, but otherwise you'll have to dig out the random-effects
variance yourself).

  However, these approaches are quite crude -- they don't deal with
uncertainty in the variance parameters.

  Ben Bolker
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.10 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org/

iEYEARECAAYFAk3n3a8ACgkQc5UpGjwzenPseACfT/mNEdtqVWCG/1sbe5cYgQIx
8JMAn3JugkkWzYISPFQfmloI4Ev2jQKp
=ElSL
-----END PGP SIGNATURE-----