zibinomial: Zero-Inflated Binomial Distribution Family Function

Description

Fits a zero-inflated binomial distribution by maximum likelihood estimation.

Usage

zibinomial(lphi="logit", lmu="logit", ephi=list(), emu=list(),
           iphi=NULL, zero=1, mv=FALSE)

Arguments

lphi, lmu

Link functions for the parameter $\phi$ and the usual binomial probability $\mu$ parameter. See Links for more choices.

ephi, emu

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 recyling is used.

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 the first. Setting zero=NULL enables both $\phi$ and $\mu$ to be modelled as a function

Logical. Currently it must be FALSE to mean the function does not handle multivariate responses. This is to remain compatible with the same argument in binomialff.

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. Half-stepping is not uncommon. If failure to converge occurs, make use of the argument iphi.

Details

This function uses Fisher scoring and is based on

P (Y = 0) = ϕ + (1 - ϕ) (1 - μ)^{N},

for $y=0$, and

P (Y = y) = (1 - ϕ) (\binom{N}{N y}) μ^{N y} (1 - μ)^{N (1 - y)} .

for $y=1/N,2/N,\ldots,1$. That is, the response is a sample proportion out of $N$ trials, and the argument size in rzibinom is $N$ here. The parameter $\phi$ satisfies $0 < \phi < 1$. The mean of $Y$ is $E(Y)=(1-\phi) \mu$ and these are returned as the fitted values. By default, the two linear/additive predictors are $(logit(\phi), logit(\mu))^T$.

Examples

Run this code

size = 10  # number of trials; N in the notation above
nn = 200
zibdata = data.frame(phi = logit(0, inv=TRUE), # 0.50
                     mubin = logit(-1, inv=TRUE), # Mean of usual binomial
                     sv = rep(size, len=nn))
zibdata = transform(zibdata, 
                    y = rzibinom(n=nn, size=sv, prob=mubin, phi=phi)/sv)
with(zibdata, table(y))
fit = vglm(y ~ 1, zibinomial, weight=sv, data=zibdata, trace=TRUE)
coef(fit, matrix=TRUE)
Coef(fit) # Useful for intercept-only models
fit@misc$p0  # Estimate of P(Y=0)
head(fitted(fit))
with(zibdata, mean(y)) # Compare this with fitted(fit)
summary(fit)

Run the code above in your browser using DataLab

Last chance! 50% off unlimited learning