[R] glm.poisson.disp versus glm.nb

[Hanke, Alex]

In response to my own query (see below),
The estimate theta from glm.nb is actually 1/phi or 1/alpha in some texts,
where phi  is the dispersion parameter for the negative binomial
distribution. However, the dispersion estimate from glm.poisson.disp does
not equal 1/theta because 1/theta is a maximum likelihood (ML) estimate
whereas the dispersion estimate from glm.poisson.disp is based on the method
of moments (MM). Correct me if I'm wrong!
The two estimates for dispersion are fairly different.
Paul and Banerjee (1998) observe that the test statistic discussed below
(TNBI) is more conservative when one uses the moment generated estimates for
phi, and that phi is underestimated by ML and overestimated by MM. However,
the overestimation by MM is slight compared to the underestimation by ML.
Does that then mean MM estimates of phi are to be preferred?

Alex

Dear list,
This is a question about overdispersion and the ML estimates of the
parameters returned by the glm.poisson.disp (L. Scrucca) and glm.nb
(Venables and Ripley) functions. Both appear to assume a negative binomial
distribution for the response variable.

Paul and Banerjee (1998) developed C(alpha) tests for "interaction and main
effects, in an unbalanced two-way layout of counts involving two fixed
factors, when data are Poisson distributed, and when data are extra
dispersed". In R I coded their C(alpha) test statistic (called TNBI) for
interaction for the case where the counts are extra-dispersed, as well as
their test for extra-dispersion (called T_a). Using the Quine data set
(Quine, 1975) the authors collapse the orginal 4x2x2x2 study design into a
2x4 table and obtained estimates of TNBI=10.36 and T_a=90.81.
Using the dispersion estimate from glm.poisson.disp and the estimates for mu
I got exactly the same results for TNBI and T_a. This made me happy. Now I
thought to try the ML estimates from glm.nb to see if the results would be
the same but I am having difficulty relating the dispersion phi from
glm.poisson.disp to theta estimated by glm.nb.
According to the R help for glm.poisson.disp " Var(y_i) =  mu_i(1+mu_i*phi)
". The help for glm.nb lead me to a book by V&R (1994) which indicates that
Var(y)=mu+mu^2/theta. From this I gathered that phi=1/theta but the
estimates do not relate to each other in this way unless one is in error. In
a document by L.P. Ammann he says a "negative binomial model can be
specified with mean mu and dispersion phi by taking theta=mu/(phi-1)". I had
a problem implementing this because in my mind mu is a vector whereas phi
and theta are scalars.



Alex Hanke
Department of Fisheries and Oceans
St. Andrews Biological Station
531 Brandy Cove Road
St. Andrews, NB
Canada
E5B 2L9