Why does AIC model rank order change in lme models with standardization of predictor variables?

Question

I can't figure this out. The AIC/AICc rank of my mixed effect models are different whether or not I standardize my predictor values using rescale.

Note, I'm not concerned that AICc is changing, as that is expected, but I am concerned that the model rank order is changing.

The standardization method being used is subtracting the mean and dividing by 2 standard deviations.

As expected, p-values and t-values do not vary between standardized and unstandardized models, just AICc, and most importantly, AICc model rank.

All variables a-k are positive continuous numeric (0-∞) and a1 and b1 are proportions with range 0-1.

# unstandardized
data_n <- data

# standardize variables
data$k <- arm::rescale(data$k)
data$a <- arm::rescale(data$a)
data$b <- arm::rescale(data$b)
data$c <- arm::rescale(data$c)
data$d <- arm::rescale(data$d)
data$e <- arm::rescale(data$e)
data$f <- arm::rescale(data$f)
data$g <- arm::rescale(data$g)
data$a1 <- arm::rescale(data$a1)
data$b1 <- arm::rescale(data$b1)

Testing AICc Rank with and without standardization

M_17 <- lme(fixed = log(y) ~ k + d + a1, data=data, random = ~ 1|Plot_No)
AICc(M_17)

## [1] 174.8177

M_17_n <- lme(fixed = log(y) ~ k + d + a1, data=data_n, random = ~ 1|Plot_No)
AICc(M17_n)

## [1] 184.1404

M_31 <- lme(fixed = log(y) ~ k + f + d, data=data, random = ~ 1|Plot_No)
AICc(M_31)

## [1] 178.2103

M_31_n <- lme(fixed = log(y) ~ k + f + d, data=data_n, random = ~ 1|Plot_No)
AICc(M31_n)

## [1] 191.7922

M_71 <- lme(fixed = log(y) ~ k + b, data=data, random = ~ 1|Plot_No)
AICc(M_71)

## [1] 179.2528

M_71_n <- lme(fixed = log(y) ~ k + b, data=data_n, random = ~ 1|Plot_No)
AICc(M71_n)

## [1] 190.5091

As shown, the ranking of M_31 and M_71 changes whether I use standardization or not.

I've been told that this isn't supposed to happen, but it is. Please help!

Answer 1

I had a very similar problem that a scatterplot showed a negative trend between my variable of interest Y and elevation but my lme model Y ~ elev... had a higher AIC than the null model: Y ~ 1.... I recalled an earlier error message with lme4 (don't think I've ever received it with using nlme package) suggesting that I should rescale the predictor variable if the order of magnitude varies considerably. When I divided elevation by 1000 and reran the model the AIC/BIC/logLik changed considerably and AIC was well below that of the null model, as I had anticipated from the scatterplot. When I repeated this example with using lmer in lme4 I got identical AIC/BIC/logLik values whether I used elev or elev/1000, as one would expect.

Answer 2

Beats me. Your codes are OK and (like you said) it shouldn't happen. I have a suggestion. Why don't you try "lmer" (instead of lme), i.e., with your data and using everything you wrote above, try the following codes, and see if you still have AICc values come out differently. Who knows you may have found some serious problems. Before panicking, let's check the calculations first.

library(lme4) m_17 <- lmer(log(y) ~ k + d + a1 + (1|Plot_No), data=data) m_17_n <- lmer(log(y) ~ k + d + a1 + (1|Plot_No), data=data_n) m_31 <- lmer(log(y) ~ k + f + d + (1|Plot_No), data=data) m_31_n <- lmer(log(y) ~ k + f + d + (1|Plot_No), data=data_n) m_71 <- lmer(log(y) ~ k + b + (1|Plot_No), data=data) m_71_n <- lmer(log(y) ~ k + b + (1|Plot_No), data=data_n) AICc(m17); AICc(m17_n) AICc(m31); AICc(m31_n) AICc(m71); AICc(m71_n)