VGAM (version 0.7-10)

posnegbinomial: Positive Negative Binomial Distribution Family Function

Description

Maximum likelihood estimation of the two parameters of a positive negative binomial distribution.

Usage

posnegbinomial(lmunb = "loge", lk = "loge", emunb =list(), ek = list(),
               ik = NULL, zero = -2, cutoff = 0.995, shrinkage.init=0.95,
               method.init=1)

Value

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

Warning

The Poisson model corresponds to k equalling infinity. If the data is Poisson or close to Poisson, numerical problems may occur. Possibly a loglog link could be added in the future to try help handle this problem.

Details

The positive negative binomial distribution is an ordinary negative binomial distribution but with the probability of a zero response being zero. The other probabilities are scaled to sum to unity.

This family function is based on negbinomial and most details can be found there. To avoid confusion, the parameter munb here corresponds to the mean of an ordinary negative binomial distribution negbinomial. The mean of posnegbinomial is $$\mu_{nb} / (1-p(0))$$ where $p(0) = (k/(k + \mu_{nb}))^k$ is the probability an ordinary negative binomial distribution has a zero value.

The parameters munb and k are not independent in the positive negative binomial distribution, whereas they are in the ordinary negative binomial distribution.

This function handles multivariate responses, so that a matrix can be used as the response. The number of columns is the number of species, say, and setting zero=-2 means that all species have a k equalling a (different) intercept only.

References

Barry, S. C. and Welsh, A. H. (2002) Generalized additive modelling and zero inflated count data. Ecological Modelling, 157, 179--188.

Williamson, E. and Bretherton, M. H. (1964) Tables of the logarithmic series distribution. Annals of Mathematical Statistics, 35, 284--297.

See Also

rposnegbin, pospoisson, negbinomial, zanegbinomial, rnbinom.

Examples

pndat = data.frame(x = runif(nn <- 2000))
pndat = transform(pndat, y1 = rposnegbin(nn, munb=exp(0+2*x), size=exp(1)),
                         y2 = rposnegbin(nn, munb=exp(1+2*x), size=exp(3)))
fit = vglm(cbind(y1,y2) ~ x, posnegbinomial, pndat, trace=TRUE)
coef(fit, matrix=TRUE)
dim(fit@y)


# Another artificial data example
pndat2 = data.frame(munb = exp(2), size = exp(3)); nn = 1000
pndat2 = transform(pndat2, y = rposnegbin(nn, munb=munb, size=size))
with(pndat2, table(y))
fit = vglm(y ~ 1, posnegbinomial, pndat2, trace=TRUE)
coef(fit, matrix=TRUE)
with(pndat2, mean(y))    # Sample mean
head(with(pndat2, munb/(1-(size/(size+munb))^size)), 1) # Population mean
head(fitted(fit), 3)
head(predict(fit), 3)


# Example: Corbet (1943) butterfly Malaya data
corbet = data.frame(nindiv = 1:24,
                    ofreq = c(118, 74, 44, 24, 29, 22, 20, 19, 20, 15, 12,
                              14, 6, 12, 6, 9, 9, 6, 10, 10, 11, 5, 3, 3))
fit = vglm(nindiv ~ 1, posnegbinomial, weights=ofreq, data=corbet)
coef(fit, matrix=TRUE)
Coef(fit)
(khat = Coef(fit)['k'])
pdf2 = dposnegbin(x=with(corbet, nindiv), mu=fitted(fit), size=khat)
print( with(corbet, cbind(nindiv, ofreq, fitted=pdf2*sum(ofreq))), dig=1)