I have a dataset which seems to have a lot of zeroes. I have already fit a poisson regression model as well as a negative binomial model. I would like to fit zero-inflated and hurdle models as well.

Before I do I would like to run a test to investigate whether my data really is zero inflated. What test(s) is/are there to determine whether my data are zero-inflated?

  • $\begingroup$ Why did you add the SAS and Stata tags if in your question you only ask about R? $\endgroup$ – Maarten Buis Oct 8 '14 at 11:19
  • $\begingroup$ I would look at the zero-inflated Poisson mixture model paper by Dianne Lambert, and use something like AIC, BIC, or -Log-likelihood to compare inflated and non-inflated models. (tandfonline.com/doi/abs/10.1080/00401706.1992.10485228) $\endgroup$ – EngrStudent Feb 20 '17 at 14:25

The score test (referenced in the comments by Ben Bolker) is performed by first calculating the rate estimate $\hat{\lambda}= \bar{x}$. Then count the number of observed 0s denoted $n_0$ and the total number of observations $n$. Calculate $\tilde{p}_0=\exp[-\hat{\lambda}]$. Then the test statistic is calculated by the formula: $\frac{(n_0 - n\tilde{p}_0 )^2}{n\tilde{p}_0(1-\tilde{p}_0) - n\bar{x}\tilde{p}_0^2}$. This test statistic has a $\chi^2_1$ distribution which can be looked up in tables or via statistical software.

Here is some R code that will do this:

pois_data <-rpois(100,lambda=1)
lambda_est <- mean(pois_data)

p0_tilde <- exp(-lambda_est)
n0 <- sum(1*(!(pois_data >0)))
n <- length(pois_data)

# number of observtions 'expected' to be zero

#now lets perform the JVDB score test 
numerator <- (n0 -n*p0_tilde)^2
denominator <- n*p0_tilde*(1-p0_tilde) - n*lambda_est*(p0_tilde^2)

test_stat <- numerator/denominator

pvalue <- pchisq(test_stat,df=1, ncp=0, lower.tail=FALSE)
  • $\begingroup$ This only works in the non-regression setting, i.e. when there are no explanatory variables and the poisson mean is constant. $\endgroup$ – Josef Dec 10 '20 at 16:24
  • $\begingroup$ @Josef I guess that was were the approach started. The result would be asymptotic but I don't see why you couldn't put in predictions from GLMs if you're using the Lambert model, logistic regression for $\{0,1\}$s and Poisson regression for (non-structural 0) counts. Do you see a reason-besides needing large enough data size-why that approach would not work? $\endgroup$ – Lucas Roberts Dec 11 '20 at 2:52
  • $\begingroup$ p_0 depends on regression parameters and varies across observations and variance of the test statistic is more complicated taking estimation of regression parameters into account. I found this because I'm currently working on it for statsmodels e.g. github.com/statsmodels/statsmodels/issues/7131 $\endgroup$ – Josef Dec 11 '20 at 7:52
  • $\begingroup$ In R I only found van der Broeck test in countreg and vcdExtra, A moment test similar to your score_test works better when poisson means can be large. my notes github.com/statsmodels/statsmodels/issues/7131 $\endgroup$ – Josef Dec 11 '20 at 7:57
  • $\begingroup$ @Josef that makes sense, in one local space of the coefficient (vector) space you might get a significant test result whereas in another you might not. For categorical variables the approach should still work (locally) but for continuous I see where the issues arise-you'd need some sort of binning. From the linked issue it looks like (I didn't read beyond titles) there are some papers which you've found to give approaches for regression models. One thing I'd wonder about is criteria for identifiability with the logit and poisson simultaneous models. $\endgroup$ – Lucas Roberts Dec 12 '20 at 17:11

I think there are different ways to do this. One thing you can do is to compare a zero-inflated negative binomial/Poisson model with its regular binomial/Poisson counter part without the zero-inflation component. It would look like this in R:

zinb <- read.csv("http://www.ats.ucla.edu/stat/data/fish.csv")
zinb <- within(zinb, {
  nofish <- factor(nofish)
  livebait <- factor(livebait)
  camper <- factor(camper)


## fit a negative binomial model
m1 <- glm.nb(count ~ child + camper, data = zinb)

## fit a zero-inflated negative binomial model
m1_zi <- zeroinfl(count ~ child + camper| persons,
               data = zinb, dist = "negbin", EM = TRUE)
## compare 2 models
vuong(m1, m1_zi)

For more information, see this ever useful tutorial.

  • 2
    $\begingroup$ but OP said they wanted to test for zero-inflation before fitting ZI model. i.e. something like a score test, e.g. Broek, Jan van den. 1995. “A Score Test for Zero Inflation in a Poisson Distribution.” Biometrics 51 (2): 738–43. doi:10.2307/2532959. Yang, Zhao, James W. Hardin, and Cheryl L. Addy. 2010. “Score Tests for Zero-Inflation in Overdispersed Count Data.” Communications in Statistics - Theory and Methods 39 (11): 2008–30. doi:10.1080/03610920902948228./03610920902948228#.VDVM81QUBAs $\endgroup$ – Ben Bolker Oct 8 '14 at 14:43
  • 3
    $\begingroup$ I'm not familiar with score-test type approaches, but examining zero-inflation beforehand should have merits. Still, the zero-inflation component sometimes works in conjunction with the regression model. For example, if you fit zeroinfl(count ~ 1|persons, ...), persons is nowhere near significant. So I believe zero-inflation should be investigated given the model under hypothesis. $\endgroup$ – Masato Nakazawa Oct 8 '14 at 14:55
  • 2
    $\begingroup$ Use of the Vuong test for this purpose may not be appropriate: sciencedirect.com/science/article/pii/S016517651400490X $\endgroup$ – Patrick B. Oct 15 '16 at 23:01
  • 1
    $\begingroup$ Sorry, I hit Enter too soon. I may be butchering the argument in the above link but I think one major problem is that the ZI model reduces to the non-ZI model at the edge of the parameter space (i.e., when the zero-inflation parameter is zero), which violates one of the assumptions of the Vuong test of non-nested models. $\endgroup$ – Patrick B. Oct 15 '16 at 23:10
  • $\begingroup$ I agree with @PatrickB, both models are nested. The score test mentioned is very slow to converge to its asymptotic distribution, if it is not too difficult to fit both model I would do that and use a likelihood-ratio test. Also remember a dataset is never zero-inflated $an$ $sich$, but always with respect to some proposed distribution, as Masato Nakazawa rightly says. $\endgroup$ – Knarpie Feb 20 '17 at 15:08

Consider some model $f(x)$. If we want to turn $f(x)$ into a zero-inflated model, then we define $g(x)$ to equal $f(x)$ with proportion $p$ and to equal $0$ with proportion $1-p$.

In this case, there are two processes at work here. One process generates only zeros and one process generates results from $f(x)$. My understanding is that a zero-inflated model is only appropriate when there is an alternate process that generates only zeros. For example, if you are attempting to estimate the number of widgets different stores sell, but some stores do not have widgets for sale, then it seems like two processes are at work here: one process that generates only zeros (those stores that cannot sell widgets because they do not ever stock widgets for sale) and another process that generates different values (those stores that do stock widgets and therefore can sell some).

Rather than having a "test" to determine whether the data are zero-inflated, I would suggest determining whether it is plausible that there are two processes at work - one being a zero-generating process at work and another process that generates non-zero numbers. If it seems reasonable given the context of your data, then use a zero-inflated model. If it doesn't seem reasonable given the context of your data, then a zero-inflated model is probably inappropriate even though it may appear to fit your data better.

(It might not be clear from what I've written above, but I want to articulate the fact that both processes can generate zeros. One process generates only zeros and the other process can generate different values which may be zero. For example, a store can stock widgets and happen to sell zero widgets. This is different from a store that does not stock widgets and therefore must sell zero widgets by default.)

  • 1
    $\begingroup$ This is a reasonable comment, but I don't think it really answers the question. I can think of several scenarios where zero-inflation is mechanistically plausible, but may or may not be necessary to model in practice. This is important because there are also costs associated to switching to a ZI model: fitting zero-inflated count models is usually much more computationally intensive than fitting a non-inflated count model, and hurdle methods are faster but the zero-inflation parameter can be hard to interpret. $\endgroup$ – Patrick B. Oct 15 '16 at 23:02

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.