glmmTMB and ar1
Without looking too closely at this, I think you also need to include (1|Common.Name) in the model; otherwise the assumption is that there is no within-group variance? This formula is taken from http://bbolker.github.io/mixedmodels-misc/notes/corr_braindump.html glmmTMB_simple_fit <- glmmTMB(y~1 + (1|f) + ar1(tt-1|f), data=d,family=gaussian) (you used +0 rather than -1 to suppress the intercept in the ar1() term, and you have a non-trivial fixed-effect term, but otherwise these are similar). Suggestions for doc improvements/pull requests welcome ... Ben Bolker
On 17-11-22 05:26 AM, Jarrod Hadfield wrote:
Hi, It is *really* great that glmmTMB allows ar1 models. However, I'm having some trouble understanding the output and reconciling the estimates with asreml. The data consist of the number of birds censused each year for 34 years. In 13 years the birds were censused twice. The model I would like to fit has year as a continuous fixed effect, and then an ar1 process across years. The residual variance should pick up the within-year variance. ?m1.glmmTMB<-glmmTMB(log(pop)~year+ar1(year.factor+0|Common.Name), data=shag_data) However, this gives one fewer parameters than I was expecting: ?summary( m1.glmmTMB) ?Family: gaussian? ( identity ) Formula:????????? log(pop) ~ year + ar1(year.factor + 0 | Common.Name) Data: shag_data ???? AIC????? BIC?? logLik deviance df.resid ??? 12.5???? 21.7???? -1.2????? 2.5?????? 42 Random effects: Conditional model: ?Groups????? Name??????????? Variance Std.Dev. Corr ?Common.Name year.factor1973 0.214434 0.46307?? (ar1) ?Residual??????????????????? 0.004099 0.06403 Number of obs: 47, groups:? Common.Name, 1 Dispersion estimate for gaussian family (sigma^2): 0.0041 Conditional model: ??????????? Estimate Std. Error z value Pr(>|z|) (Intercept) 55.73041?? 29.70301?? 1.876?? 0.0606 . year??????? -0.02464??? 0.01493? -1.650?? 0.0989 . --- Signif. codes:? 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 Is 0.214434 the process variance for the ar1?? But then where is the autocorrelation parameter? What I hoped was the equivalent model in asreml gives different answers m1.asreml<-asreml(log(pop)~year, random=~ar1v(year.factor), data=shag_data) The estimate are 0.74 (autocorrelation), 0.30 (process variance) and 0.0041 (the residual variance). Asreml uses REML not ML so this might explain some of the discrepancy but I'd be surprised if it explained all. Cheers, Jarrod