[R-sig-ME] likelihood-ratio tests in conflict with coefficiants in maximal random effect model
Levy, Roger
rlevy at ucsd.edu
Wed Mar 5 17:25:58 CET 2014
On Mar 5, 2014, at 3:15 AM, Emilia Ellsiepen <emilia.ellsiepen at gmail.com> wrote:
> 2014-03-04 18:39 GMT+01:00 Levy, Roger <rlevy at ucsd.edu>:
>>
>>
>>
>> On Mar 4, 2014, at 3:10 AM, Emilia Ellsiepen <emilia.ellsiepen at gmail.com> wrote:
>>
>>> Thank you for your support, Roger - actually the problem vanished when
>>> I used sum-coding (after installing the packages and reloading the
>>> data, I accidentally skipped this step yesterday). So, in fact with
>>> the old lme4.0 version, the model did converge. I pasted in all model
>>> summaries and LRTs for dummy and sum coding below.
>>> Would you recommend to always use the old version then?
>>
>> Thanks so much, Emilia, for following up. This is very interesting. Something doesn't quite add up, though, in the model summaries and LRTs that you pasted below, because the log-likelihoods in the model summaries are not the same as the log-likelihoods reported in the results of anova(). I notice that you're not using REML=F as an argument in your call to lmer(). Could you perhaps rerun all four model-fitting functions and the two anova() calls with REML=F, and show us the complete results?
>>
>> Best & thanks again,
>>
>> Roger
>>
> Rerunning with REML=F led to a singular convergence in the
> null-hypothesis model with dummy coding. The sum coded variant still
> converges without warnings. While the two LRTs give in general the
> same result for both coding schemes, as I would have expected, there
> are minimal differences in the logLik, AIC and BIC values.
Thanks so much once again, Emilia, for continuing to follow up. These results clear everything up. As you note, there are now only small differences between the log-likelihoods of equivalent models with the two different contrast codings. (You could probably get the models to agree even more by tightening the tolerance criterion for convergence.) Note also that the likelihood-ratio test under either parameterization gives you a result highly consistent with the t-statistic results for the maximal model, addressing one of the concerns that you originally voiced. I don’t see any evidence that the maximal model is overparameterized (though I would be happy to hear any arguments to the contrary from other list users!). So I would trust the LRT results you’re seeing now — that the interaction you care about is highly significant.
It would be worthwhile to check whether you get the same results running the same code using the newest version of lme4 with the nlminb optimizer — and I’m sure a number of list members would be interested in knowing!
Best and I’m glad that we’ve gotten to the bottom of things!
Roger
> Here are the results in detail:
> **** dummy coding:
>
>> m1 = lmer(score ~ Order*Voice + (Order*Voice|subject)+(Order*Voice|sentence), REML=F ,rawdata)
>> summary(m1)
> Linear mixed model fit by maximum likelihood
> Formula: score ~ Order * Voice + (Order * Voice | subject) + (Order *
> Voice | sentence)
> Data: rawdata
> AIC BIC logLik deviance REMLdev
> 1612 1726 -780.8 1562 1576
> Random effects:
> Groups Name Variance Std.Dev. Corr
> subject (Intercept) 0.092788 0.30461
> OrderDat second 0.051950 0.22793 -0.689
> VoicePassive 0.152335 0.39030 -0.954 0.875
> OrderDat second:VoicePassive 0.171303 0.41389 0.736 -0.991 -0.903
> sentence (Intercept) 0.013700 0.11705
> OrderDat second 0.060511 0.24599 -0.652
> VoicePassive 0.021227 0.14570 -0.215 -0.067
> OrderDat second:VoicePassive 0.045796 0.21400 0.016 -0.587 -0.404
> Residual 0.466165 0.68276
> Number of obs: 720, groups: subject, 36; sentence, 20
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 0.16341 0.07650 2.136
> OrderDat second -0.01835 0.09823 -0.187
> VoicePassive -0.05604 0.10234 -0.548
> OrderDat second:VoicePassive 0.54250 0.13194 4.112
>
> Correlation of Fixed Effects:
> (Intr) OrdrDs VcPssv
> OrderDtscnd -0.646
> VoicePassiv -0.757 0.461
> OrdrDscn:VP 0.514 -0.719 -0.730
>
>> mi = lmer(score ~ Order+Voice + (Order*Voice|subject)+(Order*Voice|sentence) ,REML=F,rawdata)
> Warning message:
> In mer_finalize(ans) : singular convergence (7)
>> summary(mi)
> Linear mixed model fit by maximum likelihood
> Formula: score ~ Order + Voice + (Order * Voice | subject) + (Order *
> Voice | sentence)
> Data: rawdata
> AIC BIC logLik deviance REMLdev
> 1622 1732 -787.2 1574 1586
> Random effects:
> Groups Name Variance Std.Dev. Corr
> subject (Intercept) 0.112099 0.33481
> OrderDat second 0.095834 0.30957 -0.743
> VoicePassive 0.206993 0.45497 -0.953 0.895
> OrderDat second:VoicePassive 0.333823 0.57777 0.765 -0.999 -0.906
> sentence (Intercept) 0.017953 0.13399
> OrderDat second 0.105527 0.32485 -0.832
> VoicePassive 0.046129 0.21478 -0.473 0.454
> OrderDat second:VoicePassive 0.183020 0.42781 0.526 -0.786 -0.819
> Residual 0.463528 0.68083
> Number of obs: 720, groups: subject, 36; sentence, 20
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 0.008662 0.065187 0.133
> OrderDat second 0.275830 0.068133 4.048
> VoicePassive 0.237267 0.069355 3.421
>
> Correlation of Fixed Effects:
> (Intr) OrdrDs
> OrderDtscnd -0.468
> VoicePassiv -0.632 -0.161
>> anova(m1,mi)
> Data: rawdata
> Models:
> mi: score ~ Order + Voice + (Order * Voice | subject) + (Order *
> mi: Voice | sentence)
> m1: score ~ Order * Voice + (Order * Voice | subject) + (Order *
> m1: Voice | sentence)
> Df AIC BIC logLik Chisq Chi Df Pr(>Chisq)
> mi 24 1622.4 1732.3 -787.21
> m1 25 1611.6 1726.1 -780.82 12.788 1 0.0003489 ***
> ---
> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
>
> **** sum coding:
>
>> contrasts(rawdata$Order) = contr.sum(2)
>> contrasts(rawdata$Voice) = contr.sum(2)
>> m1 = lmer(score ~ Order*Voice + (Order*Voice|subject)+(Order*Voice|sentence), REML=F ,rawdata)
>> summary(m1)
> Linear mixed model fit by maximum likelihood
> Formula: score ~ Order * Voice + (Order * Voice | subject) + (Order *
> Voice | sentence)
> Data: rawdata
> AIC BIC logLik deviance REMLdev
> 1612 1726 -780.9 1562 1581
> Random effects:
> Groups Name Variance Std.Dev. Corr
> subject (Intercept) 1.8907e-02 0.1375028
> Order1 1.2567e-05 0.0035449 -1.000
> Voice1 1.2530e-02 0.1119373 0.866 -0.866
> Order1:Voice1 1.0808e-02 0.1039636 0.292 -0.292 0.630
> sentence (Intercept) 2.7251e-03 0.0522022
> Order1 1.0143e-02 0.1007109 0.250
> Voice1 5.1207e-03 0.0715590 0.333 -0.465
> Order1:Voice1 2.9324e-03 0.0541514 -0.877 0.186 -0.353
> Residual 4.6608e-01 0.6826980
> Number of obs: 720, groups: subject, 36; sentence, 20
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 0.26184 0.03618 7.238
> Order1 -0.12645 0.03398 -3.721
> Voice1 -0.10761 0.03538 -3.042
> Order1:Voice1 0.13562 0.03308 4.100
>
> Correlation of Fixed Effects:
> (Intr) Order1 Voice1
> Order1 0.042
> Voice1 0.338 -0.147
> Order1:Voc1 -0.007 0.043 0.116
>> mi = lmer(score ~ Order+Voice + (Order*Voice|subject)+(Order*Voice|sentence) ,REML=F,rawdata)
>> summary(mi)
> Linear mixed model fit by maximum likelihood
> Formula: score ~ Order + Voice + (Order * Voice | subject) + (Order *
> Voice | sentence)
> Data: rawdata
> AIC BIC logLik deviance REMLdev
> 1623 1732 -787.3 1575 1589
> Random effects:
> Groups Name Variance Std.Dev. Corr
> subject (Intercept) 1.9424e-02 0.1393698
> Order1 2.6892e-06 0.0016399 -1.000
> Voice1 1.3041e-02 0.1141955 0.838 -0.838
> Order1:Voice1 2.0839e-02 0.1443579 0.290 -0.290 0.541
> sentence (Intercept) 2.4074e-03 0.0490650
> Order1 1.0471e-02 0.1023269 0.323
> Voice1 4.3300e-03 0.0658028 0.280 -0.571
> Order1:Voice1 1.1412e-02 0.1068262 -0.784 0.204 -0.300
> Residual 4.6362e-01 0.6808933
> Number of obs: 720, groups: subject, 36; sentence, 20
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 0.26513 0.03609 7.346
> Order1 -0.13474 0.03407 -3.955
> Voice1 -0.11837 0.03480 -3.401
>
> Correlation of Fixed Effects:
> (Intr) Order1
> Order1 0.063
> Voice1 0.334 -0.173
>
>> anova(m1,mi)
> Data: rawdata
> Models:
> mi: score ~ Order + Voice + (Order * Voice | subject) + (Order *
> mi: Voice | sentence)
> m1: score ~ Order * Voice + (Order * Voice | subject) + (Order *
> m1: Voice | sentence)
> Df AIC BIC logLik Chisq Chi Df Pr(>Chisq)
> mi 24 1622.5 1732.4 -787.26
> m1 25 1611.7 1726.2 -780.86 12.796 1 0.0003473 ***
> ---
> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
More information about the R-sig-mixed-models
mailing list