glmmTMB: Error when pooling models

[reiniervlinschoten]

When I try to compare fitted glmmTMB models on multiply imputed data (through MICE), I get an error when trying to use the pooled Wald test. There seems to be a problem when parameter estimates and the variance/covariance matrix are extracted:
amices/mice#435

I get the following error when trying to use anova() or D1() from the mice package on the mira objects.

all(p == dim(Uhat[[1]])) is not TRUE

I guess a workaround similarly to the one implemented for polr is necessary?

See also the reprex below.

##libraries
library(glmmTMB)
library(mice)

dff<-list()
dff[[1]] <- data.table::data.table(
  sbq_1 = c(1L,2L,1L,2L,2L,
            1L,2L,2L,4L,2L,1L,4L,2L,3L,1L,4L,2L,1L,
            3L,1L,4L,1L,3L,4L,2L,1L,2L,2L,2L,2L,
            1L,4L,3L,4L,1L,2L,3L,4L,1L,1L,3L,3L,1L,
            1L,2L,1L,1L,1L,1L),
  age = c(53L,32L,47L,58L,
          34L,29L,57L,25L,38L,29L,51L,50L,37L,33L,
          41L,36L,49L,38L,32L,51L,37L,36L,37L,39L,
          55L,55L,45L,25L,24L,27L,26L,57L,47L,30L,
          37L,39L,38L,34L,38L,56L,39L,43L,49L,39L,48L,
          44L,63L,61L,46L),
  gender = c(1L,1L,2L,1L,1L,
             1L,1L,2L,2L,1L,1L,1L,1L,1L,1L,2L,1L,1L,
             2L,1L,1L,1L,1L,1L,2L,1L,1L,2L,1L,1L,
             1L,1L,2L,1L,2L,1L,1L,2L,1L,1L,1L,2L,1L,
             1L,2L,1L,1L,2L,1L),
  ethnicity = c(2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L,1L,1L,2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L)
)

dff[[2]] <- data.table::data.table(
  sbq_1 = c(1L,2L,1L,2L,2L,
            1L,2L,2L,4L,2L,1L,4L,2L,3L,1L,4L,2L,1L,
            3L,1L,4L,1L,3L,4L,2L,1L,2L,2L,2L,2L,
            1L,4L,3L,4L,1L,2L,3L,4L,1L,1L,3L,3L,1L,
            1L,2L,1L,1L,1L,1L),
  age = c(53L,32L,47L,58L,
          34L,29L,57L,25L,38L,29L,51L,50L,37L,33L,
          41L,36L,49L,38L,32L,51L,37L,36L,37L,39L,
          55L,55L,45L,25L,24L,27L,26L,57L,47L,30L,
          37L,39L,38L,34L,38L,56L,39L,43L,49L,39L,48L,
          44L,63L,61L,46L),
  gender = c(1L,1L,2L,1L,1L,
             1L,1L,2L,2L,1L,1L,1L,1L,1L,1L,2L,1L,1L,
             2L,1L,1L,1L,1L,1L,2L,1L,1L,2L,1L,1L,
             1L,1L,2L,1L,2L,1L,1L,2L,1L,1L,1L,2L,1L,
             1L,2L,1L,1L,2L,1L),
  ethnicity = c(2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L,1L,2L,2L,2L,2L,2L,2L,
                2L,2L,2L,2L,1L,2L,2L,2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L,2L,1L,2L,2L,2L,2L,2L,
                2L,2L,2L,2L,2L,2L)
)

## backward elimination on supermodel

with_term <- list()
for (i in 1:2){
  mod <- list(glmmTMB(as.factor(sbq_1) ~ age + gender + (1 | ethnicity), data=dff[[i]]))
  with_term <- c(with_term, mod)
}

without_term <- list()
for (i in 1:2){
  mod <- list(glmmTMB(as.factor(sbq_1) ~ age + (1 | ethnicity), data=dff[[i]]))
  without_term <- c(without_term, mod)
}

fit.with <- as.mira(with_term)
fit.without <- as.mira(without_term)

D1(fit.with, fit.without)
# This gives an error
# Error in .extractParameters(model, diagonal = FALSE) : all(p == dim(Uhat[[1]])) is not TRUE

[simongrund1]

Thanks for the report. The issue was that mitml didn't support the different submodels of glmmTMB. I've added some preliminary support for glmmTMB that should allow pooling parameter estimates in the mean structure model (with testEstimates()) as well as comparisons between models that differ in their mean structure (with testModels()).

You can install the development version like this:

devtools::install_github("simongrund1/mitml")

I think this should also fix the problem with D1 in mice, but I haven't tested this. Let me know if the issue persists.

[reiniervlinschoten]

Fixed! Thank you very much :)

[minre2024]

Hi! I got the same error when trying to run following code. I, also, run this code: "devtools::install_github("simongrund1/mitml")" but it didn't work. Could you help to solve this problem?
Thank you in advance!

math_robustSE_w <- list ()
dv= cbind(data$math1, data$math2, data$math3, data$math4, data$math5)
for(i in 1:5){
math_robustSE_w <- rlmer(dv[,i] ~ ZSC183 + ZTEAFDBK + ZSC017 + ZSC155 + ZMTTRAIN + ZTEACHBEHA + ZMCLSIZE +
ZSMRATIO + ZSC018Q01TA01 + ZSC018Q09JA01 + ZSC018Q10JA01 + ZSC025Q02NA +
ZSC025Q02NA + ZSC182Q08JA01 + ZSC182Q09JA01 + (1|CNTSCHID), data=data, REML = FALSE)
}
math_robustSE_w.pooled <- mitml::testEstimates(math_robustSE_w, extra.pars = TRUE)

Error in .extractParameters(model, include.extra.pars = TRUE) :
all(p == dim(Uhat[[1]])) is not TRUE