5
$\begingroup$

I have under-dispersed count data. I do not want to transform them, and using a negative binomial error distribution (via glmer.nb) does not help.

My results are the same regardless of the distribution I use, and agree very well with the conclusion that you would draw by eyeballing the data (it is pretty obvious that there is no effect, see below). I can not seem to find any real way to deal with this problem, and I am left wondering if the most appropriate thing to do is acknowledge the under-dispersion, point to the result that is robust to every method I've thrown at it, and move on.

If you've been in this situation- could you share any tools/advice/philosophy that you found helpful?

Histogram of counts

Here is the output from my Poisson model:

mc_p<-glmer(tubes ~ status + (1|site/plant), family="poisson", data=mc) summary(mc_p) Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod'] Family: poisson ( log ) Formula: tubes ~ status + (1 | site/plant) Data: mc

 AIC      BIC   logLik deviance df.resid 
713.7    725.7   -352.8    705.7      143 

Scaled residuals: 
      Min       1Q   Median       3Q      Max 
 -2.08950 -0.38627 -0.04237  0.33830  2.87921 

Random effects:
 Groups     Name        Variance Std.Dev.
 plant:site (Intercept) 0.01958  0.1399  
 site       (Intercept) 0.01790  0.1338  
Number of obs: 147, groups:  plant:site, 20; site, 4

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  2.25665    0.08417  26.809   <2e-16 ***
statusy      0.01356    0.05301   0.256    0.798    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
        (Intr)
statusy -0.308

I am trying determine whether pollen tube counts differ between nectar-robbed and un-robbed flowers. Pollen tube counts are nested within plant (multiple flowers of each type sampled from each plant) and plants are nested within site. I've tried using a negative binomial error distribution (as I described in the first iteration of this post), but that still gave me a very high deviance relative to my residual df. It seems that there isn't a straight-forward way to apply a useful error distribution to my nested data. At the same time, the result (that there is no difference between groups) is the same regardless of how I arrange my model, and can be seen very clearly in the bar-graph above.

|cite|improve this question
$\endgroup$
  • $\begingroup$ I have a possibly silly question: why do you think your data are under-dispersed? $\endgroup$ – atiretoo - reinstate monica Dec 10 '15 at 1:51
2
$\begingroup$

If the data are overdispersed, you can estimate the same relative rates and calculate prediction intervals as in a poisson model using a quasipoisson model. Even if the data are not overdispersed, a quasipoisson model is valid and fairly efficient. A quasipoisson model just extends the poisson by estimating a dispersion parameter.

It is, of course, important to think about the originating nature of the data, the question you're trying to answer, and to ask about why this dispersion comes about, and how it might affect your interpretation of the results.

|cite|improve this answer
$\endgroup$
  • $\begingroup$ Thanks @AdamO. My data are actually under-dispersed, does this advice still apply? I've tried fitting them using glmer.nb, but it always hits the iteration limit, and I'm not sure what (if anything) I can do about that. IN regard to the nature of the data- I feel that it makes perfect sense that the data are distributed in the way that they are, and that there is no difference between the groups. $\endgroup$ – JKO Dec 9 '15 at 19:02
  • $\begingroup$ @BioUser Do you have any rationale for choosing a negative binomial probability model for these data? There are a lot of assumptions behind such a model, many idiosyncrasies compared to exponential family GLMs, and few tests of those assumptions. Have you tried fitting a poisson model? $\endgroup$ – AdamO Dec 9 '15 at 19:25
  • 1
    $\begingroup$ @BioUser Over and under dispersion are basically the same thing: they show that the mean-variance relationship is not entirely correct. Allowing for dispersion means assuming that the mean-variance relationship is at least proportional to the one assumed under the usual fully parametric GLM model. In Poisson models: var=mean. Negative binomial models are more complicated. Note well: what you are calling underdispersion may not necessarily be addressed by using a model that allows for dispersion. It could be evidence of several things. $\endgroup$ – AdamO Dec 9 '15 at 19:28
  • $\begingroup$ AdamO, thanks again! Yes- I've tried fitting a Poisson model (which looks the most appropriate to me). How would I sort out whether I truly have underdispersion, and whether or not I should allow for dispersion? I've edited my question to include the output for my Poisson model. $\endgroup$ – JKO Dec 9 '15 at 19:54
  • $\begingroup$ @BioUser Think of it this way: not all eggs contain salmonella, but you rarely worry about testing raw eggs for salmonella, you just use them in recipes where they are cooked. If the data are not over/under dispersed then, quasipoisson is a recipe which cooks the eggs. It tastes the same either way. $\endgroup$ – AdamO Dec 9 '15 at 20:31
2
$\begingroup$

A descriptively adequate treatment of under-dispersion can apparently be gotten from regression using a Conway Maxwell Poisson (COM) distribution or Consul's Generalized Poisson (GP) regression. It seems that just modeling with a Normality assumption is inefficient.

A review paper and responses for the COM model are here and an R package to implement it is here. A Stata package with references to the GP regression model is described here.

Full disclosure: I've never used these things in anger.

|cite|improve this answer
$\endgroup$
  • $\begingroup$ Thanks @conjugateprior. However, I'm under the impression that COM-Poisson can't handle nested data. Is this incorrect? $\endgroup$ – JKO Dec 9 '15 at 20:05
  • $\begingroup$ The COMPoisson is just being used as a conditional distribution, so it'd be the regression part of the package that failed to deal with nesting. It's possible there's another implementation that does that. The alternatives would be to be a) write out the likelihood function with nesting and maximize it with an optimizer (rather you than me), or b) write out the same model, add some sensible priors and run a multilevel model in JAGS (possibly more me than you), or c) Block bootstrap the standard errors to deal with the nesting (perhaps the more straightforward option). $\endgroup$ – conjugateprior Dec 9 '15 at 20:23

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.