I want to be able to pass a named list of models (merMod objects) to anova() and preserve the model names in the output. This is particularly useful in the context of using mclapply() to run a batch of slow models like glmers more efficiently in parallel. The best I've come up with is to use do.call on a de-named version of the model list, but that's not ideal because I might have models named (say) "mod12", "mod17", and "mod16" and these model names get translated to "MODEL1", "MODEL2", and "MODEL3" in the output. (That might seem trivial when looking at a single batch, but over the course of a long modeling session with dozens of models it's a surefire recipe for confusion.)
Note that this isn't the same issue as Create and Call Linear Models from List because I'm not trying to compare pairs of models across lists. It's also more complex than Using lapply on a list of models because I'm using anova() in a non-unary way.
Here's a minimal reprex:
library(lme4)
formList <- c(mod12 = angle ~ recipe + temp + (1|recipe:replicate),
mod17 = angle ~ recipe + temperature + (1|recipe:replicate),
mod16 = angle ~ recipe * temperature + (1|recipe:replicate))
modList <- lapply(formList, FUN=lmer, data=cake)
# Fails because modList is named so it's interpreted as arg-name:arg pairs
do.call(anova, modList)
# Suboptimal because model names aren't preserved
do.call(anova, unname(modList))
# Fails because object isn't merMod (or some other class covered by methods("anova"))
do.call(anova, list(object=modList[1], ...=modList[-1], model.names=names(modList)))
The second do.call returns this:
Data: ..1
Models:
MODEL1: angle ~ recipe + temp + (1 | recipe:replicate)
MODEL2: angle ~ recipe + temperature + (1 | recipe:replicate)
MODEL3: angle ~ recipe * temperature + (1 | recipe:replicate)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
MODEL1 6 1708.2 1729.8 -848.08 1696.2
MODEL2 10 1709.6 1745.6 -844.79 1689.6 6.5755 4 0.1601
MODEL3 20 1719.0 1791.0 -839.53 1679.0 10.5304 10 0.3953
Ideally, the output would look like this:
Data: ..1
Models:
mod12: angle ~ recipe + temp + (1 | recipe:replicate)
mod17: angle ~ recipe + temperature + (1 | recipe:replicate)
mod16: angle ~ recipe * temperature + (1 | recipe:replicate)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
mod12 6 1708.2 1729.8 -848.08 1696.2
mod17 10 1709.6 1745.6 -844.79 1689.6 6.5755 4 0.1601
mod16 20 1719.0 1791.0 -839.53 1679.0 10.5304 10 0.3953
How do I do this? I'm more than happy with an ugly wrapper around anova() if it means I get a more intelligible output.
You can pass a list of symbols like this:
with(modList,
do.call(anova,
lapply(names(modList), as.name)))
#refitting model(s) with ML (instead of REML)
#Data: ..1
#Models:
#mod12: angle ~ recipe + temp + (1 | recipe:replicate)
#mod17: angle ~ recipe + temperature + (1 | recipe:replicate)
#mod16: angle ~ recipe * temperature + (1 | recipe:replicate)
# Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
#mod12 6 1708.2 1729.8 -848.08 1696.2
#mod17 10 1709.6 1745.6 -844.79 1689.6 6.5755 4 0.1601
#mod16 20 1719.0 1791.0 -839.53 1679.0 10.5304 10 0.3953
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With