There are both statistical and computational issues involved here (I would not recommend doing both p-value and AIC comparisons - stick to one or the other!).
The first thing to do is to refit the original models with
REML=FALSE; I'm pretty certain that the AIC comparison does not prevent you doing the silly thing of comparing the 'log-likelihoods' (more properly 'REML criteria', as what you get when
REML=TRUE is not a log-likelihood), as
anova() does ...
We really should have considered making
NA for REML-fitted models to help people avoid this trap ...
Here's an example.
Fit full and reduced models with ML and REML:
fm1_REML <- lmer(Reaction~Days+(Days|Subject),sleepstudy)
fm0_REML <- update(fm1_REML,. ~ . - Days)
fm1_ML <- update(fm1_REML,REML=FALSE)
fm0_ML <- update(fm0_REML,REML=FALSE)
anova(fm1_ML,fm0_ML) ## p-value: 1.226e-06
anova(fm1_REML,fm0_REML) ## 'refitting models with ML': same p-value
AIC comparisons (I like
bbmle::AICtab for readability)
AICtab(fm1_ML,fm0_ML) ## delta-AIC=21.5
AICtab(fm1_REML,fm0_REML) ## delta-AIC=24.2
Here's how you would override the
anova() safeguards to do a bogus likelihood ratio test on the REML-fitted models ...
df=1,lower.tail=FALSE) ## 3.05e-07
This p-value is an order of magnitude smaller than the correct LRT, which matches up well with the larger delta-AIC in the bogus AIC comparison based on the REML-fitted models ...