7
$\begingroup$

I have two LME models with the same interaction, one containing both main effects and one containing only one main effect, say :

$$ H\_CE = Season + Crownlevel + Season:Crownlevel , random = 1|CollectorID $$ and $$ H\_CE = Season + Season:Crownlevel , random = 1|CollectorID $$

There are 4 levels of each, and every combination of Season, Crownlevel and CollectorID The AIC, BIC and log likelihood of both models are completely equal. Given the formula for AIC being

$$ \mathit{AIC} = 2k - 2\ln(L)\ $$

one would expect this to be different, even if the likelihoods are exactly the same. In the end, they have a different number of parameters. Or so I thought...

Trying this toy example in R :

library(nlme)

Season <- rep(as.factor(rep(letters[1:4],each=4)),4)
Crownlevel <-rep(as.factor(rep(letters[11:14],4)),4)
CollectorID <-rep(letters[20:23],each=16)
X <-  model.matrix(~Season+Crownlevel+Season:Crownlevel)
B <- c(1,1,-2,2,0.3,0.4,0.4,2,3,1,-2,-3,-4,2,1,2)
H_CE <- X %*% B + rnorm(16*4)
KBM <- data.frame(Season,Crownlevel,H_CE,CollectorID)

model1 <- lme(H_CE~Season+Crownlevel+Season:Crownlevel,data=KBM,
       method="ML",random=~1|CollectorID)
model1e <- lme(H_CE~Season+Season:Crownlevel,data=KBM,
       method="ML",random=~1|CollectorID)

I get :

anova(model1,model1e)
        Model df      AIC      BIC    logLik
model1      1 18 174.1834 213.0433 -69.09168
model1e     2 18 174.1834 213.0433 -69.09168

What am I missing here? Why are the numbers completely equal? It has to do something with the model specification, but I can't really see what.

The model specification in itself is faulty, I know that. But I can't explain what makes it return a different set of parameters, but exactly the same residuals, likelihood and degrees of freedom :

> all.equal(residuals(model1),residuals(model1e))
[1] TRUE

As fabians rightfully pointed out, both models are perfectly equivalent. Yet, I fail to see why in the AIC calculation the same value for the number of parameters k is used.

The k in AIC uses the df, which explains everything.

$\endgroup$
6
$\begingroup$

The models are exactly equivalent. In both models you effectively specify one parameter for each combination of levels of Season and Crownlevel - the only difference is the parameterization:

In the first model, you fit main effects for Season and Crownlevel and an interaction effect to capture the combination-specific deviations from the main effects.

In the second model, you specify only the main effect of season, and the interaction effect then captures the deviations for each crownlevel within a season.

H_CE~Season:Crownlevel

would also yield an equivalent model, with one parameter for each combination of season and crownlevel (minus one that is non-identifiable because of the intercept, i.e. constitutes the reference category).

BTW: I don't think your model specification is faulty, which specification is better depends on the inference you want to do with your model.

$\endgroup$
  • $\begingroup$ thx for pointing out the obvious. I had my aha-moment when I calculated the df. I should have more coffee before I post a question. :-) The only thing I fail to see is how the AIC ends up equal, as k is the number of parameters and not the df. Regarding the model specification: collectorID is actually nested in Crownlevel, and as far as I'm concerned one shouldn't omit a main effect when the term is included in an interaction, unless we're talking about purely nested factors. But that's not the case for Season and Crownlevel. $\endgroup$ – Joris Meys Nov 23 '10 at 10:51
  • 1
    $\begingroup$ i'd think you have 16 fixed effects, one for each Season:Crownlevel combination (i.e. df=16) for both models plus two variance parameters (error + ID), so k has to be the same for both (they ARE exactly the same models, just with different parameterizations...). $\endgroup$ – fabians Nov 23 '10 at 11:04
  • $\begingroup$ but isn't k the number of parameters, and not the number of df? The different parametrization should give a different number of parameters, no? $\endgroup$ – Joris Meys Nov 23 '10 at 16:48
  • 2
    $\begingroup$ @Joris The df is the correct value to use, not the nominal number of parameters. The former has mathematical and statistical meaning whereas the latter does not. If you use more parameters than there are df, then one or more of those parameters will not be identifiable. $\endgroup$ – whuber Nov 23 '10 at 17:19
  • $\begingroup$ @whuber: see your point. Thx for the clarification $\endgroup$ – Joris Meys Nov 23 '10 at 21:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.