[R-sig-ME] glmmTMB predicted values and AIC, vs glmer
Ben Bolker
bbo|ker @end|ng |rom gm@||@com
Fri May 21 15:27:42 CEST 2021
I don't know yet, I will dig in and try to see what's going on. An
infinite predicted value certainly seems like an issue.
The most obvious difference is that nAGQ=0 is actually doing
something slightly different (it's fitting the fixed-effect parameters
as part of the PIRLS loop rather than maximizing over them explicitly);
I would rather compare the nAGQ=1 case, just to minimize the number of
differences.
On 5/21/21 6:30 AM, John Maindonald wrote:
> The code that follows below demonstrates what I find a very odd
> issue with r-sig-mixed-models using r-project.org<mailto:r-sig-mixed-models using r-project.org>, and with the comparison between
> glmmTMB::glmmTMB and lme4::glmer, for models that have the
> same model formulae.
>
> 1) In the glmmTMB model as fitted to `ffly`, one predicted value is
> infinite. This, even though all coefficients and standard errors
> have finite values. This surely indicates that there is something
> odd in the way that predictions on the scale of the linear
> predictor are calculated.
>
> 2) Initially, I'd thought that the `Inf` value explained the
> difference between the two AIC values. But the difference
> pretty much stays the same when the "offending" point is
> removed.
>
> I've not checked what happens when there are no observation level
> random effects. Is there some difference in the way that the
> model formulae are interpreted in the two cases?
>
>
> ffly <- read.csv('https://github.com/jhmaindonald/dataR/files/6521483/ffly.csv')
> ffly$obs <- factor(ffly$obs)
>
> form1 <- cbind(Dead,Live)~0+trtGp/TrtTime+(1|obs)+(1|trtGpRep)
>
> library(lme4); library(glmmTMB)
> ObsTMB.cll <- glmmTMB(form1,
> family=binomial(link="cloglog"), data=ffly)
> Obsglmer.cll <- glmer(form1, nAGQ=0,
> family=binomial(link="cloglog"), data=ffly)
>
> round(AIC(Obsglmer.cll, ObsTMB.cll), 2)
> ## df AIC
> ## Obsglmer.cll 14 639.44
> ## ObsTMB.cll 14 636.02
>
> round(c(max(predict(ObsTMB.cll)),max(predict(Obsglmer.cll))),3)
> ## [1] Inf 3.833
>
> range(fixef(ObsTMB.cll))
> ## [1] -5.28780064553 1.15952249272
> range(vcov(ObsTMB.cll))
> ## [1] -0.0546986905338 0.2887049942215
>
> ## Try also; there are small deviations from the line
> plot(predict(Obsglmer.cll)~predict(ObsTMB.cll))
> abline(0,1)
>
> ## Remove offending point
> ffly1 <- ffly[-which.max(predict(ObsTMB.cll)),]
> ObsTMB1.cll <- glmmTMB(form1,
> family=binomial(link="cloglog"), data=ffly1)
> Obsglmer1.cll <- glmer(form1, nAGQ=0,
> family=binomial(link="cloglog"), data=ffly1)
>
> cbind(AIC(Obsglmer.cll, ObsTMB.cll), AIC1=AIC(Obsglmer1.cll, ObsTMB1.cll)[,2])
> ## df AIC AIC1
> ## Obsglmer.cll 14 639.441888969 639.441889196
> ## ObsTMB.cll 14 636.016597730 636.016597723
> ## Observe that changes are in the final four decimal places.
>
>
> round(rbind(glmer=fixef(Obsglmer1.cll),
> TMB=unclass(fixef(ObsTMB1.cll)$cond)),3)
>
> trtGpspAEgg trtGpspAL2 trtGpspAL3 trtGpspBEgg trtGpspBL2 trtGpspBL3
> glmer 0.772 -2.255 -3.299 -1.803 -2.632 -5.114
> TMB 0.790 -2.367 -3.441 -1.889 -2.757 -5.288
> trtGpspAEgg:TrtTime trtGpspAL2:TrtTime trtGpspAL3:TrtTime
> glmer 0.231 0.372 0.563
> TMB 0.288 0.398 0.602
> trtGpspBEgg:TrtTime trtGpspBL2:TrtTime trtGpspBL3:TrtTime
> glmer 0.278 0.517 1.112
> TMB 0.299 0.554 1.160
>
>
> John Maindonald email: john.maindonald using anu.edu.au<mailto:john.maindonald using anu.edu.au>
>
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models using 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