R help - Apr 2009 - question related to fitting overdispersion count data using lmer quasipoisson

Dear R-helpers:
I have a question related to fitting overdispersed count data using lmer.
Basically, I simulate an overdispsed data set by adding an observation-level
normal random shock
into exp(....+rnorm()).
Then I fit a lmer quasipoisson model.
The estimation results are very off (see model output of fit.lmer.over.quasi
below).
Can someone kindly explain to me what went wrong?

Many thanks in advance.

-Sean

#data simulation (modified from code at
http://markmail.org/message/j3zmgrklihe73p4p)
set.seed(100)
m <- 5
n <- 100
N <- n*m
#X <- cbind(1,runif(N))
X <- cbind(1,rnorm(N))
X <- cbind(runif(N),rnorm(N))
id <- rep(1:n,each=m)
#
Z <- kronecker(diag(n),rep(1,m))
#Possion with group level heterogeneity
z <- rpois(N, exp(X%*%matrix(c(1,2)) + Z%*%matrix(rnorm(n))))
#2*rnorm(n*m) is added to each observation to create overdispersion
z.overdis <- rpois(N, exp(X%*%matrix(c(1,2)) + Z%*%matrix(rnorm(n)) +
2*rnorm(n*m)))

#without observation-level random shock i.e., 2*rnorm(n*m), estimate results
are very accurate
(fit.lmer <- lmer(z ~ X + (1|id), family="poisson",verbose=F))
 #Generalized linear mixed model fit by the Laplace approximation
#Formula: z ~ X + (1 | id)
# AIC BIC logLik deviance
# 851 868   -422      843
#Random effects:
# Groups Name        Variance Std.Dev.
# id     (Intercept) 0.977    0.988
#Number of obs: 500, groups: id, 100
#
#Fixed effects:
#            Estimate Std. Error z value Pr(>|z|)
#(Intercept)  -0.0128     0.1116    -0.1      0.9
#X1            1.0615     0.0601    17.7   <2e-16 ***
#X2            2.0236     0.0214    94.7   <2e-16 ***
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
#Correlation of Fixed Effects:
#   (Intr) X1
#X1 -0.349
#X2 -0.270  0.258




#Now you can see the results are very off
(fit.lmer.over.quasi <- lmer(z.overdis ~ X + (1|id),
family=quasipoisson(link="log"),verbose=F))
#Generalized linear mixed model fit by the Laplace approximation
#Formula: z.overdis ~ X + (1 | id)
#   AIC   BIC logLik deviance
# 41867 41888 -20929    41857
#Random effects:
# Groups   Name        Variance Std.Dev.
# id       (Intercept) 175.8    13.26
# Residual              72.9     8.54
#Number of obs: 500, groups: id, 100
#
#Fixed effects:
#            Estimate Std. Error t value
#(Intercept)   1.3530     1.3492    1.00
#X1            1.0834     0.2273    4.77
#X2            1.3501     0.0783   17.25
#
#Correlation of Fixed Effects:
#   (Intr) X1
#X1 -0.099
#X2 -0.055  0.070

	[[alternative HTML version deleted]]