I am trying to fit a one-inflated Beta regression mixed model using the gamlss package using the following code:
# Minimal reproducible example:
efficiency <- rbeta(n = 85, shape1 = 2, shape2 = 2, ncp = 2) # My actual variable contains several "ones", but it does not matter here
distance <- rpois(n = 85, lambda = 4)
manager_id <- c(rep(1, 10), rep(2, 2), 3:75)
manager_id <- as.factor(manager_id)
eff <- data.frame(efficiency, distance, manager_id)
mod1 <- gamlss::gamlss(efficiency~distance + gamlss::re(random = ~1|manager_id, method = "ML"),
data = eff, family = gamlss.dist::BEOI())
Yet I get the following error:
Error in as.data.frame.default(x[[i]], optional = TRUE) :
cannot coerce class ‘"function"’ to a data.frame
I believe the error comes from the nlme::lme() function that is called by gamlss::re() to fit the random effect because it works fine when I do not include this part of the code. I could perhaps use gamlss::random() instead of gamlss::re() to fit the random effect but the help page of these functions states that re() should be used when the normality assumption of the response variable cannot be satisfied.
Could someone please help me figure out how to solve this problem?