MASS::glmmPQL diagnostic

Question

I am fitting models with MASS::glmmPQL of the form

MASS::glmmPQL(Y ~ X1+X2+X3+...,
               random=c(~1|ID), data=df , family = quasipoisson(link = log) )      

where X1...Xn are continuous predictors.

1. what model diagnostics should I look at?
2. how do I implement these diagnostics in R? I understand there is a package DHARMa that seems the perfect panacea but it doesn't work with QuasiLikelihood.
3. If you don't want to give me suggestions... Why is nobody replying to similar questions?

thanks.

Answer 1

There are problems with the MASS library's glmmPQL. It does not return a log-likelihood, so model selection will be difficult.

As for what to look at: check for warnings(), look at fixed effect estimates and significance, look at random effects variance estimates (and confidence intervals for those if you can), and look at the overdispersion parameter estimate. I would say to also look at the log-likelihood so you can compare nested models, but... I've already said that could be an issue.

There are, however, a few other options for you.

Older versions of lme4 had a method parameter for lmer and glmer which allowed you to specify PQL. While the error reporting was sometimes scant (you needed to check warnings about separability to see if your variance estimates were nonsense). However, apart from that issue the code was better in that it returned a log-likelihood which you could use.

You could also try glmer's Laplace or quadrature (say, for 2-3 points) methods. Those are typically slower than PQL, but they are similar in that they approximate the integrand and may be more useful than an old version of glmer or not getting a log-likelihood back from glmmPQL.