canderson156 canderson156 - 2 months ago 16
R Question

R: lme4::gmler.nb function produces "Error in family$family : $ operator not defined for this S4 class" depending on the order I run models

library(lme4)

dummy <- as.data.frame(cbind(speed = rpois(100, 10), pop = rep(1:4, each = 25), season = rep(1:2, each = 50), id = seq(1, 100, by = 1)))
dummy2 <- as.data.frame(cbind(speed = c(rnbinom(50, 10, 0.6), rnbinom(50, 10, 0.1)), pop = rep(1:4, each = 25), season = rep(1:2, each = 50), id = seq(1, 100, by = 1)))


poisson <- glmer(speed~pop*season + (1|id),
data=dummy, family="poisson")

neg.bin <- glmer.nb(speed ~ pop*season + (1|id),
data=dummy2, control=glmerControl(optimizer="bobyqa"))


When I run a script creating a Poisson model before a negative binomial model using the lme4 package, I get the following error when running the neg.bin model:

Error in family$family : $ operator not defined for this S4 class


However, if I run the models in the opposite order, I don't the error message.

library(lme4)

dummy <- as.data.frame(cbind(speed = rpois(100, 10), pop = rep(1:4, each = 25), season = rep(1:2, each = 50), id = seq(1, 100, by = 1)))
dummy2 <- as.data.frame(cbind(speed = c(rnbinom(50, 10, 0.6), rnbinom(50, 10, 0.1)), pop = rep(1:4, each = 25), season = rep(1:2, each = 50), id = seq(1, 100, by = 1)))



neg.bin <- glmer.nb(speed ~ pop*season + (1|id),
data=dummy2, control=glmerControl(optimizer="bobyqa"))


poisson <- glmer(speed~pop*season + (1|id),
data=dummy, family="poisson")


The neg.bin model example does have convergence warnings, but the same pattern is happening with my actual models which are converging fine. How is running the Poisson model first affecting the neg.bin model?

Answer

Because you have masked R function poisson. The following would work fine (except that there is convergence warning for neg.bin):

library(lme4)
set.seed(0)
dummy <- as.data.frame(cbind(speed = rpois(100, 10), pop = rep(1:4, each = 25), season = rep(1:2, each = 50), id = seq(1, 100, by = 1)))
dummy2 <- as.data.frame(cbind(speed = c(rnbinom(50, 10, 0.6), rnbinom(50, 10, 0.1)), pop = rep(1:4, each = 25), season = rep(1:2, each = 50), id = seq(1, 100, by = 1)))

## use a different name for your model, say `poisson_fit`
poisson_fit <- glmer(speed~pop*season + (1|id),
         data=dummy, family="poisson")

negbin_fit <- glmer.nb(speed ~ pop*season + (1|id),
            data=dummy2, control=glmerControl(optimizer="bobyqa"))

Note, perhaps you wonder why masking poisson can still be an issue, even if you use family = "poisson" rather than family = poisson. Because if you feed family a string, glmer will later call match.fun to find the corresponding function poisson. This happens too for glm. Either way, we need to ensure that model fitting routines have no trouble getting the correct function.

It is always safe to use family = stats::poisson, by the way.