7
$\begingroup$

I am modelling a zero-truncated process with a count model, and am trying to determine whether the data are overdispersed. The Poisson distribution has a variance equal to its mean, $$\newcommand{\Var}{\operatorname{Var}} \Var(y) = E(y) = \lambda $$ The negative binomial model relaxes this assumption by estimating an overdispersion parameter $\alpha$. There are two regimens for this method, NB1 $$\Var(y) = \lambda \alpha$$ and NB2 $$\Var(y) = \lambda (1 + \lambda/\alpha)$$ I have seen the second version more often. I have estimated three models with the same covariates but the three different link functions (Poisson, NB1, NB2).

Model   logLik   alpha  p(alpha)
Poisson -7942    1      NA
NB1     -6399    1.001  0
NB2     -7944    403.4  0

Clearly the overdispersion parameter is significant and I should therefore use one of the NB models. But how can the log-likelihood for NB2 be lower than for Poisson when there is an additional parameter in the model?

When I estimate all three models, I get a warning that my variance-covariance matrix may not be positive semi-definite. But I don't have highly collinear variables, and all of the relevant matrices have strictly positive eigenvalues, so I don't see why I'm getting this warning. Could this be related to the likelihood problem?

I am using the glmmadmb function from the R package of the same name.

$\endgroup$
1
  • 3
    $\begingroup$ I would recommend that you also try to estimate the variance-covariance matrix of the Poisson estimates using the Huber/White/Sandwich linearized estimator. That estimator of the variance-covariance matrix does not assume $E(y_i) = Var(y_i)$, nor does it even require that $Var(y_i)$ be constant across $i$. That would be less of apples-to-orangutangs comparison with the NB model. I can't tell if glmmadmb allows you to calculate the errors that way or not. $\endgroup$
    – dimitriy
    Commented Mar 29, 2013 at 20:45

1 Answer 1

5
$\begingroup$

You are right, the NB2 model should not have a worse fit than the poisson model, because the latter is a special case of the former, where $\alpha=0$. There is definitely a problem. It sounds to me like your numerical maximum likelihood maximizer had problems with your specification and therefore didn't really hit the global maximum.

Now I don't use R, but in Stata the NB2 model often has convergence problems. Stata has an option "difficult", which uses a different algorithm to jump from nonconcave regions of the likelihood function. Using this often solves the problem. Maybe R too has different options for your maximizer; perhaps, you could try different starting values and see whether that makes a difference.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.