[R-sig-ME] anova in lmer modifies REML to ML

[Ken Knoblauch]

Hi,

I noticed that the AIC, BIC, etc. returned in printing
an object obtained by fitting with lmer with the default
REML method differ from those returned by  comparing two "mer"
objects with anova.  It looks as if the anova method
automatically uses the ML rather than the REML estimates.
In addition, in Ch. 2 of the online draft of
"lme4: Mixed-effects Modeling with R" on R-Forge,
all of the models that are compared using anova are fit
with the argument "REML = 0", even though only the random
effects vary across models.  I understand from reading
Pinheiro & Bates that fits using the REML criterion with
different fixed effects  structure cannot be compared because
of the  lack of invariance across one-to-one transformations,
but that REML fits can be compared if the fixed effects
structure does not change (p. 83 of P&B).  If I haven't
gotten anything wrong above, my question is whether
the behavior of anova, transforming to ML, is on purpose,
for example, to safeguard against those of us who might
try to compare REML fits with different fixed effects,
or perhaps the comparisons of the ML fits are preferred
even when only the random effects have been modified?
I didn't see (haven't quite finished it yet) any mention
of this behavior in the book draft and didn't find this
described on the help pages (but I may have missed something,
too), which is why I'm wondering if this is the intended
behavior.  Thank you.

A short example demonstrating this follows, using an online
data set that I have been working.

library(lme4)

###get the data and put it into long format
site <- "http://vision.arc.nasa.gov/modelfest/data/modelfestbaselinedata.csv"
  ModelFest <- read.csv(site, header = FALSE)
  Obs <- ModelFest$V1
  ModelFest.df <- data.frame(LContSens = c(t(ModelFest[, -1])),
	Obs = rep(Obs, each = ncol(ModelFest) - 1),
	Stim = rep(paste("Stim", 1:43, sep = ""), each = 4)
	)
ModelFest.df$Stim <- with(ModelFest.df, factor(Stim,
		levels = unique(Stim)))

mf1.lmer <- lmer(LContSens ~ Stim + (1 | Obs),
	ModelFest.df)
mf2.lmer <- update(mf1.lmer, . ~ . + (1 | Obs:Stim))

summary(mf1.lmer)@AICtab
summary(mf2.lmer)@AICtab
anova(mf1.lmer, mf2.lmer)
summary(update(mf1.lmer, REML = 0))@AICtab
summary(update(mf2.lmer, REML = 0))@AICtab


Ken


sessionInfo()
R version 2.10.1 Patched (2010-03-10 r51250)
i386-apple-darwin9.8.0

locale:
[1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods
[7] base

other attached packages:
[1] lme4_0.999375-32   Matrix_0.999375-38 lattice_0.18-3

loaded via a namespace (and not attached):
[1] grid_2.10.1  tools_2.10.1




-- 
Ken Knoblauch
Inserm U846
Stem-cell and Brain Research Institute
Department of Integrative Neurosciences
18 avenue du Doyen Lépine
69500 Bron
France
tel: +33 (0)4 72 91 34 77
fax: +33 (0)4 72 91 34 61
portable: +33 (0)6 84 10 64 10
http://www.sbri.fr/members/kenneth-knoblauch.html