[R-sig-ME] AIC and BIC with ML and REML

[Douglas Bates]

On Mon, Mar 30, 2009 at 4:35 PM, Gorjanc Gregor
<Gregor.Gorjanc at bfro.uni-lj.si> wrote:
> Hi!

> A time ago I forgot the details about the calculation about how AIC and BIC are calculated and I did
> a search. All went fine that day, i.e., AIC = - 2 * l + 2 * p, where l is log-likelihood value and p is a number
> of parameters in the model. Similar for BIC. However, I have read a paper today where they state that when
> we fit a model using ML we should use

> AIC = - 2 * l + 2 * (p + q),

> while the following should be used in the case of REML

> AIC = - 2 * l + 2 * q,

> where l is as above, p is a number of free parameters in the fixed part of the model and q is the number of
> variance components for the random part of the model. Is this correct? I checked the calculations
> with lmer() and it seems that only the first approach is taken. Can anyone comment on this?

In a way this question relates to the earlier discussion on whether to
regard the random effects as parameters or as unobserved random
variables.  The difference between -2 * l + 2 * p and -2 * l + 2 * (p
+ q) can be considered to be a question of how many parameters there
are in the model.  In particular, do the random effects count as
parameters?  I had an interesting discussion with Georges Monette
about this a few days ago and both of us feel the saying the random
effects don't affect the parameter count is underestimating the
complexity of the model but saying they should add q to the number of
parameters is overestimating the complexity.

I don't know a good answer to the question of how to "count" the
number of parameters in a mixed model (but I also don't feel that AIC
or BIC should be taken too seriously - these quantities are, at best,
a guide for model comparisons).


> mu <- 10
> a <- c(0, 3)
> b <- rnorm(n=20, sd=3)
>
> facA <- gl(n=2, k=100)
> facB <- gl(n=20, k=10)
>
> y <- rnorm(n=length(facA), mean=(mu + a[facA] + b[facB]), sd=3)
>
> library(package="lme4")
>
> lmer(y ~ facA + (1 | facB))
> ##  AIC  BIC logLik deviance REMLdev
> ## 1054 1067 -522.9     1050    1046
> ##
> ## (1054 - 2 * 522.9) / 2 = 4.1 --> 4 parameters, mu, a2, sigma^2_b, and sigma^2_e --> p = 2, q = 2
>
> lmer(y ~ facA + (1 | facB), REML=FALSE)
> ##   AIC  BIC logLik deviance REMLdev
> ##  1058 1071 -524.9     1050    1046
> ##
> ## (1058 - 2 * 524.9) / 2 = 4.1  --> 4 parameters, mu, a2, sigma^2_b, and sigma^2_e --> p = 2, q = 2
>
> Lep pozdrav / With regards,
>    Gregor Gorjanc
> ----------------------------------------------------------------------
> University of Ljubljana       PhD student
> Biotechnical Faculty          www: http://gregor.gorjanc.googlepages.com
> Department of Animal Science  blog: http://ggorjan.blogspot.com
> Groblje 3                     mail: gregor.gorjanc <at> bfro.uni-lj.si
> SI-1230 Domzale               fax: +386 (0)1 72 17 888
> Slovenia, Europe              tel: +386 (0)1 72 17 861
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>