Random effect equal to 0 in generalized linear mixed model [duplicate]

Sorry if I’m missing something very obvious here but I am new to mixed effect modelling.

I am trying to model a binomial presence/absence response as a function of percentages of habitat within the surrounding area. My fixed effect is the percentage of the habitat and my random effect is the site (I mapped 3 different farm sites).

glmmsetaside <- glmer(treat~setas+(1|farm),
       family=binomial,data=territory)

When verbose=TRUE:

0:     101.32427: 0.333333 -0.0485387 0.138083 
1:     99.797113: 0.000000 -0.0531503 0.148455  
2:     99.797093: 0.000000 -0.0520462 0.148285  
3:     99.797079: 0.000000 -0.0522062 0.147179  
4:     99.797051: 7.27111e-007 -0.0508770 0.145384  
5:     99.797012: 1.45988e-006 -0.0495767 0.141109  
6:     99.797006: 0.000000 -0.0481233 0.136883  
7:     99.797005: 0.000000 -0.0485380 0.138081  
8:     99.797005: 0.000000 -0.0485387 0.138083  

My output is this:

Generalized linear mixed model fit by the Laplace approximation 
Formula: treat ~ setasidetrans + (1 | farm) 

AIC   BIC logLik deviance
105.8 112.6  -49.9     99.8
Random effects:
 Groups Name        Variance Std.Dev.
farm   (Intercept)  0        0  
Number of obs: 72, groups: farm, 3

Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept)   -0.04854    0.44848  -0.108    0.914
setasidetrans  0.13800    1.08539   0.127    0.899

Correlation of Fixed Effects:
            (Intr)
setasidtrns -0.851

I basically do not understand why my random effect is 0? Is it because the random effect only has 3 levels? I don’t see why this would be the case. I have tried it with lots of different models and it always comes out as 0.

It cant be because the random effect doesn’t explain any of the variation because I know the habitats are different in the different farms.

Here is an example set of data using dput:

list(territory = c(1, 7, 8, 9, 10, 11, 12, 13, 14, 2, 3, 4, 5, 
6, 15, 21, 22, 23, 24, 25, 26, 27, 28, 16, 17, 18, 19, 20, 29, 
33, 34, 35, 36, 37, 38, 39, 40, 30, 31, 32, 41, 45, 46, 47, 48, 
49, 50, 51, 52, 42, 43, 44, 53, 55, 56, 57, 58, 59, 60, 61, 62, 
54, 63, 65, 66, 67, 68, 69, 70, 71, 72, 64), treat = c(1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0), farm = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 
3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3), 
built = c(5.202332763, 1.445026852, 2.613422283, 2.261705833, 
2.168842186, 1.267473928, 0, 0, 0, 9.362387965, 17.55433115, 
4.58020626, 4.739300829, 8.638442377, 0, 1.220760647, 7.979990338, 
13.30789514, 0, 8.685544976, 3.71617163, 0, 0, 6.802926951, 
8.925512803, 8.834006678, 4.687723044, 9.878232478, 8.097800267, 
0, 0, 0, 0, 5.639651095, 9.381654651, 8.801754791, 5.692392532, 
3.865304919, 4.493438554, 4.826277798, 3.650995554, 8.20818417, 
0, 8.169597157, 8.62030666, 8.159474015, 8.608979238, 0, 
8.588288678, 7.185700856, 0, 0, 3.089524893, 3.840381223, 
31.98103158, 5.735501995, 5.297691011, 5.17141191, 6.007539933, 
2.703345394, 4.298077606, 1.469986793, 0, 4.258511595, 0, 
21.07029581, 6.737664009, 14.36176373, 3.056631919, 0, 32.49289428, 
0)

It goes on with around 10 more columns for different types of habitat (like built, setaside is one of them) with percentages in it.

Answer

With just three farms, there is no point in trying to pretend that you can fit a Gaussian distribution to three points. Analyze this simply as lm(response~as.factor(farm) + treat+other stuff), and don’t bother with lmer; you won’t be able to do much better than ANOVA, anyway.

Generally, hitting exactly zero is not that unusual. The variance estimate is a nonlinear function of the data, the difference between the overall variance and the within-site variance. If the true variance is zero, this nonlinear statistic has a distribution that puts non-zero mass to the left of zero (this will also be true if the true value is a small positive quantity, but the sampling variability is large enough to overshoot below zero). Due to the way the estimator is programmed, however (Cholesky factorization), it can only take non-negative values. So whenever the unattainably best estimate would have been at zero (as in your balanced-by-design situation) or below it, the log-likelihood will be maximized at zero, with a negative gradient to the right of it. Self & Liang (1987) is the standard biostat reference for the problem; I better like Andrews (1999) which is even more general.

Attribution
Source : Link , Question Author : Cec.g , Answer Author : Ben Bolker

Leave a Comment