What is causing the singularity in a glmm with simple random effect?

Question

First time poster but have been very grateful over the past couple months for this forum. First and foremost, I apologize in advance if I am not following the right procedures in asking a question.

I am searching for why my models come out singular.

I have 9 farms, and am treating farm as a random effect. I have 3 treatments in each farm with 2 replicates for each treatment (for most of the variables)

My response variable is count data on worms and due to overdispersion I am using a negative binomial distribution.

I have a number of explanatory variables that I reduced be testing for collinearity via spearman coefficients (>0.50) and making decisions on which to keep via biological importance, as well as a few removed that were clearly linearly dependent on one another. I also centered and scaled them.

From a full model of these variables, I built out approx. 120 models with glmer.nb() from lme4 and used model.sel() from the MuMIn package to determine the best model(s).

Unfortunately, many of my (better) models are returning the error "?isSingular" And in checking via

tt<-getME(mod1, "theta")
ll<-getME(mod1, "lower")
min(tt[ll==0])

or

theta <- getME(mod1,"theta")
diag.element <- getME(mod1,"lower")==0
any(theta[diag.element]<1e-5)

the models seem truly singular.

In reading through some stackexchange questions/answers and as well as some good websties (including the following https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#what-methods-are-available-to-fit-estimate-glmms) I've decided for my purposes, getting rid of the only random effect term does not make sense and I will not attempt to pursue a Bayesian approach.

But I would really like to understand why some of my models are singular and others are not.

What I've noted thus far:

• Out of the 8 models with singularity problems, they all contain 4 of the same variables , one of these variables SNH_perc, has the same value throughout each farm due to it being measured at a kilometer radius from the center of the farm. I tried creating a little bit of noise in these values but this did not fix the singularity.

• In using VarCorr(mod1) those with singularity problems have an extremely small value 1e-6 or so, while those without have about 0.1, indicating these models (often with more variables?) have less variation when more dimensions (or certain dimensions) are added

• I got the eigenvalues of the Variance/Covariance Choleski Decomposition and they seem off but I am not sure how to interpret them either

an example for one model being

mod1 <- glmer.nb(abundance_SUM ~ SOM_perc + P2O5 + Cu + silt_percent + SNH_perc + (1|vineyard), data=firstrun.scale)

    eigen() decomposition
    

$values
[1] 0.07609357, 0.07371038, 0.06832976, 0.06689690, 0.06391916, 0.06232956


$vectors
           [,1]      [,2]        [,3]       [,4] [,5]         [,6]
[1,] -0.2138548 0.2237181  0.64555958 0.77197345    1  0.922543154
[2,] -0.8997568 0.9746539  0.07141780 0.04855806    0 -0.007181145
[3,] -0.3617484 0.0000000  0.75959292 0.63379737    0 -0.380728292
[4,]  0.1176516 0.0000000 -0.01901948 0.00000000    0 -0.047982131
[5,]  0.0000000 0.0000000  0.02844624 0.00000000    0  0.034387977
[6,]  0.0000000 0.0000000  0.00000000 0.00000000    0  0.020584223

a sample of my dataset with the relevant explanatory variables (I'm not sure how much of it is useful to provide).

  abundance_SUM farm treatment SOM_perc    P2O5    Mg     Cu silt_percent SNH_perc SNH_perc2
1              8 AT01        BG 7.598926 26.1835 13.40   1.66     70.65202     55.6      55.6
2              4 AT01        BG 7.964178 32.4680 14.00   1.68     71.18537     55.6      55.6
3             22 AT01        AC 7.848078 23.1355 13.40   1.66     70.92315     55.6      55.5
4             14 AT01        AC 8.192078 28.0170 15.95   1.51     71.15727     55.6      55.5
5             20 AT01        CC 8.058861 28.0720 13.65   1.46     67.84648     55.6      55.4
6             12 AT01        CC 8.239543 28.0400 13.25   2.06     70.80634     55.6      55.4
7              7 AT02        BG 6.622578 80.4045 13.35   7.04     57.07563     35.7      35.7
8             10 AT02        BG 6.962443 72.7000 14.30   4.89     57.49698     35.7      35.7
9             13 AT02        AC 7.263278 62.7370 13.95   6.71     61.04811     35.7      35.8
10            11 AT02        AC 6.953843 63.4000 12.15   7.68     61.04811     35.7      35.8
11            13 AT02        CC 7.258978 62.4965 14.40  54.15     65.65210     35.7      35.9
12            10 AT02        CC 7.074078 59.1835 12.75  52.07     65.65210     35.7      35.9
13            46 AT03        BG 4.194980 12.3990 20.05  10.71     46.27264     59.7      59.7
14            46 AT03        BG 4.704530 15.9010 20.50  10.11     46.27264     59.7      59.7
15            28 AT03        AC 5.749596 20.0230 17.70   8.04     47.82043     59.7      59.6
16            34 AT03        AC 5.874213 22.9345 18.70  10.09     47.82043     59.7      59.6
17            33 AT03        CC 5.169096 17.9850 21.50   9.59     46.08007     59.7      59.8
18            34 AT03        CC 5.267913 18.4005 23.25  10.52     46.08007     59.7      59.8
19            26 AT04        BG 5.865696 11.0740 10.15  15.10     42.56881      7.7       7.7
20            42 AT04        BG 5.753896 11.2875 10.15  12.89     42.56881      7.7       7.7
21            16 AT04        AC 3.861730 24.7900 13.45  29.98     34.74912      7.7       7.8
22            34 AT04        AC 3.268413 24.7310 11.65  30.02     34.74912      7.7       7.8
23            10 AT04        CC 4.936813 21.2350 10.70  21.65     34.88536      7.7       7.6
24             8 AT04        CC 4.807813 20.0090 11.20  21.28     34.88536      7.7       7.6
25             3 AT05        BG 4.395096 37.5385 13.70  78.91     59.70392     12.3      12.3
26             5 AT05        BG 4.863630 41.0965 13.75  75.80     59.70392     12.3      12.3
27             4 AT05        AC 5.839730 31.6870 14.70  88.75     62.55672     12.3      12.4
28             9 AT05        AC 5.237813 43.3220 12.60  80.68     62.55672     12.3      12.4
29             2 AT05        CC 5.710813 32.4700 14.95 106.10     58.84137     12.3      12.5
30             9 AT05        CC 4.885213 34.1055 12.00  99.67     58.84137     12.3      12.5

Thank you in advance for any hints in the right direction!

EDIT: Adding the summary output for both glm and glmm version of best model

Call:
MASS::glm.nb(formula = abundance_SUM ~ SOM_perc + P2O5 + Cu + 
    silt_percent + SNH_perc, data = firstrun.scale, init.theta = 6.692239628, 
    link = log)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-2.29715  -0.91392  -0.00525   0.60202   2.14098  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)   2.75125    0.06393  43.038  < 2e-16 ***
SOM_perc     -0.18855    0.07293  -2.585 0.009730 ** 
P2O5         -0.20610    0.06948  -2.966 0.003013 ** 
Cu           -0.18915    0.08093  -2.337 0.019433 *  
silt_percent -0.25421    0.07006  -3.628 0.000285 ***
SNH_perc      0.23310    0.07741   3.011 0.002603 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(6.6922) family taken to be 1)

    Null deviance: 108.025  on 53  degrees of freedom
Residual deviance:  55.154  on 48  degrees of freedom
AIC: 375.88

Number of Fisher Scoring iterations: 1


              Theta:  6.69 
          Std. Err.:  1.85 

 2 x log-likelihood:  -361.88 



Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: Negative Binomial(6.6922)  ( log )
Formula: abundance_SUM ~ SOM_perc + P2O5 + Cu + silt_percent + SNH_perc +      (1 | vineyard)
   Data: firstrun.scale

     AIC      BIC   logLik deviance df.resid 
   377.9    393.8   -180.9    361.9       46 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.66497 -0.80737 -0.00517  0.65162  2.77370 

Random effects:
 Groups   Name        Variance  Std.Dev. 
 vineyard (Intercept) 3.389e-11 5.821e-06
Number of obs: 54, groups:  vineyard, 9

Fixed effects:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)   2.75126    0.06392  43.043  < 2e-16 ***
SOM_perc     -0.18855    0.07374  -2.557 0.010564 *  
P2O5         -0.20610    0.06699  -3.077 0.002093 ** 
Cu           -0.18915    0.08376  -2.258 0.023925 *  
silt_percent -0.25421    0.07260  -3.501 0.000463 ***
SNH_perc      0.23310    0.07526   3.097 0.001952 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
            (Intr) SOM_pr P2O5   Cu     slt_pr
SOM_perc     0.030                            
P2O5         0.052 -0.006                     
Cu           0.068 -0.235 -0.332              
silt_percnt  0.083 -0.174  0.272  0.022       
SNH_perc    -0.044 -0.406 -0.182  0.437 -0.084
optimizer (Nelder_Mead) convergence code: 0 (OK)
boundary (singular) fit: see ?isSingular

Answer 1

It is evident from the ouput that there is extremely small variation within vineyard. Note that variance of the random intercepts is 3.389e-11 which is basically zero. This is why you have a singular fit.

Also note that the output from the glm.nb model is basically the same as the output for the fixed effects from the glmm.

So you do not need random effects for this model.

As for the procedure of "information-theoretic approach using AICc values" I don't know what thay is, but when you fit hundreds of models to a dataset there is a distinct possibility that you overfit, and end up with a model that does not generalise well to new data.