1
$\begingroup$

I am having trouble with analysing some of my data. I'm trying to test the effect of a treatment with two levels (Treatment/Control) on the abundance of ants belonging to different dominance categories on plants (Response = Abundance, Predictors = Site, Treatment, Dominance).

There are two sites (SiteC/SiteT) each containing 3 sets of paired plots )((3*treatment + 3*control)*2 Sites). There are 3 categories of ant dominance (D, S or O). Each plot contains 9 plants. Sample provided below.

 Site    Pair Plot Site_Plot Treatment   Check     Plant   Dominance Abundance Presence
1   SiteC SiteC_1   1C  SiteC_1C   Control Initial  Plant_73    Dominant         0        0
2   SiteC SiteC_1   1C  SiteC_1C   Control Initial  Plant_73 Subdominant         0        0
3   SiteC SiteC_1   1C  SiteC_1C   Control Initial  Plant_73 Subordinate         0        0
4   SiteC SiteC_1   1C  SiteC_1C   Control Initial  Plant_76    Dominant         0        0
5   SiteC SiteC_1   1C  SiteC_1C   Control Initial  Plant_76 Subdominant         0        0
6   SiteC SiteC_1   1C  SiteC_1C   Control Initial  Plant_76 Subordinate         1        1
308 SiteT SiteT_3   3T  SiteT_3T Treatment Initial  Plant_50    Dominant         0        0
309 SiteT SiteT_3   3T  SiteT_3T Treatment Initial  Plant_50 Subdominant         0        0
310 SiteT SiteT_3   3T  SiteT_3T Treatment Initial  Plant_50 Subordinate        10        1
311 SiteT SiteT_3   3T  SiteT_3T Treatment Initial  Plant_53    Dominant         0        0
312 SiteT SiteT_3   3T  SiteT_3T Treatment Initial  Plant_53 Subdominant         0        0
313 SiteT SiteT_3   3T  SiteT_3T Treatment Initial  Plant_53 Subordinate         1        1

The initial model was: Average_Abundance ~ Treatment * Site * Dominance + (1|Pair)

The issues are that:

  • Sites responded differently to the treatment and there are interactions.

  • Very low replicates

  • Data is zero-inflated, overdispersed and are counts

  • Low number of levels in categorical variables (Site = 2, Treatment = 2 and Dominance = 3)

One solution suggested to me was to use simple statistics, specifically paired t-tests to analyse this data. However the non-parametric equivalent Wilcoxon signed-rank tests cannot return significance for n=3.

Some solutions I am trying to use:

  • Instead of using the average abundance per plot, I have used plant as the replicate and nested plot in plot pairs to account for spatial autocorrelation.

  • use a ZINB GLMM (glmmTMB) instead of glmer.nb (lme4) to address zero-inflation and overdispersion.

However! This brings me to my problem and my questions.

The full model below will not run when ziformula=~. (when the zero-inflated model is the same as the conditional model).

zinb_p1 <- glmmTMB(Abundance ~ Site * Treatment * Dominance   + (1|Pair/Site_Plot),
                       data=ant_stats,
                       ziformula=~.,
                       family=nbinom1)

I get the following warning message

Warning message:
In fitTMB(TMBStruc) :
  Model convergence problem; non-positive-definite Hessian matrix. See vignette('troubleshooting')

If i change ziformula=~. to ziformula=~1 i get the following output, which I assume does not have a zero-inflated model fitted.

 > summary(zinb_p3)
 Family: nbinom1  ( log )
Formula:          Abundance ~ Site * Treatment * Dominance + (1 | Pair/Site_Plot)
Zero inflation:             ~1
Data: ant_stats

     AIC      BIC   logLik deviance df.resid 
   722.2    782.2   -345.1    690.2      297 

Random effects:

Conditional model:
 Groups         Name        Variance  Std.Dev. 
 Site_Plot:Pair (Intercept) 5.902e-11 7.682e-06
 Pair           (Intercept) 5.120e-02 2.263e-01
Number of obs: 313, groups:  Site_Plot:Pair, 12; Pair, 6

Overdispersion parameter for nbinom1 family (): 4.74 

Conditional model:
                                                  Estimate Std. Error z value Pr(>|z|)   
(Intercept)                                     -5.453e-01  4.745e-01  -1.149  0.25045   
SiteTWP                                          1.120e+00  5.559e-01   2.015  0.04394 * 
TreatmentTreatment                              -1.175e-01  6.301e-01  -0.186  0.85209   
DominanceSubdominant                             1.455e-01  6.271e-01   0.232  0.81657   
DominanceSubordinate                             1.412e+00  5.107e-01   2.765  0.00569 **
SiteTWP:TreatmentTreatment                      -5.310e-01  7.825e-01  -0.679  0.49741   
SiteTWP:DominanceSubdominant                    -2.991e+01  6.858e+05   0.000  0.99997   
SiteTWP:DominanceSubordinate                    -2.251e+00  6.934e-01  -3.246  0.00117 **
TreatmentTreatment:DominanceSubdominant         -5.825e-01  9.608e-01  -0.606  0.54436   
TreatmentTreatment:DominanceSubordinate          1.213e-01  7.122e-01   0.170  0.86474   
SiteTWP:TreatmentTreatment:DominanceSubdominant -1.786e+00  3.371e+06   0.000  1.00000   
SiteTWP:TreatmentTreatment:DominanceSubordinate  1.612e+00  9.631e-01   1.674  0.09421 . 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Zero-inflation model:
            Estimate Std. Error z value Pr(>|z|)
(Intercept)   -19.09    5788.09  -0.003    0.997

My Questions:

  1. Is this model overparametised? Is that why it won't work?

  2. If I want to continue to use plant as the replicate rather than plot but nesting plot in pair of plots causes working models to have convergence errors, what alternative is there to account for spatial autocorrelation?

  3. Is there a different approach I could be using that would allow me to keep all of my fixed variables?

$\endgroup$
  • $\begingroup$ Have you followed the suggestions in the vignette (troubleshooting) that is mentioned in the warning? The output of the first model that yields the warning would be helpful here, I think. $\endgroup$ – Daniel Nov 28 '18 at 12:57
0
$\begingroup$

Per the Introduction provided by the glmmTMB creators, ziformula = ~ 1 means that a zero-inflation parameter has been applied to all observations. ziformula = ~ 0 would exclude an adjustment for zero-inflation.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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