VGAM (version 0.8-1)

zipoisson: Zero-Inflated Poisson Distribution Family Function

Description

Fits a zero-inflated Poisson distribution by full maximum likelihood estimation.

Usage

zipoisson(lphi="logit", llambda = "loge", ephi=list(), elambda =list(),
          iphi = NULL, method.init=1, shrinkage.init=0.8, zero = NULL)

Arguments

lphi

Link function for the parameter $\phi$. See Links for more choices.

llambda

Link function for the usual $\lambda$ parameter. See Links for more choices.

ephi, elambda

List. Extra argument for the respective links. See earg in Links for general information.

iphi

Optional initial values for $\phi$, whose values must lie between 0 and 1. The default is to compute an initial value internally. If a vector then recycling is used.

method.init

An integer with value 1 or 2 which specifies the initialization method for $\lambda$. If failure to converge occurs try another value and/or else specify a value for shrinkage.init and/or else specify a value

shrinkage.init

How much shrinkage is used when initializing $\lambda$. The value must be between 0 and 1 inclusive, and a value of 0 means the individual response values are used, and a value of 1 means the median or mean is used. This argument is used in conju

zero

An integer specifying which linear/additive predictor is modelled as intercepts only. If given, the value must be either 1 or 2, and the default is none of them. Setting zero=1 makes $\phi$ a single parameter.

Value

• An object of class "vglmff" (see vglmff-class). The object is used by modelling functions such as vglm and vgam.

Warning

Numerical problems can occur, e.g., when the probability of zero is actually less than, not more than, the nominal probability of zero. For example, in the Angers and Biswas (2003) data below, replacing 182 by 1 results in nonconvergence. Half-stepping is not uncommon. If failure to converge occurs, try using combinations of method.init, shrinkage.init, iphi, and/or zero=1 if there are explanatory variables.

Details

The model is a mixture of a Poisson distribution and the value 0; it has value 0 with probability $\phi$ else is Poisson($\lambda$) distributed. The model can be written $$P(Y=0) = \phi + (1-\phi) \exp(-\lambda),$$ and for $y=1,2,\ldots$, $$P(Y=y) = (1-\phi) \exp(-\lambda) \lambda^y / y!.$$ Here, the parameter $\phi$ satisfies $0 < \phi < 1$. The mean of $Y$ is $(1-\phi) \lambda$ and these are returned as the fitted values. By default, the two linear/additive predictors are $(logit(\phi), \log(\lambda))^T$. This function implements Fisher scoring.

References

Thas, O. and Rayner, J. C. W. (2005) Smooth tests for the zero-inflated Poisson distribution. Biometrics, 61, 808--815.

Data: Angers, J-F. and Biswas, A. (2003) A Bayesian analysis of zero-inflated generalized Poisson model. Computational Statistics & Data Analysis, 42, 37--46.

Cameron, A. C. and Trivedi, P. K. (1998) Regression Analysis of Count Data. Cambridge University Press: Cambridge.

See Also

zapoisson, Zipois, yip88, zipebcom, rpois.

Examples

zipdat = data.frame(x = runif(nn <- 2000))
zipdat = transform(zipdat, phi    = logit(-0.5 + 1*x, inverse=TRUE),
                           lambda =  loge( 0.5 + 2*x, inverse=TRUE))
zipdat = transform(zipdat, y = rzipois(nn, lambda, phi))
with(zipdat, table(y))
fit = vglm(y ~ x, zipoisson, zipdat, trace=TRUE)
coef(fit, matrix=TRUE)  # These should agree with the above values


# Another example: McKendrick (1926). Data from 223 Indian village households
cholera = data.frame(ncases = 0:4, # Number of cholera cases,
                     wfreq = c(168, 32, 16, 6, 1)) # Frequencies
fit = vglm(ncases ~ 1, zipoisson, wei=wfreq, data=cholera, trace=TRUE)
coef(fit, matrix=TRUE)
with(cholera, cbind(actual=wfreq, fitted= round(
      dzipois(ncases, lambda=Coef(fit)[2], phi=Coef(fit)[1]) * sum(wfreq), dig=2)))


# Another example: data from Angers and Biswas (2003)
abdat = data.frame(y = 0:7, w = c(182, 41, 12, 2, 2, 0, 0, 1))
abdat = subset(abdat, w>0)
fit = vglm(y ~ 1, zipoisson(lphi=probit, iphi=0.3), abdat, wei=w, trace=TRUE)
fit@misc$prob0  # Estimate of P(Y=0)
coef(fit, matrix=TRUE)
Coef(fit)  # Estimate of phi and lambda
fitted(fit)
with(abdat, weighted.mean(y, w)) # Compare this with fitted(fit)
summary(fit)