Solved – How to fit a zero inflated poisson model with only offset (without coefficients)

poisson distributionrregressionzero inflation

I have already got a poisson estimated lambda, and actual result y, and I would like to see if the model is good.

To start with, I check if the dispersion is alright.

glm(y ~ 0, offset=log(lambda), data=data, family=quasipoisson)

And then I would like to see if there is zero inflation factor by using pscl package.

zeroinfl(y ~ 1, offset=log(lambda), data=data ,dist="poisson")

However, I am not able to do this:

zeroinfl(y ~ 0, offset=log(lambda), data=data ,dist="poisson")

By introducing an interception, it is likely to interfere with the zero estimator so it becomes unclear if there is any zero (and how big) inflation in the original model.

Is it possible to have 0 estimated parameters (besides zeros)?

Best Answer

Actually, your second line still works - it just overrides your input $\lambda$ with an estimated $\lambda$. Unless you have good reason to believe your input $\lambda$ is more accurate than the estimate that is done by zeroinfl, I'd just leave it.

On the other hand, if you're really firm about using the input $\lambda$, the way to check for a zero inflation factor is to compare the observed number of zeros with the expected number of zeros given $\lambda$. You can do a standard significance test for this, or just use the difference divided by sample size as an estimate of the zero inflation factor:

lambda <- 2
p_zero <- 0.2
x = rpois(1000, lambda) * rbinom(1000, 1, 1-p_zero)

expected_zero <- exp(-lambda)
observed_zero <- mean(x==0)
> expected_zero
[1] 0.1353353
> observed_zero
[1] 0.33

and the difference of 0.1946647 would be your estimate of the zero inflation factor given your estimate of $\lambda$.

The significance test in this case can be done using the usual Normal approximation to the distribution of proportions to compare the observed and expected frequencies of zeroes:

> zstat <- (observed_zero-expected_zero) / sqrt(expected_zero*(1-expected_zero)/length(x))
> zstat
[1] 17.99525
> pnorm(zstat)
[1] 1
Related Question