Hi, Prof. Bates (and all other interested parties), I haven't been using the lme4 package recently for binomial modeling but have recently returned for some analysis and found some puzzling problems. The current package version (0.999375-27) incorrectly estimates the dispersion parameter when using a quasibinomial. I think this is referred to indirectly in the following thread: https://stat.ethz.ch/pipermail/r-sig-mixed-models/2008q4/001501.html However, I have attached a working example that shows (I believe) the last known package where the dispersion is believable. For now, I'm using the older package to do the analysis but would be grateful for any corrections or workarounds. <code> # library(lme4, lib = "D:/R/archive") ## version 0.99875-9 library(lme4) ## version 0.999375-27 set.seed(2) ## simulate data A <- factor(rep(1:10, each = 10)) n <- length(A) rA <- rnorm(nlevels(A)) y <- rbinom(n, 1, binomial()$linkinv(0.5 + rA[A])) (f <- lmer(y ~ (1 | A), family = quasibinomial)) extractVar.lmer <- function (object, ...) { vc <- VarCorr(object) vdc <- unlist(lapply(lapply(vc, as.matrix), diag)) names(vdc) <- sub("\\.\\(.*\\)", "", names(vdc)) Var <- c(vdc, attr(vc, "sc")^2) names(Var)[length(Var)] <- "residual" Var } v <- extractVar.lmer(f) c(v[1]/v[2], v[2]) sessionInfo() </code> <output version="0.99875-9"> > (f <- lmer(y ~ (1 | A), family = quasibinomial)) Generalized linear mixed model fit using Laplace Formula: y ~ (1 | A) Family: quasibinomial(logit link) AIC BIC logLik deviance 130 135.2 -63 126 Random effects: Groups Name Variance Std.Dev. A (Intercept) 0.22075 0.46984 Residual 0.95095 0.97517 number of obs: 100, groups: A, 10 Fixed effects: Estimate Std. Error t value (Intercept) 0.7443 0.2572 2.893 > > v <- extractVar.lmer(f) > c(v[1]/v[2], v[1]) A residual 0.2321391 0.9509488 > sessionInfo() R version 2.7.2 (2008-08-25) i386-pc-mingw32 locale: LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] lme4_0.99875-9 Matrix_0.999375-15 lattice_0.17-15 loaded via a namespace (and not attached): [1] grid_2.7.2 </output> <output version="0.999375-27"> > (f <- lmer(y ~ (1 | A), family = quasibinomial)) Generalized linear mixed model fit by the Laplace approximation Formula: y ~ (1 | A) AIC BIC logLik deviance 132 139.8 -63 126 Random effects: Groups Name Variance Std.Dev. A (Intercept) 0.046833 0.21641 Residual 0.205057 0.45283 Number of obs: 100, groups: A, 10 Fixed effects: Estimate Std. Error t value (Intercept) 0.7463 0.1191 6.264 > > v <- extractVar.lmer(f) > c(v[1]/v[2], v[1]) A residual 0.2283881 0.2050567 > > sessionInfo() R version 2.7.2 (2008-08-25) i386-pc-mingw32 locale: LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] lme4_0.999375-27 Matrix_0.999375-15 lattice_0.17-15 loaded via a namespace (and not attached): [1] grid_2.7.2 </output> Thanks, --sundar
Binomial GLMM with different versions of lme4
1 message · Sundar Dorai-Raj