R: Negative-Binomial Model for Counts

negbin {aod}

R Documentation

Negative-Binomial Model for Counts

Description

The function fits a negative-binomial log linear model accounting for overdispersion in counts y.

Usage

  negbin(formula, random, data, phi.ini = NULL, warnings = FALSE, 
         na.action = na.omit, fixpar = list(),
         hessian = TRUE, control = list(maxit = 2000), ...)

Arguments

`formula`	A formula for the fixed effects. The left-hand side of the formula must be the counts `y` i.e., positive integers (`y >= 0`). The right-hand side can involve an offset term.
`random`	A right-hand formula for the overdispersion parameter(s) φ.
`data`	A data frame containing the response (`y`) and explanatory variable(s).
`phi.ini`	Initial values for the overdispersion parameter(s) φ. Default to 0.1.
`warnings`	Logical to control printing of warnings occurring during log-likelihood maximization. Default to FALSE (no printing).
`na.action`	A function name. Indicates which action should be taken in the case of missing value(s).
`fixpar`	A list with 2 components (scalars or vectors) of the same size, indicating which parameters are fixed (i.e., not optimized) in the global parameter vector (b, φ) and the corresponding fixed values. For example, `fixpar = list(c(4, 5), c(0, 0))` means that 4th and 5th parameters of the model are set to 0.
`hessian`	A logical. When set to `FALSE`, the hessian and the variances-covariances matrices of the parameters are not computed.
`control`	A list to control the optimization parameters. See `optim`. By default, set the maximum number of iterations to 2000.
`...`	Further arguments passed to `optim`.

Details

For a given count y, the model is:

y | λ ~ Poisson(λ)

with λ following a Gamma distribution Gamma(r, θ).
If G denote the gamma function, then:

P(λ) = r^{-θ} * λ^{θ - 1} * exp(-λ / r) / G(θ)

E[λ] = θ * r

Var[λ] = θ * r^2

The marginal negative-binomial distribution is:

P(y) = G(y + θ) * (1 / (1 + r))^θ * (r / (1 + r))^y / (y! * G(θ))

The function uses the parameterization μ = θ * r = exp(X b) = exp(η) and φ = 1 / θ, where X is a design-matrix, b is a vector of fixed effects, η = X b is the linear predictor and φ the overdispersion parameter.
The marginal mean and variance are:

E[y] = μ

Var[y] = μ + φ * μ^2

The parameters b and φ are estimated by maximizing the log-likelihood of the marginal model (using the function optim()). Several explanatory variables are allowed in b. Only one is allowed in φ.
An offset can be specified in the formula argument to model rates y/T. The offset and the marginal mean are log(T) and μ = exp(log(T) + η), respectively.

Value

An object of formal class “glimML”: see glimML-class for details.

Author(s)

Matthieu Lesnoff matthieu.lesnoff@cirad.fr, Renaud Lancelot renaud.lancelot@cirad.fr

References

Lawless, J.F., 1987. Negative binomial and mixed Poisson regression. The Canadian Journal of Statistics, 15(3): 209-225.

Examples

  # without offset
  data(salmonella)
  negbin(y ~ log(dose + 10) + dose, ~ 1, salmonella)
  library(MASS) # function glm.nb in MASS
  fm.nb <- glm.nb(y ~ log(dose + 10) + dose,
                  link = log, data = salmonella)
  coef(fm.nb)
  1 / fm.nb$theta # theta = 1 / phi
  c(logLik(fm.nb), AIC(fm.nb))
  # with offset
  data(dja)
  negbin(y ~ group + offset(log(trisk)), ~ group, dja)
  # phi fixed to zero in group TREAT
  negbin(y ~ group + offset(log(trisk)), ~ group, dja,
    fixpar = list(4, 0))
  # glim without overdispersion
  summary(glm(y ~ group + offset(log(trisk)),
    family = poisson, data = dja))
  # phi fixed to zero in both groups
  negbin(y ~ group + offset(log(trisk)), ~ group, dja,
    fixpar = list(c(3, 4), c(0, 0)))

[Package aod version 1.3 Index]