# AIC/BIC vs the rule of “must include lower order interaction”

I am running a series of mixed effect models, which include both linear and quadratic term of a variable T (continuous) and the main IV I (categorical), and facing a dilemma.

Model 2 include interaction terms with both the linear and quadratic T. Model 1 include only the quadratic interaction.

With Model 1, I obtained the lowest AIC and BIC. In comparison, for Model 2 the AIC is 5.8 higher and BIC is 26 higher. If I followed the rule of thumb for AIC/BIC-based selection, Model 1 is the one I should select (and the one I preferred).

On the other hand, there's also a lot of discussion regarding the necessity of including the lower order term interaction when quadratic interaction is in the model.

However, it will be very difficult for me to justify the selection of a model that's very unlikely (probability: exp(-5.8/2)=5%) to be better than the one with AICmin. Thus the dilemma.

• 1) Did you mix things up here? You say you should select model 1 which has linear and quadratic term. Note by the way that some literature has AIC and BIC as "larger is better" and some as "lower is better"; you should always say which version you are referring to. 2) Is there a subject matter reason in your specific situation to not include a linear term? Sometimes it can be justified. 3) " model that's very unlikely (probability: exp(-5.8/2)=5%) to be better" - this looks like a misinterpretation of the AIC. There is no such probability for a model to be better, unless you go Bayesian. – Lewian Oct 2 '20 at 23:47
• @Lewian Thanks for the comment. I follow the recommendation that the lower the AIC and BIC, the better. Model 2 has a AIC that's 5.8 higher and BIC that's 26 higher, compared with Model 1. And for the probability, I was using this guidance "the quantity exp((AICmin − AICi)/2) can be interpreted as being proportional to the probability that the ith model minimizes the (estimated) information loss." Finally, there is no theoretical justification for including or excluding the linear interaction term. – user6606453 Oct 3 '20 at 1:53
• @Lewian running out of space :) We have always included the linear term T and quadratic term T^2, but we are unsure if we can include the interaction term I:T^2 while leaving out the interaction term I:T. – user6606453 Oct 3 '20 at 1:58
• So your question says "With Model 2, I obtained the lowest AIC and BIC", but then that its AIC is "higher". This looks contradictory. You also use "lower is better" versions of AIC, so effectively it now reads like Model 1 is better according to AIC. Model 1 has linear and quadratic interactions as you state, so what's wrong with it? To me the things that you state look inconsistent and I strongly recommend to revise the writing of your question. – Lewian Oct 3 '20 at 9:10
• " I was using this guidance "the quantity exp((AICmin − AICi)/2) can be interpreted as being proportional to the probability that the ith model minimizes the (estimated) information loss." You don't state where this comes from and within what model this "probability" is defined. The estimated information loss is observable, so it does not make sense to talk about a probability to minimise it. The guidance looks either confused or taken out of context. – Lewian Oct 3 '20 at 9:12

It's certainly possible to have a situation in which only the interaction is significant. A trivial example is:

> set.seed(12)
> x1 <- rnorm(100,0,1)
> x2 <- rnorm(100,0,1)
> cor(x1,x2)
[1] 0.01592198
> y <- x1*x2 + rnorm(100,0,.3)
> modFull <- lm(y~x1+x2+x1:x2)
> modIntOnly <- lm(y~x1:x2)
> anova(modIntOnly,modFull)
Analysis of Variance Table

Model 1: y ~ x1:x2
Model 2: y ~ x1 + x2 + x1:x2
Res.Df    RSS Df Sum of Sq     F Pr(>F)
1     98 8.8677
2     96 8.7568  2   0.11093 0.608 0.5465
> AIC(modFull)
[1] 50.25384
> AIC(modIntOnly)
[1] 47.51263


The individual x1 and x2 terms clearly add nothing to this model, and have insignificant p-values in the summary of modFull (not shown). The AIC is lower (better) when they are omitted. But now let's say that instead of having available the actual values of x1 and x2 underlying this true model, we only have w1 and w2 that are shifted by exactly 10 units each:

> w1<-x1+10
> w2<-x2+10
> modelShiftedFull <- lm(y~w1+w2+w1:w2)
> modelShiftedIntOnly <- lm(y~w1:w2)
> anova(modelShiftedFull,modelShiftedIntOnly)
Analysis of Variance Table

Model 1: y ~ w1 + w2 + w1:w2
Model 2: y ~ w1:w2
Res.Df    RSS Df Sum of Sq      F    Pr(>F)
1     96  8.757
2     98 94.034 -2   -85.277 467.44 < 2.2e-16
> AIC(modelShiftedFull)
[1] 50.25384
> AIC(modelShiftedIntOnly)
[1] 283.6363


Now the individual w1 and w2 terms are critical to the model, with highly significant p-values in modelShiftedFull (not shown). It's the same underlying "reality" in terms of the outcome y, captured identically with respect to AIC by both "full" models. The only difference is that we now have available to us values of the independent variables that have been shifted in location from their "true" values. How sensitive do you want your model to be with respect to this type of problem?

So "must include lower order interaction terms" is more complicated than a "rule." The issue is discussed extensively on this site, for example here and here. I think of it as the generally best way to avoid unexpected problems. Furthermore, unless you are in danger of overfitting, there is little to be lost by including nominally "insignificant" predictors in a model, given the risks of omitted-variable bias. So I would turn this question around: why omit the lower-order terms at all?

One last thought: consider a flexible restricted cubic spline instead of a quadratic fit of your continuous predictor. Then you can readily evaluate the degree of flexibility (number of knots) needed to model your data without imposing a particular parametric form. I suspect that such a model would be superior, by the AIC criterion, to your quadratic fit.

• Thanks for the example. It is certainly illustrative of something I wasn't fully aware. I know nothing about flexible restricted cubic spline in modeling. Is there a good place I can start, especially its comparison with quadratic term? – user6606453 Oct 3 '20 at 21:06
• @user6606453 restricted cubic spline regression is provided by the rcs() function in Frank Harrell's rms package in R, but should be available in most serious statistical software. Harrell briefly introduces them here, describing them in more detail in section 2.4 of his course notes. They are different from penalized splines, which can lead to some confusion. There are some comparisons with polynomials here. – EdM Oct 3 '20 at 22:25