[R] How can I extract the AIC score from a mixed model object produced using lmer?

Henric Nilsson (Public) nilsson.henric at gmail.com
Sun Dec 30 02:32:36 CET 2007 

Douglas Bates wrote:
> On Dec 19, 2007 9:42 AM, David Hewitt <dhewitt at vims.edu> wrote:
>>
>> David Barron-3 wrote:
>>> You can calculate the AIC as follows:
>>>
>>> (fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
>>> aic1 <- AIC(logLik(fm1))
>>>
>>>
> 
>> Is AIC() [extractAIC()] "valid" for models with random effects? I noticed
>> that the help page for extractAIC() does not list models with random
>> effects. I think this boils down to the difference between the likelihoods
>> for models with and without random effects, and I don't know. Just
>> curious...
> 
> The log-likelihood for a linear mixed model is well-defined.  Whether
> this makes AIC valid or not depends on how comfortable you are with
> the idea of AIC in the first place.  My impression is that the
> justification for AIC is not entirely rigorous but I must admit that I
> haven't gone back to look at the original literature on it.
> 
> To the best of my knowledge and ability the log-likelihood from a
> model fit by lmer with method = "ML" is properly defined and
> accurately evaluated.  (The default estimation criterion in lmer is
> REML and models fit by REML provide a close approximation to the
> log-likelihood but not an exact result.  If you really want a
> log-likelihood and AIC value you should refit with method = "ML".)
> What is later done to the log-likelihood to obtain the AIC value is
> more problematic.  In particular, one needs to provide a value for the
> number of parameters in the model and that can be tricky.  Recently I
> was working with models for data on test scores by students over time.
>  There were over 200,000 students.  Under one way of counting
> parameters, if I incorporate a random effect for the student this
> costs me only 1 parameter, corresponding to the variance component for
> that random effect.  However, I am incorporating over 200,000 random
> effects to help model the observed responses.  So is the number of
> parameters 1 or over 200,000?  I don't know.

I can't remember the exact details, but I do remember that the issue is 
discussed in

@ARTICLE{LMM:Vaida+Blanchard:2005,
   author = {Florin Vaida and Suzette Blanchard},
   title = {Conditional Akaike information for mixed-effects models},
   journal = {Biometrika},
   year = {2005},
   volume = {92},
   pages = {351–370},
   abstract = {This paper focuses on the Akaike information criterion,
               AIC, for linear mixed-effects models in the analysis of
               clustered data. We make the distinction between questions
               regarding the population and questions regarding the
               particular clusters in the data. We show that the AIC in
               current use is not appropriate for the focus on clusters,
               and we propose instead the conditional Akaike information
               and its corresponding criterion, the conditional AIC,
               cAIC. The penalty term in cAIC is related to the effective
               degrees of freedom p for a linear mixed model proposed by
               Hodges & Sargent (2001); p reflects an intermediate level
               of complexity between a fixed-effects model with no
               cluster effect and a corresponding model with fixed
               cluster effects. The cAIC is defined for both maximum
               likelihood and residual maximum likelihood estimation. A
               pharmacokinetics data application is used to illuminate
               the distinction between the two inference settings, and to
               illustrate the use of the conditional AIC in model
               selection.},
   keywords = {Akaike information; AIC; effective degrees of freedom;
               linear mixed model}
}


HTH,
Henric



> 
> Regarding the fact the the extractAIC help page doesn't mention models
> with random effects, it can't list all the possible methods because
> any package can add methods to a generic function.
> 
>>
>>> On 12/18/07, Peter H Singleton <psingleton at fs.fed.us> wrote:
>>>> I am running a series of candidate mixed models using lmer (package lme4)
>>>> and I'd like to be able to compile a list of the AIC scores for those
>>>> models so that I can quickly summarize and rank the models by AIC. When I
>>>> do logistic regression, I can easily generate this kind of list by
>>>> creating
>>>> the model objects using glm, and doing:
>>>>
>>>>> md <- c("md1.lr", "md2.lr", "md3.lr")
>>>>> aic <- c(md1.lr$aic, md2.lr$aic, md3.lr$aic)
>>>>> aic2 <- cbind(md, aic)
>>>> but when I try to extract the AIC score from the model object produced by
>>>> lmer I get:
>>>>
>>>>> md1.lme$aic
>>>> NULL
>>>> Warning message:
>>>> In md1.lme$aic : $ operator not defined for this S4 class, returning NULL
>>>>
>>>> So... How do I query the AIC value out of a mixed model object created by
>>>> lmer?
>>>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>> -----
>> David Hewitt
>> Virginia Institute of Marine Science
>> http://www.vims.edu/fish/students/dhewitt/
>> --
>> View this message in context: http://www.nabble.com/How-can-I-extract-the-AIC-score-from-a-mixed-model-object-produced-using-lmer--tp14406832p14419438.html
>> Sent from the R help mailing list archive at Nabble.com.
>>
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

More information about the R-help mailing list