[R-sig-ME] different results across versions for glmer/lmer with the quasi-poisson or quasi-binomial families: the lattest version might not be accurate...
Gosselin Frederic
frederic.gosselin at cemagref.fr
Fri Nov 20 16:16:17 CET 2009
************** Slightly modified version of a mail Already posted in the R-help list
************* but Ben Bolker suggested me to also post it here
************* Sorry for cross-posting ******************
Dear R-helpers,
this mail is intended to mention a rather strange result and generate potential useful comments on it. I am not aware of another posts on this issue ( RSiteSearch("quasipoisson lmer version dispersion")). (But see http://r-forge.r-project.org/tracker/index.php?func=detail&aid=307&group_id=60&atid=298).
Using the exemple in the reference of the lmer function (in lme4 library) and turning it into a quasi-poisson or quasi-binomial analysis, we get different results across lme4 and R versions, especially for the dispersion parameter (section 1 in the following). For the quasi-poisson, we found even stronger differences with another, bigger database on which we worked (residual variance from approx. 0.7 to 1.2). We also tested the differences between versions for plain likelihood families (poisson and binomial): the differences were slighter, and were actually very small in the case of the poisson family.
My first reaction was that the later version should have been better. However, the quasibinomial results of this version are rather strange. Therefore, I simulated a bigger data set corresponding to Poisson and binomial mixed models (section 2 in the following). I actually found that the later version of lme4 was suspect but not the old one. My temptative conclusion is therefore to use the old versions of lme4 (here: lme4 version 0.99875-9) when using quasi-binomial and quasi-poisson methods.
Any comment/insight appreciated.
All the best,
Frédéric Gosselin
Engineer & Researcher (PhD) in Forest Ecology
Cemagref
Domaine des Barres
F-45290 Nogent sur Vernisson
France
############### 1- EXAMPLE FROM LMER HELP
############# here are the commands for the quasipoisson:
library(lme4)
gm1 <- lmer(incidence ~ period + (1 | herd), family = quasipoisson, data = cbpp)
summary(gm1)
########## here is the result under R2.5.1 and Package lme4 version 0.99875-9:
Generalized linear mixed model fit using Laplace
Formula: incidence ~ period + (1 | herd)
Data: cbpp
Family: quasipoisson(log link)
AIC BIC logLik deviance
112.2 122.3 -51.11 102.2
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.35085 0.59233
Residual 1.40470 1.18520
number of obs: 56, groups: herd, 15
Fixed effects:
Estimate Std. Error t value
(Intercept) 1.2812 0.2200 5.824
period2 -1.1240 0.3315 -3.391
period3 -1.3203 0.3579 -3.689
period4 -1.9477 0.4808 -4.051
Correlation of Fixed Effects:
(Intr) perid2 perid3
period2 -0.339
period3 -0.314 0.219
period4 -0.233 0.163 0.151
#############" here is the result on R2.7.1 and Package lme4 version 0.999375-26
Generalized linear mixed model fit by the Laplace approximation
Formula: incidence ~ period + (1 | herd)
Data: cbpp
AIC BIC logLik deviance
114.2 126.4 -51.1 102.2
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.32421 0.56939
Residual 1.29474 1.13786
Number of obs: 56, groups: herd, 15
Fixed effects:
Estimate Std. Error t value
(Intercept) 1.2762 0.2115 6.035
period2 -1.1249 0.3187 -3.530
period3 -1.3190 0.3438 -3.837
period4 -1.9450 0.4615 -4.215
Correlation of Fixed Effects:
(Intr) perid2 perid3
period2 -0.339
period3 -0.314 0.219
period4 -0.233 0.163 0.151
############# now the commands for the quasibinomial:
library(lme4)
toto<-as.double(cbpp$incidence>0)
gm2 <- lmer(toto ~ period + (1 | herd), family = quasibinomial, data = cbpp)
summary(gm2)
########## here is the result under R2.5.1 and Package lme4 version 0.99875-9:
Generalized linear mixed model fit using Laplace
Formula: toto ~ period + (1 | herd)
Data: cbpp
Family: quasibinomial(logit link)
AIC BIC logLik deviance
72.04 82.17 -31.02 62.04
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 5.0259e-10 2.2418e-05
Residual 1.0052e+00 1.0026e+00
number of obs: 56, groups: herd, 15
Fixed effects:
Estimate Std. Error t value
(Intercept) 2.662 1.048 2.540
period2 -2.078 1.188 -1.750
period3 -2.950 1.180 -2.501
period4 -3.135 1.194 -2.626
Correlation of Fixed Effects:
(Intr) perid2 perid3
period2 -0.882
period3 -0.888 0.784
period4 -0.878 0.775 0.780
############# here is the result on R2.7.1 and Package lme4 version 0.999375-26
Generalized linear mixed model fit by the Laplace approximation
Formula: toto ~ period + (1 | herd)
Data: cbpp
AIC BIC logLik deviance
74.04 86.2 -31.02 62.04
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.00000 0.00000
Residual 0.13522 0.36772
Number of obs: 56, groups: herd, 15
Fixed effects:
Estimate Std. Error t value
(Intercept) 2.6391 0.3806 6.933
period2 -2.0513 0.4324 -4.744
period3 -2.9267 0.4293 -6.817
period4 -3.1091 0.4345 -7.155
Correlation of Fixed Effects:
(Intr) perid2 perid3
period2 -0.880
period3 -0.887 0.780
period4 -0.876 0.771 0.777
############### 2- SIMULATED EXAMPLE
############# here are the commands both for the quasibinomial and quasipoisson:
library(lme4)
set.seed(1)
period<-rnorm(1000)
herd<-rep(1:50,eac=20)
herd.effect<-rnorm(50)
toto<-rbinom(1000,1,exp(period+herd.effect[herd])/(1+exp(period+herd.effect[herd])))
gm2s <- lmer(toto ~ period + (1 | herd), family = quasibinomial)
summary(gm2s)
gm2sL <- lmer(toto ~ period + (1 | herd), family = binomial)
summary(gm2sL)
#now poisson with similar data:
#essai avec jeu de données simulé:
set.seed(2)
toto<-rpois(1000,exp(period+herd.effect[herd]))
gm1s <- lmer(toto ~ period + (1 | herd), family = quasipoisson)
summary(gm1s)
gm1sL <- lmer(toto ~ period + (1 | herd), family = poisson)
summary(gm1sL)
########## here is the quasi-poisson result under R2.5.1 and Package lme4 version 0.99875-9: (all is OK)
Generalized linear mixed model fit using Laplace
Formula: toto ~ period + (1 | herd)
Family: quasipoisson(log link)
AIC BIC logLik deviance
1151 1166 -572.7 1145
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.83605 0.91436
Residual 0.97980 0.98985
number of obs: 1000, groups: herd, 50
Fixed effects:
Estimate Std. Error t value
(Intercept) 0.02775 0.13419 0.21
period 0.97424 0.02291 42.52
Correlation of Fixed Effects:
(Intr)
period -0.157
############# here is the quasipoisson result on R2.7.1 and Package lme4 version 0.999375-26: the variance estimates are suspect
Generalized linear mixed model fit by the Laplace approximation
Formula: toto ~ period + (1 | herd)
AIC BIC logLik deviance
1153 1173 -572.7 1145
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.35305 0.59418
Residual 0.41370 0.64320
Number of obs: 1000, groups: herd, 50
Fixed effects:
Estimate Std. Error t value
(Intercept) 0.02757 0.08720 0.32
period 0.97423 0.01489 65.44
Correlation of Fixed Effects:
(Intr)
period -0.157
########## here is the quasi-binomial result under R2.5.1 and Package lme4 version 0.99875-9: (all is OK)
Generalized linear mixed model fit using Laplace
Formula: toto ~ period + (1 | herd)
Family: quasibinomial(logit link)
AIC BIC logLik deviance
1145 1160 -569.6 1139
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.85697 0.92573
Residual 0.94584 0.97254
number of obs: 1000, groups: herd, 50
Fixed effects:
Estimate Std. Error t value
(Intercept) -0.05493 0.14949 -0.367
period 1.02879 0.08360 12.306
Correlation of Fixed Effects:
(Intr)
period -0.001
############# here is the quasi-binomial result on R2.7.1 and Package lme4 version 0.999375-26: the variance estimates are very suspect
Generalized linear mixed model fit by the Laplace approximation
Formula: toto ~ period + (1 | herd)
AIC BIC logLik deviance
1147 1167 -569.6 1139
Random effects:
Groups Name Variance Std.Dev.
herd (Intercept) 0.12520 0.35384
Residual 0.13819 0.37173
Number of obs: 1000, groups: herd, 50
Fixed effects:
Estimate Std. Error t value
(Intercept) -0.05503 0.05714 -0.96
period 1.03117 0.03199 32.24
Correlation of Fixed Effects:
(Intr)
More information about the R-sig-mixed-models
mailing list