Solved – lme4: Problems trying to run negative binomial glmer.nb

binomial distributionlme4-nlmenegative-binomial-distributionproportion;r

this is my first post and I am very new to the world of R computing so I hope I have asked everything in the right format, included enough information etc.

My data set consists of proportion data for behavioural observations of beetles on a carcass therefore are between 0 and 1. I want to test these against my fixed effects (Call, Population, Pronotum width) as well as random effects (Focal Male Id and Mouse order)

model1 <- glmer(cbind(Total.ON, Not.ON) ~ Call + Population + Pronotum.width.mm  + (1|Focal.Male.ID) 
           + (1|Mouse.Order), family = binomial, data = Beetle.data)

Generalized linear mixed model fit by maximum likelihood (Laplace 
Approximation) ['glmerMod']
Family: binomial  ( logit )
Formula: cbind(Total.ON, Not.ON) ~ Call + Population + Pronotum.width.mm +      
(1 | Focal.Male.ID) + (1 | Mouse.Order)
Data: Beetle.data

 AIC      BIC   logLik deviance df.resid 
299.2    329.2   -141.6    283.2      306 

Scaled residuals: 
Min      1Q  Median      3Q     Max 
-2.4854 -0.2559 -0.1704  0.4027  1.8366 

Random effects:
 Groups        Name        Variance Std.Dev.
 Focal.Male.ID (Intercept) 1.67079  1.29259 
 Mouse.Order   (Intercept) 0.00104  0.03224 
Number of obs: 314, groups:  Focal.Male.ID, 152; Mouse.Order, 2

Fixed effects:
                Estimate Std. Error z value Pr(>|z|)    
(Intercept)          -5.3314     2.8099  -1.897 0.057779 .  
CallNo-competitor     1.4068     0.3818   3.685 0.000229 ***
PopulationGlos        0.1858     0.5889   0.316 0.752354    
PopulationLab         0.2099     0.6192   0.339 0.734691    
PopulationNewbottle   0.2181     0.5811   0.375 0.707435    
Pronotum.width.mm     0.4534     0.5551   0.817 0.413996    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
        (Intr) CllN-c PpltnG PpltnL PpltnN
CllN-cmpttr -0.160                            
PopulatnGls -0.155  0.021                     
PopulatinLb -0.276  0.022  0.487              
PpltnNwbttl -0.050  0.012  0.504  0.473       
Prntm.wdth. -0.982  0.053  0.047  0.177 -0.055
convergence code: 0
Model failed to converge with max|grad| = 1.61379 (tol = 0.001, component 1)

After running this model, i tested the overdispersion using the overdisp_fun function from the GLMM page and got these ratios

> overdisp_fun(model1)
 chisq      ratio        rdf          p  
87.310363   0.285328 306.000000   1.000000 

After this I thought I would try fitting a negative binomial model to my data to see if this improves the dispersion ratios but it is here I run into problems.

I have tried using the glmer.nb in the lme4 package but keep getting error messages

model1 <- glmer.nb(cbind(Total.ON, Not.ON) ~ Call + Population + Pronotum.width.mm  + (1|Focal.Male.ID) 
           + (1|Mouse.Order), data = Beetle.data)

Error in (function (fr, X, reTrms, family, nAGQ = 1L, verbose = 0L, maxit = 100L,  : updateXwts: dimension mismatch

If anyone has any ideas or recommendations for my model in general it would be very appreciated.

Thanks in advance

Best Answer

There a several things going on here.

  • It's not worth worrying about any aspects of diagnostics or model fit until you sure you've specified your model (and in this particular case your response variable) correctly. For binomial GLM(M)s in R, the response matrix should be cbind(num_successes,num_failure). Is your Total.ON the total number of cases measured for each observation (i.e. successes + failures)? If so, you should use cbind(Total.ON-Not.ON,Not.ON) ...

    if (as mentioned in comments) Total.ON is a proportion successful and Not.ON is a total number of observations, you should use either cbind(Total.ON*Not.ON,(1-Total.ON)*Not.ON) or (my preference) use Total.ON as your response and include weights=Not.ON as an additional argument.

  • the results of your overdispersion test:
chisq      ratio        rdf          p  
87.310363   0.285328 306.000000   1.000000 

actually suggest underdispersion; the ratio is less than one and the p-value is 1.00

  • It's a very common (and understandable) confusion to think that the negative binomial is an extension of the binomial distribution that's appropriate for accounting for overdispersion. Actually, the negative binomial extends the Poisson distribution. The section on overdispersion in the GLMM FAQ suggests various methods for dealing with overdispersion in binomial models: observation-level random effects (= logit-Normal-binomial models), beta-binomial models (e.g. in glmmTMB), quasi-likelihood, ...
  • if you do end up having underdispersion, you can use quasi-likelihood (at present I don't know of an underdispersed binomial analogue); there is a new section in the GLMM FAQ on underdispersion