fitting a GLMER model with glmmTMB

Question

I'm trying to fit a generalized linear mixed model with glmmTMB

summary(glmmTMB(cbind(SARA_ph58, 1)~ `Melk T`+VetT+EiwT+
    `VET EIWIT ratio`+LactT+CelT+UrmT+vetg+eiwitg+lactg+
    `DS opname`+`boluses per day`+`chewings per bolus`+
    `rumination (min/d)`+ Activiteit + (1|experiment),
     data = dataset1geenNA, family = binomial()))

When i run this code i get some output but is also get the next warning message:

1: In fitTMB(TMBStruc) : Model convergence problem; non-positive-definite Hessian matrix. See vignette('troubleshooting')
2: In sqrt(diag(vcov)) : NaNs produced

does anybody know how to solve this problem?

Output:

 Family: binomial  ( logit )
Formula:          cbind(SARA_ph58, 1) ~ `Melk T` + VetT + EiwT + `VET EIWIT ratio` +      LactT + CelT + UrmT + vetg + eiwitg + lactg + `DS opname` +  
    `boluses per day` + `chewings per bolus` + `rumination (min/d)` +      Activiteit + (1 | experiment)
Data: dataset1geenNA

     AIC      BIC   logLik deviance df.resid 
      NA       NA       NA       NA       79 

Random effects:

Conditional model:
 Groups     Name        Variance  Std.Dev. 
 experiment (Intercept) 5.138e-08 0.0002267
Number of obs: 96, groups:  experiment, 3

Conditional model:
                       Estimate Std. Error z value Pr(>|z|)  
(Intercept)          -1.595e+01  1.605e+01  -0.994   0.3202  
`Melk T`             -2.560e-01  1.330e-01  -1.925   0.0542 .
VetT                 -7.499e+00  3.166e+00  -2.369   0.0178 *
EiwT                  8.353e+00  4.885e+00   1.710   0.0872 .
`VET EIWIT ratio`     2.100e+01  1.545e+01   1.359   0.1742  
LactT                -2.086e+00  8.571e-01  -2.434   0.0149 *
CelT                 -1.430e-04  6.939e-04  -0.206   0.8367  
UrmT                  1.300e-02  3.978e-02   0.327   0.7438  
vetg                  1.166e-03         NA      NA       NA  
eiwitg               -2.596e-03  5.180e-03  -0.501   0.6162  
lactg                 7.862e-03         NA      NA       NA  
`DS opname`          -1.882e-02  8.416e-02  -0.224   0.8231  
`boluses per day`    -3.200e-02  1.226e-01  -0.261   0.7940  
`chewings per bolus`  1.758e-02  6.688e-02   0.263   0.7927  
`rumination (min/d)` -1.468e-03  3.145e-03  -0.467   0.6408  
Activiteit            4.265e-03  4.625e-03   0.922   0.3564
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

We need a lot more information to help you solve this. The development version of glmmTMB has a diagnose function that might help you: see glmmtmb.github.io/glmmTMB for instructions about how to diagnose. There is a basic troubleshooting vignette (the one pointed to in the warning message) here. Can you at least show us the output? — Ben Bolker
@BenBolker i added a picture of the whole output, let me know if it readable for you. — PietroFleer
Thanks, but if you can include the output as text (i.e. cut-and-paste from your R session) rather than as a screenshot that would be much better. Put it in a code block; if you don't know how to do that, just put the text in and I or someone else will edit. — Ben Bolker

Ben Bolker Ben Bolker · Accepted Answer · 2021-05-27T16:21:15

There are a number of issues here.

The proximal problem is that you have a (near) singular fit: glmmTMB is trying to make the variance zero (5.138e-08 is as close as it can get). Because it fits on a log-variance (actually log-standard-deviation) scale, this means that it's trying to go to -∞, which makes the covariance matrix of the parameters impossible to estimate.

The main reason this is happening is that you have a very small number of groups (3) in your random effect (experiment).

These are extremely common issues with mixed models: you can start by reading ?lme4::isSingular and the relevant section of the GLMM FAQ.

The simplest solution would be to treat experiment as a fixed effect, in which case you would no longer have a mixed model and you could back to plain glm().

Another slightly worrying aspect of your code is the response variable cbind(SARA_ph58, 1). If SARA_ph58 is a binary (0/1) variable you can use just SARA_ph58. If you pass a two-column matrix as you are doing, the first column is interpreted as the number of successes and the second column as the number of failures; it looks like you may have been trying to specify that the total number of trials for each observation is 1 (again, if this is the case you can just use SARA_ph58 as the response).

One final note is that lme4::glmer is a little more tolerant of singular fits than glmmTMB.

fitting a GLMER model with glmmTMB

1 Answers