Suppose I have four linear mixed effects models as shown:

mod1 <- lmer(X ~ (1|D), data = df)
mod2 <- lmer(X ~ A + (1|D), data = df)
mod3 <- lmer(X ~ A + B + (1|D), data = df)
mod4 <- lmer(X ~ A + B + C + (1|D), data = df)

All models have the same random effects and each model just adds another new variable to the previous model, so I believe they fit the criteria for nested models? I want to compare the variance explained in the fixed effects of the response (X) between all the four models. My initial approach is to just compare all the models via an ANOVA as such:

anova(mod1, mod2, mod3, mod4)

However, from my understanding this would just compare mod1 vs mod2, then mod2 vs mod3, and mod3 vs mod4, however if I want to compare all combinations of models, i.e.:

• mod1 vs mod2
• mod1 vs mod3
• mod1 vs mod4
• mod2 vs mod3
• mod2 vs mod4
• mod3 vs mod4

Do I need to create separate ANOVAs for each of these comparisons, or is there a more efficient way to do this in R? Is there anything that needs to be taken into account when doing these comparisons when using the lmer package vs lm package (linear models)?

Thank you in advance.