Dear all, I would like to ask, if there is a way to make the variance / dispersion parameter $\theta$ (referring to MASS, 4th edition, p. 206) in the function glm.nb dependent on the data, e.g. $1/ \theta = exp(x \beta)$ and to estimate the parameter vector $\beta$ additionally. If this is not possible with glm.nb, is there another function / package which might do that? Thank you very much for your answer in advance! Best, Martin
Negative Binomial Regression - glm.nb
3 messages · Martin Spindler, Brian Ripley, Ben Bolker
On 28/02/2013 07:27, Martin Spindler wrote:
Dear all, I would like to ask, if there is a way to make the variance / dispersion parameter $\theta$ (referring to MASS, 4th edition, p. 206) in the function glm.nb dependent on the data, e.g. $1/ \theta = exp(x \beta)$ and to estimate the parameter vector $\beta$ additionally.
That is no longer a glm, so no.
If this is not possible with glm.nb, is there another function / package which might do that?
You can maximize the likelihood directly. How to do that is exemplified in the optimization chapter of MASS.
Thank you very much for your answer in advance! Best, Martin
Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595
Martin Spindler <Martin.Spindler <at> gmx.de> writes:
Dear all, I would like to ask, if there is a way to make the variance / dispersion parameter $\theta$ (referring to MASS, 4th edition, p. 206) in the function glm.nb dependent on the data, e.g. $1/ \theta = exp(x \beta)$ and to estimate the parameter vector $\beta$ additionally. If this is not possible with glm.nb, is there another function / package which might do that?
As Brian Ripley says, that's outside the scope of glm.nb(),
and a later chapter in MASS tells you how to do it.
The mle2 function in the bbmle package offers a possibly
convenient shortcut. Something like
mle2(response~dnbinom(mu=exp(logmu),size=exp(logtheta)),
parameters=list(logmu~[formula for linear predictor of log(mu)],
logtheta~[formula for linear predictor of log(beta)]),
start=list(logmu=[starting value for logmu intercept],
logbeta=[starting value for logbeta intercept]),
data=...)
should work ...
Ben Bolker