Charlie_J
Charlie_J

Reputation: 109

Error: requires numeric/complex matrix/vector arguments for %*%; cross validating glmmTMB model

I am adapting some k-fold cross validation code written for glmer/merMod models to a glmmTMB model framework. All seems well until I try and use the output from the model(s) fit with training data to predict and exponentiate values into a matrix (to then break into quantiles/number of bins to assess predictive performance). I can get get this line to work using glmer models, but it seems when I run the same model using glmmTMB I get Error in model.matrix: requires numeric/complex matrix/vector arguments There are many other posts out there discussing this error code and I have tried converting the data frame into matrix form and changing the class of the covariates with no luck. Separately running the parts before and after the %*% works but when combined I get the error. For context, this code is intended to be run with use/availability data so the example variables may not make sense, but the problem gets shown well enough. Any suggestions as to what is going on?

library(lme4)
library(glmmTMB)

# Example with mtcars dataset
data(mtcars)

# Model both with glmmTMB and lme4
m1 <- glmmTMB(am ~ mpg + wt + (1|carb), family = poisson, data=mtcars)
m2 <- glmer(am ~ mpg + wt + (1|carb), family = poisson, data=mtcars)

#--- K-fold code (hashed out sections are original glmer version of code where different)---

# define variables
k <- 5
mod <- m1 #m2
dt <- model.frame(mod) #data used

reg.list <- list() # initialize object to store all models used for cross validation
  
# finds the name of the response variable in the model dataframe
  resp <- as.character(attr(terms(mod), "variables"))[attr(terms(mod), "response") + 1]
    
# define column called sets and populates it with character "train"
  dt$sets <- "train"
    
# randomly selects a proportion of the "used"/am records (i.e. am = 1) for testing data 
    dt$sets[sample(which(dt[, resp] == 1), sum(dt[, resp] == 1)/k)] <- "test"
    
# updates the original model using only the subset of "trained" data
    reg <- glmmTMB(formula(mod), data = subset(dt, sets == "train"), family=poisson,
                   control = glmmTMBControl(optimizer = optim, optArgs=list(method="BFGS")))

    #reg <- glmer(formula(mod), data = subset(dt, sets == "train"), family=poisson,
    #            control = glmerControl(optimizer = "bobyqa", optCtrl=list(maxfun=2e5)))
    
 reg.list[[i]] <- reg # store models
    
# uses new model created with training data (i.e. reg) to predict and exponentiate values 
      predall <- exp(as.numeric(model.matrix(terms(reg), dt) %*% glmmTMB::fixef(reg)))
     #predall <- exp(as.numeric(model.matrix(terms(reg), dt) %*% lme4::fixef(reg)))

Upvotes: 1

Views: 244

Answers (1)

Ben Bolker
Ben Bolker

Reputation: 226761

Without looking at the code too carefully: glmmTMB::fixef(reg) returns a list (with elements cond (conditional model parameters), zi (zero-inflation parameters), disp (dispersion parameters) rather than a vector.

If you replace this bit with glmmTMB::fixef(reg)[["cond"]] it will probably work.

Upvotes: 2

Related Questions