VGAM (version 0.8-7)

Posbinom: Positive-Binomial Distribution

Description

Density, distribution function, quantile function and random generation for the positive-binomial distribution.

Usage

dposbinom(x, size, prob, log = FALSE)
pposbinom(q, size, prob)
qposbinom(p, size, prob)
rposbinom(n, size, prob)

Arguments

x, q

vector of quantiles.

p

vector of probabilities.

n

number of observations. Fed into runif.

size

number of trials. It is the $N$ symbol in the formula given in posbinomial.

prob

probability of success on each trial.

log

See dbinom.

Value

• dposbinom gives the density, pposbinom gives the distribution function, qposbinom gives the quantile function, and rposbinom generates random deviates.

Details

The positive-binomial distribution is a binomial distribution but with the probability of a zero being zero. The other probabilities are scaled to add to unity. The mean therefore is $$\mu / (1-(1-\mu)^N)$$ where $\mu$ is the argument prob above. As $\mu$ increases, the positive-binomial and binomial distributions become more similar. Unlike similar functions for the binomial distribution, a zero value of prob is not permitted here.

See Also

posbinomial, zabinomial, zibinomial, rbinom.

Examples

prob = 0.2; size = 10
table(y <- rposbinom(n = 1000, size, prob))
mean(y)  # Sample mean
size * prob / (1-(1-prob)^size) # Population mean

(ii = dposbinom(0:size, size, prob))
cumsum(ii) - pposbinom(0:size, size, prob)  # Should be 0s
table(rposbinom(100, size, prob))

table(qposbinom(runif(1000), size, prob))
round(dposbinom(1:10, size, prob) * 1000) # Should be similar

barplot(rbind(dposbinom(x = 0:size, size, prob),
                           dbinom(x = 0:size, size, prob)),
        beside = TRUE, col = c("blue", "green"),
        main=paste("Positive-binomial(", size, ",", prob, ") (blue) vs",
        " Binomial(", size, ",", prob, ") (green)", sep = ""),
        names.arg = as.character(0:size), las = 1)

# Simulated data example
nn = 1000; sizeval1 = 10; sizeval2 = 20
pdat <- data.frame(x2 = seq(0, 1, length = nn))
pdat <- transform(pdat, prob1  = logit(-2 + 2 * x2, inverse = TRUE),
                        prob2  = logit(-1 + 1 * x2, inverse = TRUE),
                        sizev1 = rep(sizeval1, len = nn),
                        sizev2 = rep(sizeval2, len = nn))
pdat <- transform(pdat, y1 = rposbinom(nn, size = sizev1, prob = prob1),
                        y2 = rposbinom(nn, size = sizev2, prob = prob2))
with(pdat, table(y1))
with(pdat, table(y2))
# Multivariate response
fit2  = vglm(cbind(y1, y2) ~ x2, posbinomial(mv = TRUE),
             trace  = TRUE, pdat, weight = cbind(sizev1, sizev2))
coef(fit2, matrix = TRUE)