Hi all,
I'm using lme4::glmer to estimate Poisson mixed models in a very simple
context (single random effect). I'm interested in the model likelihood/AIC
across many simulated datasets.
To investigate whether the Laplace approximation was appropriate for my
data context, I explored using the argument nAGQ to improve the accuracy of
the likelihood estimation. When I changed nAGQ to a value > 1, I saw an
unexpectedly huge change in the likelihood; log-likelihoods tended to be
off by ~200. Other statistics packages (e.g. GLMMadaptive) yield estimates
that agree with lme4's Laplace approximation, as did a manual likelihood
estimate, and not with the nAGQ > 2 estimate.
The following code reproduces the problem I'm encountering.
*# r-sig-mixed-models GLMM question*
library(lme4)
set.seed(51)
*# Simulate some random effect-driven Poisson data*
random_effects <- rnorm(10, 0, 2)
group <- rep(1:10, 10)
simulated_data <- data.frame(y = rpois(n = 100, lambda = exp(3 +
random_effects[group])),
group = group)
*# Fit models with Laplace (nAGQ = 1) and nAGQ = 11*
fit_Laplace <- glmer(y ~ (1|group), data = simulated_data, family =
poisson())
fit_AGQ <- glmer(y ~ (1|group), data = simulated_data, family = poisson(),
nAGQ = 11)
logLik(fit_Laplace)
logLik(fit_AGQ)
logLik(fit_Laplace) - logLik(fit_AGQ) *# Huge difference!*
When I execute the above code, I see a difference in likelihood of
-218.8894. I've tested across many simulations and on 2 different machines
(Mac and Linux). My version of lme4 is up to date.
Has anyone run into this issue before? Am I using the glmer function wrong,
or is it possible there's something going on under the hood?
Thanks,
Ben
[[alternative HTML version deleted]]