negbin {aod} | R Documentation |
The function fits a negative-binomial log linear model accounting for overdispersion in counts y.
negbin(formula, random, data, phi.ini = NULL, warnings = FALSE, na.action = na.omit, fixpar = list(), hessian = TRUE, control = list(maxit = 2000), ...)
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) \phi. |
data |
A data frame containing the response (y ) and explanatory variable(s). |
phi.ini |
Initial values for the overdispersion parameter(s) \phi. 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, \phi) 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 . |
For a given count y, the model is:
y | \lambda ~ Poisson(\lambda)
with \lambda following a Gamma distribution Gamma(r, \theta).
If G denote the gamma function, then:
P(\lambda) = r^{-\theta} * \lambda^{\theta - 1} * exp(-\lambda / r) / G(\theta)
E[\lambda] = \theta * r
Var[\lambda] = \theta * r^2
The marginal negative-binomial distribution is:
P(y) = G(y + \theta) * (1 / (1 + r))^\theta * (r / (1 + r))^y / (y! * G(\theta))
The function uses the parameterization \mu = \theta * r = exp(X b) = exp(\eta) and \phi = 1 / \theta,
where X is a design-matrix, b is a vector of fixed effects, \eta = X b is the linear predictor and
\phi the overdispersion parameter.
The marginal mean and variance are:
E[y] = \mu
Var[y] = \mu + \phi * \mu^2
The parameters b and \phi 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 \phi.
An offset can be specified in the formula
argument to model rates y/T. The offset and the marginal mean
are log(T) and \mu = exp(log(T) + \eta), respectively.
An object of formal class “glimML”: see glimML-class
for details.
Matthieu Lesnoff matthieu.lesnoff@cirad.fr, Renaud Lancelot renaud.lancelot@cirad.fr
Lawless, J.F., 1987. Negative binomial and mixed Poisson regression. The Canadian Journal of Statistics, 15(3): 209-225.
glimML-class
, glm
and optim
,
glm.nb
in the recommended package MASS,
gnlr
in package gnlm available at www.luc.ac.be/~jlindsey/rcode.html.
# 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)))