[R-sig-ME] glmmTMB and ar1

[Jarrod Hadfield]

Hi Ben,

Common.Name is a dummy variable, it only has one level:

https://cran.r-project.org/web/packages/glmmTMB/vignettes/covstruct.html

the model doesn't converge when adding 1|common.Name.

The vignette is a little ambiguous since the process variance is set to 
one in the example, but even then I see no estimate of it in the summary 
even though the surrounding text suggests it should be there.

Cheers,

Jarrod
||
||

||



On 22/11/2017 15:51, Ben Bolker wrote:
>    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
>>
>>
>>
>>
>>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20171122/623a600b/attachment.ksh>