BASICS OF COUNT MODEL FIT STATISTICS

Most statisticians consider overdispersion the key problem when considering count model fit.

That is, an analyst typically attempts to evaluate whether a count model is extradispersed – which usually means overdispersed. If there is evidence of overdispersion in a Poisson model, the problem then is to determine what gives rise to it. If we can determine the cause, we can employ the appropriate model to use on the data.

The deviance goodness-of-fit (GOF) test is based on the view that the deviance is distributed as Chi2. The Chi2 distribution has two parameters – the mean and scale. For the deviance GOF, this is the deviance statistic and residual degrees of freedom, which is the number of observations less predictors in the model, including the intercept and interactions. If the resulting Chi2 p-value is less than 0.05, the model is considered well fit.

Employing the deviance GOF test on the modeled German health data yields the following model:

library(COUNT) 
## Loading required package: msme
## Loading required package: MASS
## Loading required package: lattice
## Loading required package: sandwich
data(rwm5yr) 
rwm1984 <- subset(rwm5yr, year==1984)
modx<-glm(docvis ~ outwork + female + age, family=poisson, data=rwm1984)
modx
## 
## Call:  glm(formula = docvis ~ outwork + female + age, family = poisson, 
##     data = rwm1984)
## 
## Coefficients:
## (Intercept)      outwork       female          age  
##    -0.15269      0.26508      0.28650      0.02267  
## 
## Degrees of Freedom: 3873 Total (i.e. Null);  3870 Residual
## Null Deviance:       25790 
## Residual Deviance: 24000     AIC: 31090
dev<-deviance(modx) 
df<-df.residual(modx)
p_value<-1-pchisq(dev,df)
cbind.data.frame(c("Deviance GOF", "D", "df", "p-value"),c(" ", round(dev,4), df, p_value))
##   c("Deviance GOF", "D", "df", "p-value") c(" ", round(dev, 4), df, p_value)
## 1                            Deviance GOF                                   
## 2                                       D                         24003.9403
## 3                                      df                               3870
## 4                                 p-value                                  0

The deviance is 24190.3581 and residual degrees of freedom is 3871.

These are both saved postestimation statistics that can be used for just this sort of purpose.

For the deviance GOF test, as well as the Pearson GOF tests that follow, a Chi2 p < 0.05 indicates that the model is considered well fit. More accurately, statistically speaking,

with a p < 0.05, the deviance GOF test indicates that we can reject the hypothesis that the model is not well fitted.

We also assume that the model has not been misspecified and that the model may appear to fit well by chance alone, as is the case for all such tests. Also, since it was first proposed, statisticians have discovered that many models appearing to be well fitted on the basis of the deviance test in fact poorly fit the data. Recall that a p-value of 0.05 accepts that we will be mistaken 1 out of 20 times on average.

If the value of D is very large, then we can generally be safe in rejecting the goodness of the model fit:

Deviance is in effect a measure of the distance between the most full or complete (saturated) model we can fit and the proposed model we are testing for fit.

The smaller the distance, or deviance, between them, the better the fit. The test statistic evaluates whether the value of the deviance, for a specific size of model, is close enough to that of the saturated model that it cannot be rejected as not fitting.

Some statisticians have used the deviance dispersion as the basis for scaling standard errors. However, when scaling standard errors, the Pearson dispersion statistic better captures the excess variability in the data, adjusting model standard errors in such a manner as to reflect what the standard errors would be if the excess variability were not present in the data.

Poisson models with scaled standard errors are called quasipoisson:

A Pearson dispersion in excess of 1.0 indicates likely Poisson model overdispersion. Whether the overdispersion is significant depends on (1) the value of the dispersion statistic, (2) the number of observations in the model, and (3)the structure of the data; for example, if the data are highly unbalanced.

Function to Calculate Pearson Chi2 and Dispersion Statistics

A code exists in R to calculate the deviance and Pearson Chi2 statistics and associated p-values for the associated GOF tests.

Pearson goodness-of-fit

pr <- sum(residuals(modx, type="pearson")^2) # get Pearson Chi2
pr
## [1] 44628.56
pchisq(pr, modx$df.residual, lower=F) # calc p-value
## [1] 0

Deviance goodness-of-fit

dv<-modx$deviance
dv
## [1] 24003.94
pchisq(dv, modx$df.residual, lower= F) # calc p-vl
## [1] 0

CALCULATING THE DISPERSION STATISTIC

P_disp <- function(x) {
  pr <- sum(residuals(x, type="pearson")^2)
  dispersion <- pr/x$df.residual
  cat("\n Pearson Chi2 = ", pr , "\n Dispersion = ", dispersion, "\n")
}

P_disp.r calculates the Pearson Chi2 and associated dispersion statistics.

P_disp(modx)
## 
##  Pearson Chi2 =  44628.56 
##  Dispersion =  11.53193

Not all overdispersion is real; apparent overdispersion can sometimes be identified and the model amended to eliminate it.

Overdispersion: What, Why, and How

TESTING OVERDISPERSION

The first set of tests to give a model when there is evidence of overdispersion are those listed previously.

These are possible remedies for apparent overdispersion, not specific tests of overdispersion. If equidispersion can result in a Poisson model by simply adding an interaction term, or a predictor that was previously neglected, then we do not need to test themodel for future overdispersion.

The important thing to remember, though, is that for the Poisson model, almost all violations of distributional assumptions, or other situations that result in a poorly fitted model, also result in overdispersion, and in some cases underdispersion:

Score Test

A score test to evaluate whether the amount of overdispersion in a Poisson model is sufficient to violate the basic assumptions of the model may be defined as \[z=\dfrac{(y-\mu)^2-y}{\mu\sqrt{2}}\]

library(COUNT) 
data(rwm5yr) 
rwm1984 <- subset(rwm5yr, year==1984)
summary(poi <- glm(docvis ~ outwork + age, family=poisson, data=rwm1984))
## 
## Call:
## glm(formula = docvis ~ outwork + age, family = poisson, data = rwm1984)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.0335170  0.0391815  -0.855    0.392    
## outwork      0.4079314  0.0188447  21.647   <2e-16 ***
## age          0.0220842  0.0008377  26.362   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 25791  on 3873  degrees of freedom
## Residual deviance: 24190  on 3871  degrees of freedom
## AIC: 31279
## 
## Number of Fisher Scoring iterations: 6
mu <-predict(poi, type="response")
z <- ((rwm1984$docvis - mu)^2 - rwm1984$docvis)/ (mu * sqrt(2))
summary(zscore <- lm(z ~ 1))
## 
## Call:
## lm(formula = z ~ 1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
##   -8.11   -7.40   -5.92   -5.10 2874.56 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    7.308      1.033   7.074 1.78e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 64.3 on 3873 degrees of freedom

The z score test is 7.3, with a t-probability of < 0.0005.

z tests the hypothesis that the Poisson model is overdispersed; it evaluates whether the data are Poisson or negative binomial.

This example indicates that the hypothesis of no overdispersion is rejected (i.e., that it is likely that real overdispersion exists in the data).

The test is based on two assumptions: * The data set on which the test is used is large. * z is t-distributed.

The preceding z-test result indicates that there is overdispersion in the data. Other z tests have been constructed, and you will find a variety of alternatives in the literature. The most popular alternative is perhaps an auxiliary regression test.

Lagrange Multiplier Test

The Lagrange multiplier test is a Chi2 test. It can be defined as \[\chi^2=\dfrac{\displaystyle\left(\displaystyle\sum_{i=1}^{n}\mu_i^2-n\bar{y}_i\right)^2}{2\displaystyle\sum_{i=1}^{n}\mu_i^2}\] with one degree of freedom.

obs <- nrow(rwm1984) 
mmu <- mean(mu) 
nybar <- obs*mmu 
musq <- mu*mu
mu2 <- mean(musq)*obs
chival <- (mu2 - nybar)^2/(2*mu2) 
chival
## [1] 11526.78
pchisq(chival,1,lower.tail = FALSE)
## [1] 0

With one degree of freedom, the test appears to be significant – the hypothesis of no overdispersion is again rejected.