Fitting Quasi-Poisson Models with lme4 or glmmTMB: A Comparative Analysis

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(&gt;|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