Since "quasi" families cannot be used in glmer, there is an efficient approach to adjust the standard errors of the parameters and the associated statistics post-fitting (https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#fitting-models-with-overdispersion). But as a result, I only have the updated summary of the model with the fitted parameters. How do I proceed if I want to perform post-hoc tests on the model with the quasi-likelihood adjusted parameters, such as pairwise comparisons of user-defined contrasts with emmeans?
I tried to override the glmer models parameter with the adjusted parameters and then use the manipulated model with emmeans, but I don't know how to manipulate the glmer models parameter or if this even works.
Thank you!
You can compute the dispersion parameter yourself (as the sum of squared Pearson residuals divided by the residual df) and use the vcov. argument to emmeans to scale the covariance matrix accordingly ...
library(lme4)
library(emmeans)
gm1 <- glmer(cbind(incidence, size-incidence) ~ period + (1|herd),
family = binomial,
data = cbpp)
phi <- sum(residuals(gm1, type = "pearson")^2)/df.residual(gm1)
em1 <- emmeans(gm1, ~ period)
em2 <- emmeans(gm1, ~ period, vcov. = vcov(gm1)*phi)
pairs(em1)
contrast estimate SE df z.ratio p.value
period1 - period2 0.992 0.303 Inf 3.272 0.0059
period1 - period3 1.128 0.323 Inf 3.495 0.0027
period1 - period4 1.580 0.422 Inf 3.743 0.0010
period2 - period3 0.136 0.376 Inf 0.363 0.9837
period2 - period4 0.588 0.464 Inf 1.267 0.5843
period3 - period4 0.452 0.478 Inf 0.945 0.7807
pairs(em2)
contrast estimate SE df z.ratio p.value
period1 - period2 0.992 0.337 Inf 2.942 0.0172
period1 - period3 1.128 0.359 Inf 3.143 0.0091
period1 - period4 1.580 0.469 Inf 3.366 0.0042
period2 - period3 0.136 0.418 Inf 0.326 0.9880
period2 - period4 0.588 0.516 Inf 1.139 0.6654
period3 - period4 0.452 0.531 Inf 0.850 0.8307
Note that the parameters are not adjusted by quasi-likelihood adjustments, only their (co)variances/standard deviations.