[R-sig-ME] generalized linear mixed models: large differences when using glmmPQL or lmer with laplace approximation
Ken Beath
ken at kjbeath.com.au
Tue Oct 7 23:07:50 CEST 2008
There is a large difference between the estimated std of the random
effect, usually a sign that the glmmPQL approximation isn't working.
I would try using adaptive Gauss-Hermite. Set the nAGQ parameter to
increasing values until the results stabilise.
Ken
On 08/10/2008, at 12:21 AM, Martijn Vandegehuchte wrote:
> Dear list,
>
> First of all, I am a mere ecologist, trying to get the truth out of
> his data, and not a statistician, so forgive me my lack of
> statistical background and possible conceptual misunderstandings.
>
> I am currently comparing generalized linear mixed models in glmmPQL
> and lmer, with a quasipoisson family, and have found out that
> parameter estimates are quite different for both methods. I read
> some of the discussions on the R-forum and it seems that the Laplace
> approximation used in the current version of lmer is generally
> preferred to the PQL method. I am an ex-SAS user, and in proc
> glimmix in SAS the default is PQL, and the estimates and p-values
> are almost exact the same as with glmmPQL in R. But lmer gives quite
> different results, and now I am wondering what would be the best
> option for me.
>
> First of all, parameter estimates of a same model can be somewhat
> different in lmer or glmmPQL. Second of all, in lmer, I only get t-
> values but no associated p-values (apparently they are omitted
> because of the uncertainty about the df). But if I compare the t-
> values generated by glmmPQL with those of a same model in lmer, the
> differences are substantial. My dataset consists of 120
> observations, so basically you could guess the order of magnitude of
> the p-values in lmer based on the t-value and a "large" df.
>
> First example:
> In lmer:
>
>> model<-lmer(schirufu~diameter+leafvit+densroot+cover+nemcm+(1|
>> site),family=quasipoisson)
>> summary (model)
> Generalized linear mixed model fit by the Laplace approximation
> Formula: schirufu ~ diameter + leafvit + densroot + cover + nemcm +
> (1 | site)
> AIC BIC logLik deviance
> 2045 2068 -1015 2029
> Random effects:
> Groups Name Variance Std.Dev.
> site (Intercept) 12.700 3.5638
> Residual 15.182 3.8964
> Number of obs: 120, groups: site, 6
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 1.31017 1.47249 0.890
> diameter -0.24799 0.29180 -0.850
> leafvit 1.29007 0.21041 6.131
> densroot 0.31024 0.04939 6.281
> cover -0.24544 0.22179 -1.107
> nemcm 0.24817 0.12028 2.063
>
> Correlation of Fixed Effects:
> (Intr) diamtr leafvt densrt cover
> diameter 0.031
> leafvit -0.083 0.321
> densroot 0.011 -0.017 -0.202
> cover 0.021 -0.448 0.016 0.214
> nemcm -0.014 0.114 0.114 0.310 -0.017
>>
>
> Although no p-values are given, it suggests that fixed effects
> leafvit, densroot and nemcm would be significant.
> In glmmPQL:
>
>> model<-glmmPQL(schirufu~diameter+leafvit+densroot+cover
>> +nemcm,random=~1|site,family=quasipoisson)
> iteration 1
> iteration 2
> iteration 3
> iteration 4
> iteration 5
>> summary(model)
> Linear mixed-effects model fit by maximum likelihood
> Data: NULL
> AIC BIC logLik
> NA NA NA
>
> Random effects:
> Formula: ~1 | site
> (Intercept) Residual
> StdDev: 0.7864989 4.63591
>
> Variance function:
> Structure: fixed weights
> Formula: ~invwt
> Fixed effects: schirufu ~ diameter + leafvit + densroot + cover +
> nemcm
> Value Std.Error DF t-value p-value
> (Intercept) 1.4486735 0.4174843 109 3.470007 0.0007
> diameter -0.2600504 0.3477017 109 -0.747913 0.4561
> leafvit 1.2236406 0.2489291 109 4.915619 0.0000
> densroot 0.3236446 0.0596342 109 5.427164 0.0000
> cover -0.2523163 0.2698555 109 -0.935005 0.3519
> nemcm 0.2336305 0.1451751 109 1.609301 0.1104
> Correlation:
> (Intr) diamtr leafvt densrt cover
> diameter 0.130
> leafvit -0.335 0.313
> densroot 0.027 -0.022 -0.203
> cover 0.090 -0.463 0.015 0.214
> nemcm -0.056 0.097 0.107 0.301 -0.014
>
> Standardized Within-Group Residuals:
> Min Q1 Med Q3 Max
> -2.4956188 -0.4154369 -0.1333850 0.1724601 4.7355928
>
> Number of Observations: 120
> Number of Groups: 6
>>
>
> Note the difference in parameter estimates. Also, the fixed effect
> nemcm now is not significant any more.
>
> Second example,now with an offset:
> In lmer:
>
>> model<-lmer(nemcm~diameter+leafvit+densroot+rootvit+cover+schirufu
>> +(1|site), offset= loglength, family=quasipoisson)
>> summary (model)
> Generalized linear mixed model fit by the Laplace approximation
> Formula: nemcm ~ diameter + leafvit + densroot + rootvit + cover +
> schirufu + (1 | site)
> AIC BIC logLik deviance
> 1593 1618 -787.4 1575
> Random effects:
> Groups Name Variance Std.Dev.
> site (Intercept) 21.522 4.6392
> Residual 173.888 13.1867
> Number of obs: 120, groups: site, 6
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 0.06733 1.92761 0.0349
> diameter 0.14665 0.60693 0.2416
> leafvit -0.19902 0.48802 -0.4078
> densroot -0.49178 0.64221 -0.7658
> rootvit 0.37699 0.46810 0.8054
> cover -0.23545 0.57896 -0.4067
> schirufu 0.23226 0.46866 0.4956
>
> Correlation of Fixed Effects:
> (Intr) diamtr leafvt densrt rootvt cover
> diameter -0.016
> leafvit 0.015 0.396
> densroot 0.055 -0.233 -0.291
> rootvit -0.038 -0.251 -0.629 0.277
> cover 0.024 -0.796 -0.133 0.253 0.117
> schirufu -0.032 0.137 -0.029 -0.505 -0.078 -0.121
>>
>
> This suggests no significant effects at all.
> In glmmPQL:
>
>> model<-glmmPQL(nemcm~diameter+leafvit+densroot+rootvit+cover
>> +schirufu+offset(loglength),random=~1|site, family=quasipoisson)
> iteration 1
> iteration 2
> iteration 3
>> summary (model)
> Linear mixed-effects model fit by maximum likelihood
> Data: NULL
> AIC BIC logLik
> NA NA NA
>
> Random effects:
> Formula: ~1 | site
> (Intercept) Residual
> StdDev: 0.2684477 4.507758
>
> Variance function:
> Structure: fixed weights
> Formula: ~invwt
> Fixed effects: nemcm ~ diameter + leafvit + densroot + rootvit +
> cover + schirufu + offset(loglength)
> Value Std.Error DF t-value p-value
> (Intercept) 0.1131898 0.1656949 108 0.6831220 0.4960
> diameter 0.1225231 0.1976568 108 0.6198779 0.5366
> leafvit -0.2191361 0.1697784 108 -1.2907181 0.1996
> densroot -0.4733839 0.2221562 108 -2.1308604 0.0354
> rootvit 0.3858120 0.1615706 108 2.3878846 0.0187
> cover -0.2075038 0.1922054 108 -1.0795940 0.2827
> schirufu 0.2028444 0.1633954 108 1.2414323 0.2171
> Correlation:
> (Intr) diamtr leafvt densrt rootvt cover
> diameter -0.050
> leafvit 0.077 0.360
> densroot 0.217 -0.168 -0.262
> rootvit -0.163 -0.202 -0.632 0.257
> cover 0.084 -0.772 -0.098 0.200 0.073
> schirufu -0.103 0.099 -0.050 -0.483 -0.068 -0.075
>
> Standardized Within-Group Residuals:
> Min Q1 Med Q3 Max
> -1.1146287 -0.5208003 -0.1927005 0.2462878 7.9755368
>
> Number of Observations: 120
> Number of Groups: 6
>>
>
> Again some differences in parameter estimates, but now the two fixed
> effects densroot and rootvit turn out to be significant.
> So my questions are:
> - what would you recommend me to use? lmer or glmmPQL (laplace
> approximation or penalized quasi-likelihood)?
> - if lmer is the better option, is there a way to get a reliable p-
> value for the fixed effects?
> I have experienced that deleting a term and comparing models using
> anova() always overestimates the significance of that term, probably
> because the quasipoisson correction for overdispersion is not taken
> into account.
>
> Thank you very much beforehand,
>
> Martijn.
>
> --
> Martijn Vandegehuchte
> Ghent University
> Department Biology
> Terrestrial Ecology Unit
> K.L.Ledeganckstraat 35
> B-9000 Ghent
> telephone: +32 (0)9/264 50 84
> e-mail: martijn.vandegehuchte at ugent.be
>
> website TEREC: www.ecology.ugent.be/terec
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
More information about the R-sig-mixed-models
mailing list