Fitting a Quasi-Poisson Model with lme4 or glmmTMB
=====================================================
In this post, we’ll explore how to fit a mixed-effects quasi-poisson model using the lme4 package in R. We’ll also cover how to do it with the glmmTMB package, which is known for its flexibility and accuracy.
What is a Quasi-Poisson Model?
A quasi-poisson model is an extension of the Poisson distribution that accounts for overdispersion, or excessive variation in the data. This type of model is commonly used when the response variable has a high variance-to-mean ratio, which can occur in situations like counting data (e.g., number of ticks on a chick) or time-to-event data.
Why Use Quasi-Poisson Models?
Quasi-poisson models offer several advantages over traditional Poisson models:
- They can handle overdispersion more effectively.
- They provide accurate estimates of variance components, which is essential for understanding the structure of the data.
- They can be used with any distribution that can be transformed into a Poisson-like distribution.
Fitting Quasi-Poisson Models with lme4
Step 1: Load Libraries and Data
First, we need to load the required libraries and data. In this example, we’ll use the grouseticks dataset from Elston et al (2001).
library(lme4)
data(grouseticks)
Step 2: Prepare Data
Next, we prepare our data by transforming the HEIGHT variable using the scale function.
g <- transform(grouseticks, sHEIGHT = drop(scale(HEIGHT)))
form <- TICKS ~ YEAR + sHEIGHT + (1|BROOD) + (1|LOCATION)
Step 3: Fit Model
Now we fit our quasi-poisson model using the glmer function.
full_mod1 <- glmer(form, family = "poisson", data = g)
Step 4: Compute Deviance and Residuals
We compute the deviance and residuals of the fitted model.
deviance(full_mod1) / df.residual(full_mod1)
residuals(full_mod1, type = "pearson")
Fitting Quasi-Poisson Models with glmmTMB
Step 1: Load Libraries and Data
First, we load the required libraries and data. In this example, we’ll use the grouseticks dataset from Elston et al (2001).
library(glmmTMB)
data(grouseticks)
Step 2: Prepare Data
Next, we prepare our data by transforming the HEIGHT variable using the scale function.
g <- transform(grouseticks, sHEIGHT = drop(scale(HEIGHT)))
form <- TICKS ~ YEAR + sHEIGHT + (1|BROOD) + (1|LOCATION)
Step 3: Fit Model
Now we fit our quasi-poisson model using the glmerMB function.
full_mod1 <- glmerMB(form, family = "poisson")
Quasi-likelihood Adjustment Function
To apply a quasi-likelihood adjustment to our standard errors and p-values, we can use the following function:
quasi_table <- function(model, ctab = coef(summary(model))) {
phi <- sum(residuals(model, type = "pearson")^2) / df.residual(model)
qctab <- within(as.data.frame(ctab), {
Std. Error <- Std. Error * sqrt(phi)
z value <- Estimate / Std. Error
Pr(>|z|) <- 2*pnorm(abs(z value), lower.tail = FALSE)
})
return(qctab)
}
Applying Quasi-likelihood Adjustment
We apply this function to our fitted model.
printCoefmat(quasi_table(full_mod1), digits = 2)
Comparing Results with lme4 and glmmTMB
As shown in the example, applying the quasi-likelihood adjustment function yields identical results for both lme4 and glmmTMB. The estimates are unchanged, but the standard errors and p-values have been adjusted accordingly.
The broom.mixed package can also be used to tidy the output of glmer models. This provides a more concise and readable format for the results.
library(tidyverse)
library(broom.mixed)
phi <- sum(residuals(full_mod1, type = "pearson")^2) / df.residual(full_mod1)
(full_mod1 %>%
tidy(effect = "fixed") %>%
transmute(term = term, estimate = estimate,
std.error = std.error * sqrt(phi),
statistic = estimate/std.error,
p.value = 2*pnorm(abs(statistic), lower.tail = FALSE)))
Conclusion
In this post, we explored how to fit a mixed-effects quasi-poisson model using lme4 and glmmTMB. We also covered the concept of quasi-likelihood adjustment and provided an example function that can be used to adjust standard errors and p-values. The results obtained from both packages are identical when applying this adjustment.
We hope that this post has been informative and helpful in understanding how to fit quasi-poisson models using R packages like lme4 and glmmTMB.
Last modified on 2024-04-13