Learn R Programming

VGAMdata (version 1.1-12)

oiposbinomial: One-Inflated Positive Binomial Distribution Family Function

Description

Fits a one-inflated positive binomial distribution by maximum likelihood estimation.

Usage

oiposbinomial(lpstr1 = "logitlink", lprob = "logitlink",
    type.fitted = c("mean", "prob", "pobs1", "pstr1", "onempstr1"),
    iprob = NULL, gpstr1 = ppoints(9), gprob  = ppoints(9),
    multiple.responses = FALSE, zero = NULL)

Value

An object of class "vglmff"

(see vglmff-class). The object is used by modelling functions such as vglm

and vgam.

Arguments

lpstr1, lprob

Link functions for the parameter \(\phi\) and the positive binomial probability \(\mu\) parameter. See Links for more choices. See CommonVGAMffArguments also. For the one-deflated model see below.

type.fitted

See CommonVGAMffArguments and fittedvlm.

iprob, gpstr1, gprob

For initial values; see CommonVGAMffArguments.

multiple.responses

Logical. See binomialff and posbinomial.

zero

See CommonVGAMffArguments for information.

Author

T. W. Yee

Details

These functions are based on $$P(Y=y) = \phi + (1-\phi) N \mu (1-\mu)^N / (1-(1-\mu)^N),$$ for \(y=1/N\), and $$P(Y=y) = (1-\phi) {N \choose Ny} \mu^{Ny} (1-\mu)^{N(1-y)} / (1-(1-\mu)^N).$$ for \(y=2/N,\ldots,1\). That is, the response is a sample proportion out of \(N\) trials, and the argument size in roiposbinom is \(N\) here. Ideally \(N > 2\) is needed. The parameter \(\phi\) is the probability of a structural one, and it satisfies \(0 < \phi < 1\) (usually). The mean of \(Y\) is \(E(Y)=\phi + (1-\phi) \mu / (1-(1-\mu)^N)\) and these are returned as the default fitted values. By default, the two linear/additive predictors for oiposbinomial() are \((logit(\phi), logit(\mu))^T\).

See Also

roiposbinom, posbinomial, binomialff, rbinom.

Examples

Run this code
size <- 10  # Number of trials; N in the notation above
nn <- 200
odata <- data.frame(pstr1  = logitlink( 0, inv = TRUE),  # 0.50
              mubin1 = logitlink(-1, inv = TRUE),  # Binomial mean
              svec   = rep(size, length = nn),
              x2     = runif(nn))
odata <- transform(odata,
                   mubin2 = logitlink(-1 + x2, inv = TRUE))
odata <- transform(odata,
              y1 = roiposbinom(nn, svec, pr = mubin1, pstr1 = pstr1),
              y2 = roiposbinom(nn, svec, pr = mubin2, pstr1 = pstr1))
with(odata, table(y1))
fit1 <- vglm(y1 / svec ~  1, oiposbinomial, data = odata,
             weights = svec, trace = TRUE, crit = "coef")
fit2 <- vglm(y2 / svec ~ x2, oiposbinomial, data = odata,
             weights = svec, trace = TRUE)

coef(fit1, matrix = TRUE)
Coef(fit1)  # Useful for intercept-only models
head(fitted(fit1, type = "pobs1"))  # Estimate of P(Y = 1)
head(fitted(fit1))
with(odata, mean(y1))  # Compare this with fitted(fit1)
summary(fit1)

Run the code above in your browser using DataLab