本篇介绍一种常见的广义线性模型:泊松回归。泊松分布是离散型分布,它的概率分布函数如下:
写成指数族分布的形式如下:
对照指数族分布的通式:
可得,
广义线性模型假设 与解释变量存在线性关系,即
又因为泊松分布的数学期望就是 ,因此泊松回归可写成如下形式:
泊松回归是计数模型,因变量必须为自然数。一般针对具有稳定发生频数/率的事件进行建模。在glm
函数中,泊松回归的family
参数设置为poisson(link = "log")
或poisson()
。
示例
示例数据是来自dlnm
工具包的chicagoNMMAPS
数据集,它记录的是芝加哥逐日的死亡人口数和温度、环境污染物浓度等信息。
library(dlnm)
names(chicagoNMMAPS)
## [1] "date" "time" "year" "month" "doy" "dow" "death" "cvd" "resp"
## [10] "temp" "dptp" "rhum" "pm10" "o3"
在一段时间内,自然状态下地区的死亡人口可以看作是等频数发生的,因此可以使用泊松回归进行建模。
model <- glm(death ~ 1, family = poisson(),
data = chicagoNMMAPS)
summary(model)
##
## Call:
## glm(formula = death ~ 1, family = poisson(), data = chicagoNMMAPS)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.6735 -0.9850 -0.1323 0.7891 21.2791
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.748568 0.001302 3648 <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: 9873.8 on 5113 degrees of freedom
## Residual deviance: 9873.8 on 5113 degrees of freedom
## AIC: 43525
##
## Number of Fisher Scoring iterations: 4
在模型中添加温度变量temp
作为解释变量:
model.2 <- glm(death ~ temp, family = poisson(),
data = chicagoNMMAPS)
summary(model.2)
##
## Call:
## glm(formula = death ~ temp, family = poisson(), data = chicagoNMMAPS)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.3097 -0.8563 -0.0719 0.7530 22.5216
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.7927697 0.0017346 2762.97 <2e-16 ***
## temp -0.0044903 0.0001197 -37.51 <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: 9873.8 on 5113 degrees of freedom
## Residual deviance: 8471.8 on 5112 degrees of freedom
## AIC: 42125
##
## Number of Fisher Scoring iterations: 4
在模型中添加温度的时间滞后项作为解释变量:
library(dplyr)
model.3 <- glm(death ~ temp + lag(temp, 1) + lag(temp, 2) +
lag(temp, 3) + lag(temp, 4) + lag(temp, 5),
family = poisson(),
data = chicagoNMMAPS)
summary(model.3)
##
## Call:
## glm(formula = death ~ temp + lag(temp, 1) + lag(temp, 2) + lag(temp,
## 3) + lag(temp, 4) + lag(temp, 5), family = poisson(), data = chicagoNMMAPS)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.8876 -0.8192 -0.0514 0.7029 22.6409
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.8028333 0.0017758 2704.577 < 2e-16 ***
## temp 0.0019664 0.0003720 5.286 1.25e-07 ***
## lag(temp, 1) -0.0011023 0.0005198 -2.121 0.034 *
## lag(temp, 2) -0.0021813 0.0005301 -4.115 3.87e-05 ***
## lag(temp, 3) -0.0005632 0.0005295 -1.064 0.287
## lag(temp, 4) -0.0008187 0.0005193 -1.577 0.115
## lag(temp, 5) -0.0028462 0.0003720 -7.651 1.99e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 9856.6 on 5108 degrees of freedom
## Residual deviance: 7768.2 on 5102 degrees of freedom
## (5 observations deleted due to missingness)
## AIC: 41398
##
## Number of Fisher Scoring iterations: 4
偏移项
在很多情况下,保持恒定的并非计数变量(如示例中的死亡人口),而是比率变量。在上述示例中,死亡人口受到总人口或老年人口(如60岁以上人口)数量的影响:
泊松回归要求因变量必须为非负整数,比率变量显然不能满足这一条件。
利用对数的运算法则,有如下变形:
在上式中,log(pop)
作为系数恒为1的解释变量,即偏移项。因此可以设置如下形式的模型:
model.4 <- glm(death ~ temp, offset = log(pop),
family = poisson(),
data = chicagoNMMAPS)
chicagoNMMAPS
数据集中无pop
变量,上式只作为演示,不能运行。
准泊松回归
泊松分布的均值和方差相等,均为 :
## (Dispersion parameter for poisson family taken to be 1)
但很多情况下该条件并满足:
deviance(model.2)/df.residual(model.2)
## [1] 1.657238
这种情况可以使用准泊松分布族quasipoisson()
:
model.22 <- glm(death ~ temp, family = quasipoisson(),
data = chicagoNMMAPS)
summary(model.22)
结果对比:
coef(summary(model.2))
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.792769652 0.0017346450 2762.9687 0.000000e+00
## temp -0.004490284 0.0001196991 -37.5131 5.633927e-308
coef(summary(model.22))
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.792769652 0.0023079821 2076.60605 0.000000e+00
## temp -0.004490284 0.0001592622 -28.19428 1.102289e-162
泊松回归和准泊松回归的系数估计值是一样的,但后者的标准误(Std. Error)较大。
相关阅读:
往期推荐阅读:
