My study has a complicated design and I am not sure if I am modeling my zero-inflated data correctly. I have seed abundances and seedling abundances for 11 species. I have one main "treatment" with four levels (C, P, I, R). For the sampling design: we sampled at 5 sites. At each site, each of the 4 treatments were replicated 4 times (e.g. PA, PB, PC, PD). So there are 4 replicates of 4 treatments nested in 5 sites (so Site is my random factor). The response variable (one species' seed abundances) are zero-inflated (count data). I checked for overdispersion by running a glm and dividing the residual deviance by the degrees of freedom (I believe this is the right approach?). All variables were overdispersed (ratio greater than 1); so I went with negative binomial, instead of Poisson, distributions. I've tried two ways of modeling this: zeroinfl and glmmadmb. But after learning zeroinfl is not useful for mixed models, I am trying to model with glmmADMB:

glmmNB <- glmmadmb(CON_XAL~Treatment+(1|Site), data = SR.year.raw, 
                   zeroInflation=TRUE, family="nbinom")

My outcome looks like this:

enter image description here

Summary of my response variable:

 Min.   1st Qu.    Median      Mean   3rd Qu.      Max.   
(0.0)     (0.0)     (0.0)   (232.1)     (7.5)  (6245.0)    

My questions are:

  1. Am I specifying the random effect correctly in the glmmadmb model?
  2. What post-hoc tests can I do with glmmADMB to look at how treatments differ from one another? For example, I have a treatment called "control". But in the model output, there is no significance value for the control. So I suppose the output is saying that the "island" treatment is significantly different from the control. But how do I know if the "Island" treatment is significantly different from the "plantation" treatment?
  3. Can you recommend some ways of graphing the results, i.e., something I ultimately could include in my paper? I ask this because from everything I have read, all of the example graphs are for comparing different models (Poisson vs ngbinom), but I can't seem to find code for a good final conclusions graph.
  4. Is the warning message a problem? Can I ignore it?
  • $\begingroup$ Your question may be flagged as off-topic since it involves program specific questions, but let me add: a quick look at the zeroinfl documentation shows that it is not designed for mixed model/longitudinaly analysis. The way you think you are specifying the random effects is actually specifying a different model for each component of the mixture model. It looks like glmmADMB is your best choice. There is a vignette for that package showing model testing and plots. Not an expert on this, but hope it helps. $\endgroup$
    – Moose
    Apr 26, 2016 at 15:59
  • $\begingroup$ If you use the save.dir option of glmmadmb and make the files glmmadmb.dat and glmmadmb.pin available, I'll take a look and see what if anything useful I can say. In general, zero inflaiton and overdispersion are often confounded for these models, so that when used together it might appear that neither is significant, while if you use either one by itself it may be that one or both is significant -- that's life! I would use likelihood ratio tests to investigate this. $\endgroup$ Apr 26, 2016 at 16:07
  • $\begingroup$ Moose, thank you. Yes, the zeroinfl models were recommended by a supervisor and I see now that they do not handle mixed models. So I will focus now on glmmadmb. Dave Fournier, I have the glmmadmb.dat and glmmadmb.pin files, but how do I make them available here in this forum? Thank you for the help everyone. $\endgroup$
    – Sylvia
    May 1, 2016 at 18:54

1 Answer 1

  1. No, zeroinfl() currently does not support random effects. So the formula you specified actually means something different: You use a fixed treatment effect in the count part and a fixed site effect in the zero-inflation part. See vignette("countreg", package = "pscl") for more details.

  2. If you want random effects, then no. If you use fixed interaction effects instead, you could still try to find a suitable model with zeroinfl(). But with your number of observations this is probably not the best solution.

  3. As the model is not the one you would want to fit, this is not relevant here.

  4. For zeroinfl() there would be and I suppose that for glmmADMB there are as well. But I'm not an expert on that.

  5. You could employ effect plots for the covariate effects or rootograms for the goodness of fit. It depends on what you really want to show.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.