[R-sig-ME] random effects with nbinom glmmadmb

[Ben Bolker]

On 13-12-20 07:52 AM, Jurij Diaci BF wrote:
> Hi,
> 
> I model emergence of seedling counts (sp: spruce) in four canopy openings
> created by a bark beetle calamity (random factor opening). Two openings were
> fenced and two were not fenced (factor fence). Within openings plots were
> placed in a stratified random manner to sample two different orographic
> features: plateau and sinkhole (factor relief) and two positions within
> canopy opening  (factor location). Sampling was repeated after seven years
> (factor year). I sampled also several covariates (c1: rockiness; c2: ground
> vegetation coverage).
> 
> The total dataset consist of: 2(fence)*2(relief)*2(location)*2(year)*30 =
> 480 plots.
> I applied a negative binomial model:
> 
> nbm.sp <- glm.nb(sp  ~  fence  + relief + location + year+c1 +c2, link =
> "log", data = krese2) Results seem OK, but I'm aware of pseudoreplication.

> 
> Therefore, I applied a mixed negative binomial model:
> 
> mnbm.sp <- glmmadmb(sp ~ fence  + relief + location + year + location:year +
>                       (1 | opening), family="nbinom", data = krese2) 

> Results are similar compared to the first model, with slightly larger std.
> errors.
> 
> I understand that 4 openings are at the lower limit for calculation of
> random effects.
> 
> Do the results below look OK?
> 
> Is this a correct approach? Namely, there is a fixed factor "fence"
> connected to openings:
> 
> table(krese2$fence, krese2$opening)

>          opening
>          a   b   c  d
> fenced   0   0 120 120
> nonfenc 40 200   0   0
> 
> Call:

> glmmadmb(formula = sp ~ fence  + relief + location + year + location:year +
>                       (1 | opening), data = krese2, family = "nbinom")
> AIC: 2933.7
> 
> Coefficients:
> 
>                 Estimate Std. Error z value Pr(>|z|)
> (Intercept)   1.588      0.194    8.20  2.4e-16 ***
> reliefsink    0.243      0.170    1.43  0.15350
> fencefe       0.488      0.211    2.31  0.02082 *
> locationedge  0.785      0.172    4.56  5.0e-06 ***
> year2013     -0.408      0.118   -3.46  0.00054 ***
> 
> reliefsink:locationedge   -0.541      0.240   -2.26  0.02413 *
> Number of observations: total=480, opening=4 Random effect variance(s):
> Group=Odd
>             Variance StdDev
> (Intercept)  0.02702 0.1644
> Negative binomial dispersion parameter: 0.66844 (std. err.: 0.048389)
> Log-likelihood: -1458.86


  This looks reasonable overall.  Some comments:

 * Why do your fixed effects differ between models?  You have the
location:year interaction in the GLMM and not the GLM, and the covariates
in the GLM but not the GLMM.  (But I don't see the location:year parameter
in the output ... is that really location:relief ??)
 * It would be nice to compare the GLM with 'opening' as a fixed
effect, but you can't quite do this -- the fence effect will be partially
confounded with opening ...
 * Be aware that the p-values quoted are Wald effects, so don't take
account of finite-size effects. It's a bit hard to get the values
exactly, but your fence effect would require 8 df to be significant
at the 0.05 level (qt(0.975,8)=2.306) ... so I would be quite careful
with that p-value.  (You might be able to use MCMC to get more reliable
values -- I believe the glmmADMB vignette has an example.)