library(tidyverse)
library(broom)
library(knitr)
# add other packages as needed
Lecture 05 AE: Poisson Regression
Goodness-of-fit and overdispersion
<- read_csv("data/fHH1.csv") hh_data
Household vs. age, age^2, and location
<- glm(total ~ age + I(age^2) + location,
hh_age_loc data = hh_data, family = poisson)
tidy(hh_age_loc, conf.int = T) |>
kable(digits = 4)
term | estimate | std.error | statistic | p.value | conf.low | conf.high |
---|---|---|---|---|---|---|
(Intercept) | -0.3843 | 0.1821 | -2.1107 | 0.0348 | -0.7444 | -0.0306 |
age | 0.0704 | 0.0069 | 10.1900 | 0.0000 | 0.0569 | 0.0840 |
I(age^2) | -0.0007 | 0.0001 | -10.9437 | 0.0000 | -0.0008 | -0.0006 |
locationDavaoRegion | -0.0194 | 0.0538 | -0.3605 | 0.7185 | -0.1250 | 0.0859 |
locationIlocosRegion | 0.0610 | 0.0527 | 1.1580 | 0.2468 | -0.0423 | 0.1641 |
locationMetroManila | 0.0545 | 0.0472 | 1.1542 | 0.2484 | -0.0378 | 0.1473 |
locationVisayas | 0.1121 | 0.0417 | 2.6853 | 0.0072 | 0.0308 | 0.1945 |
Calculating residuals in R
You can get the residuals in two ways:
residuals
function withtype = "pearson"
ortype = "deviance"
.augment
function withtype.residuals = "pearson"
ortype.residuals = deviance
.
Calculate model deviance & Goodness-of-fit
Ex 1
Use the deviance residuals to calculate the model deviance.
# Calculate model deviance
Ex 2
Check your work above using the glance
function.
# Get deviance using glance function
Ex 3
Conduct a goodness-of-fit test
# Conduct goodness-of-fit test
Quasi-Poisson
<- glm(total ~ age + I(age^2) + location,
hh_age_loc_q data = hh_data, family = quasipoisson)
Ex 4
What do we expect the dispersion parameter to be if there is no overdispersion?
Ex 5
Use the Pearson residuals to calculate the dispersion parameter.
# Calculate dispersion parameter
Ex 6
View the dispersion parameter in the model output
#summary(hh_age_loc_q)
Negative binomial regression model
<- MASS::glm.nb(total ~ age + I(age^2) + location, data = hh_data)
hh_age_loc_nb tidy(hh_age_loc_nb) |>
kable(digits = 4)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | -0.3753 | 0.2076 | -1.8081 | 0.0706 |
age | 0.0699 | 0.0079 | 8.8981 | 0.0000 |
I(age^2) | -0.0007 | 0.0001 | -9.5756 | 0.0000 |
locationDavaoRegion | -0.0219 | 0.0625 | -0.3501 | 0.7262 |
locationIlocosRegion | 0.0577 | 0.0615 | 0.9391 | 0.3477 |
locationMetroManila | 0.0562 | 0.0551 | 1.0213 | 0.3071 |
locationVisayas | 0.1104 | 0.0487 | 2.2654 | 0.0235 |