[R] Negative Binomial Model

[Achim Zeileis]

On Wed, 14 May 2008, Mike Ryckman wrote:

> Hello,
>
> Thank you for your replies. I think that information largely focuses on
> robust standard errors... maybe I wasn't clear.

No, you were not and, as pointed out previously (and described in the 
posting guide), a reproducible example would clearly help.

> By default Stata parameterizes the dispersion. From the stata help files 
> it is calculated as 1 + alpha*exp(xb + offset).

If I interpret this correctly, this is would be a vector of length n, not 
a scaler dispersion.

> As I understand it, this process should happen outside of a decision to use
> robust errors. I think the normal glm.nb function in MASS uses some kind of
> constant in estimating the error process in a manner more analogous to how
> they would be estimated in any normal glm model.

There are different parametrizations of the negative binomial model but 
afaik both have dispersion = 1 (in the GLM sense) but have a variance 
function that can account for over-dispersion (compared to the Poisson 
variance function).

However, if you want to use a different dispersion, you can plug it into 
the summary method
   summary(fm, dispersion = ...)

Z


> I think the problem is that stata is using a different method to calculate
> the standard errors - and I'm not sure what that method is or how to do 
> it in R.
> The gamlss package produces the closest results - but even those are still
> not
> exactly the same.
>
> Mike
>
> -----Original Message-----
> From: Ted Harding [mailto:Ted.Harding at manchester.ac.uk]
> Sent: Wednesday, May 14, 2008 2:38 AM
> To: r-help at r-project.org
> Cc: Mike Ryckman
> Subject: Re: [R] Negative Binomial Model
>
> On 14-May-08 09:11:31, Achim Zeileis wrote:
>> On Wed, 14 May 2008, Mike Ryckman wrote:
>>
>>> Hello,
>>>
>>> I am trying to run a negative binomial regression model in R and
>>> can't get the standard errors to match the output I get from the
>>> Stata nbreg command. I've tried a few different options but
>>> haven't had much luck. The closest I've found is:
>>>
>>> gamlss(formula, family = NBI, sigma.formula = ~ 1,data=dataframe)
>>>
>>> ...But this is still a little off most of the time and pretty far off
>>> at
>>> other
>>> times (compared with the Stata output). The glm.nb from the MASS
>>> package
>>> produces the correct coefficients, but different (usually very
>>> different)
>>> standard errors.
>>>
>>> Could anybody explain this and point me in the right direction? I'd
>>> really
>>> appreciate it.
>>
>> Well, you don't really give us enough information to know (a
>> reproducible
>> R example and the desired standard errors from Stata would have been
>> helpful). My guess is that Stata uses some sort of "robust" standard
>> errors, i.e., sandwich standard errors. Try something like:
>>    library("MASS")
>>    library("sandwich")
>>    library("lmtest")
>>    fm <- glm.nb(...)
>>    coeftest(fm, vcov = sandwich)
>>
>> See also the following thread for some more discussion for count data
>> regression in R and Stata:
>>    https://stat.ethz.ch/pipermail/r-help/2008-May/161640.html
>
> Better links for following this thread are:
>
> Start:
>  https://stat.ethz.ch/pipermail/r-help/2008-May/161591.html
>
> Then click on "Next message:" for the second. Unfortunately,
> at that point the thread broke (Paul next responded to the list
> as a reply to an off-list message I sent him). So to continue,
> next take Achim's URL above:
>
>  https://stat.ethz.ch/pipermail/r-help/2008-May/161640.html
>
> and thereafter continue to click on "Next message:" until the
> thread runs out.
>
> That was a very helpful thread, for me!
> Best wishes,
> Ted.
>
> --------------------------------------------------------------------
> E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk>
> Fax-to-email: +44 (0)870 094 0861
> Date: 14-May-08                                       Time: 10:38:24
> ------------------------------ XFMail ------------------------------
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>