[R] AICs from lmer different with summary and anova

Michael Dewey info at aghmed.fsnet.co.uk
Sun Apr 19 17:36:11 CEST 2009


At 16:22 15/04/2009, Jonathan Williams wrote:
>Dear R Helpers,
>
>I have noticed that when I use lmer to analyse data, the summary function
>gives different values for the AIC, BIC and log-likelihood compared with the
>anova function.

I do not think I have seen a reply to this.

What happens if you fit the models originally using ML rather than REML?


>Here is a sample program
>
>#make some data
>set.seed(1);
>datx=data.frame(array(runif(720),c(240,3),dimnames=list(NULL,c('x1','x2','y'
>))))
>id=rep(1:120,2); datx=cbind(id,datx)
>
>#give x1 a slight relation with y (only necessary to make the random effects
>non-zero in this artificial example)
>datx$x1=(datx$y*0.1)+datx$x1
>
>library(lme4)
>
>#fit the data
>fit0=lmer(y~x1+x2+(1|id), data=datx); print(summary(fit0),corr=F)
>fit1=lmer(y~x1+x2+(1+x1|id), data=datx); print(summary(fit1),corr=F)
>
>#compare the models
>anova(fit0,fit1)
>
>
>Now, look at the output, below. You can see that the AIC from
>"print(summary(fit0))" is 87.34, but the AIC for fit0 in "anova(fit0,fit1)"
>is 73.965. There are similar changes for the values of BIC and logLik.
>
>Am I doing something wrong, here? If not, which are the real AIC and logLik
>values for the different models?
>
>Thanks for your help,
>
>Jonathan Williams
>
>
>Output:-
>
> > fit0=lmer(y~x1+x2+(1|id), data=datx); print(summary(fit0),corr=F)
>Linear mixed model fit by REML
>Formula: y ~ x1 + x2 + (1 | id)
>    Data: datx
>    AIC   BIC logLik deviance REMLdev
>  87.34 104.7 -38.67    63.96   77.34
>Random effects:
>  Groups   Name        Variance Std.Dev.
>  id       (Intercept) 0.016314 0.12773
>  Residual             0.062786 0.25057
>Number of obs: 240, groups: id, 120
>
>Fixed effects:
>             Estimate Std. Error t value
>(Intercept)  0.50376    0.05219   9.652
>x1           0.08979    0.06614   1.358
>x2          -0.06650    0.06056  -1.098
> > fit1=lmer(y~x1+x2+(1+x1|id), data=datx); print(summary(fit1),corr=F)
>Linear mixed model fit by REML
>Formula: y ~ x1 + x2 + (1 + x1 | id)
>    Data: datx
>    AIC   BIC logLik deviance REMLdev
>  90.56 114.9 -38.28    63.18   76.56
>Random effects:
>  Groups   Name        Variance  Std.Dev. Corr
>  id       (Intercept) 0.0076708 0.087583
>           x1          0.0056777 0.075351 1.000
>  Residual             0.0618464 0.248689
>Number of obs: 240, groups: id, 120
>
>Fixed effects:
>             Estimate Std. Error t value
>(Intercept)  0.50078    0.05092   9.835
>x1           0.09236    0.06612   1.397
>x2          -0.06515    0.06044  -1.078
> > anova(fit0,fit1)
>Data: datx
>Models:
>fit0: y ~ x1 + x2 + (1 | id)
>fit1: y ~ x1 + x2 + (1 + x1 | id)
>      Df     AIC     BIC  logLik  Chisq Chi Df Pr(>Chisq)
>fit0  5  73.965  91.368 -31.982
>fit1  7  77.181 101.545 -31.590 0.7839      2     0.6757

Michael Dewey
http://www.aghmed.fsnet.co.uk




More information about the R-help mailing list