[R-sig-ME] Prediction, Poisson and Intercepts

Manuel Morales Manuel.A.Morales at williams.edu
Fri Oct 23 18:37:21 CEST 2009


Dear list,

I have a question about how to efficiently generate predictions from a
poisson mixed model in lmer. A previous post by Bert Gunter
(http://markmail.org/message/22ncrgbt76bj5bmo) suggested the following:

model.matrix(terms(a.lmer.model),new.data) %*% fixef(a.lmer.model)

However, the back-transformed values from these predictions do not match
the original means because the fixed effect estimate of the intercept
does not include the random effects (see code below).

One alternative is:

new.Intercept <- log(mean(exp(coef(a.lmer.model)[[1]][,1])))
new.fixed <- fixef(a.lmer.model)
new.fixed[1] <- new.Intercept
model.matrix(terms(a.lmer.model),newdata) %*% new.fixed

But this seems awfully clunky ... Any suggestions for a more elegant
solution?

Demonstration code for mean model below:
library(reshape)
library(lme4)

seed=1
data <- melt(sapply(c(.1,.5,1,5),function(x) rpois(50,x)))

mean(data$value)  #Overall mean

m.glm <- glm(value~1,data=data, family=poisson)
exp(coef(m.glm))     #Exponent of intercept = population mean

m.lmer <- lmer(value~(1|X2),data=data, family=poisson, nAGQ=20)
ranef(m.lmer)$X2+fixef(m.lmer)==coef(m.lmer)
mean(exp(unlist(coef(m.lmer))))   #Mean of exp of coefs = overall mean


-- 
http://mutualism.williams.edu




More information about the R-sig-mixed-models mailing list