# Parametric modelling of variance of count data

I am looking to model some data, but I am not sure what type of model I can use. I have count data, and I want a model that will give parametric estimates of both the mean and the variance of the data. That is, I have various predictive factors and I want to determine if any of them influence the variance (not just the group mean).

I know that Poisson regression will not work because the variance is equal to the mean; this assumption is not valid in my case, so I know there is overdispersion. However, a negative binomial model only generates a single overdispersion parameter, not one that is a function of the predictors in the model. What model can do this?

Additionally, a reference to a book or paper which discusses the model and/or an R package which implements the model would be appreciated.

• How do you know there is overdispersion without first doing the Poisson regression? After all, comparing the variance of the raw (response) values to their mean isn't relevant: what matters is the goodness of fit of the Poisson model (this is the analog of evaluating the distribution of residuals in a linear model compared to evaluating the distribution of the response variable). Another way to put this is that the link between the independent variables and the response can create the appearance of overdispersion even in a beautifully accurate Poisson model.
– whuber
Oct 19, 2011 at 21:31
• @whuber That's a fair point. For a single categorical predictor looking at the variance and mean of the sub-groups would be sufficient to detect overdispersion, but for a multivariate Poisson regression, it is not. For the sake of argument, let's assume both a Poisson and negative binomal regression have been done and the negative binomial shows a better fit via anova model comparison. That should indicate overdispersion. Given that, how could the variance/overdispersion be modeled parametrically rather than as a constant? Oct 19, 2011 at 21:46
• I think there's a chapter in McCullagh and Nelder, Generalized linear models, 2nd edition, that covers this (but my copy's at work)...there won't be a real likelihood, but you can use quasi-likelihood, and so that may be the title of the chapter. You apply iteratively reweighted least squares even though there's no probability model that corresponds.
– Karl
Oct 19, 2011 at 23:30
• Chapter 10 of McCullagh and Nelder discusses Joint Modelling of Mean and Dispersion, i.e. parameterizing both the mean and the variance. Extended quasi-likelihood is the main tool, but in some situations there can be concerns about that method Dec 2, 2011 at 18:33

You can model the negative binomial dispersion parameter itself as a function of variables and parameters using the gamlss package in R. I provide an excerpt from an introduction to it:

Why should I use GAMLSS

If your response variable is count (discrete) data it is very likely that the Poisson distribution will not fit well. GAMLSS provides a variety of discrete distributions (including the negative binomial) that you can try out. The dispersion parameter can be also modelled as a function of explanatory variables.

The www.gamlss.org website has documentation and links to several papers on the approaches used in the package.

• Both replies are helpful and provide good references. I'm awarding the bounty to this one because (a) it preceded the other by four minutes and (b) the gamlss solution is new to me (I'm familiar with nbreg). But hats off to @timbp for providing a good reply; I hope you'll continue to contribute to our site.
– whuber
Dec 4, 2011 at 21:56
• @whuber, I was also torn as to which to accept as "the" answer because both were very helpful. I went with this one because it included an R package reference which I can use; the book reference in the other answer has been good reading and should not be discounted. Thank you for offering the bounty which prompted these two good answers. Dec 6, 2011 at 17:52

Stata provides the -gnbreg- command, which allows you to model the dispersion parameter. You can view Stata help for the command at http://www.stata.com/help.cgi?nbreg

Stata calls this the generalised negative binomial model. Joseph Hilbe discusses it in his book "Negative Binomial Regression", section 10.4, as "NB-H: Heterogeneous negative binomial regression".