熱中症搬送人員数のオフセット項を65歳以上人口に変更してみる

　以前の記事では都道府県人口の対数をオフセット項として一般化線形回帰分析を行った．実際のところ，年代別の搬送人員としては65歳以上の高齢者が圧倒的に多い．そのため，東京など労働人口の多いところでは予測性能が悪化する可能性がある．今回はオフセット項の都道府県人口を3区分に分け，65歳以上人口の対数をオフセット項として採用してみたところ予測性能が改善したと思われたので記事とした．

データの準備

データのダウンロード

　例によってe-Statからダウンロードする．詳細はこの記事にある通りで，15歳未満人口，15歳以上65歳未満人口，65歳以上人口を追加する．

Power Queryでデータを整形する

　必要な列は年，都道府県コード，都道府県名，総人口，15歳未満人口，15歳以上64歳未満人口，65歳以上人口である．

SQL Serverにインポート

　ウィザードを使用して熱中症データベースにインポートする．

クエリ

　下記クエリを発行して抽出されたデータをEXCELにコピペする．

USE HeatStrokeDB;
GO
SELECT	Year(H.日付)	AS Year
,	MONTH(H.日付)	AS Month
,	H.日付	AS Date
,	P.都道府県コード	AS PrefCode
,	P.都道府県	AS Pref
,	T.日最高気温	AS Temp
,	V.日平均蒸気圧	AS Vapor
,	P.[15歳未満人口]	AS Pop14
,	P.[15～64歳人口]	AS Pop14_65
,	P.[65歳以上人口]	AS Pop65
,	H.[搬送人員（計）]	AS Num
FROM	dbo.T_HeatStroke	AS H
INNER	JOIN	dbo.T_MaxTemperature	AS T
ON	H.都道府県コード = T.都道府県コード
AND	H.日付 = T.年月日
INNER	JOIN	dbo.T_VaporPressure	AS V
ON	H.都道府県コード = V.都道府県コード
AND	H.日付 = V.年月日
INNER	JOIN	dbo.T_Population	AS P
ON	H.都道府県コード = P.都道府県コード
AND	YEAR(H.日付) = P.調査年
ORDER	BY	H.都道府県コード, H.日付

Rでの処理

　下記コードを実行する．

> y = HeatStroke5$Num
> x1 = HeatStroke5$Temp
> x2 = HeatStroke5$Vapor
> x5 = HeatStroke5$Month
> p1 = HeatStroke5$Pop14
> p2 = HeatStroke5$Pop14_65
> p3 = HeatStroke5$Pop65

> z1 = glm(y ~ x1 + x2 + x5, offset = log(p1), family = poisson(link = "log"))
> z2 = glm(y ~ x1 + x2 + x5, offset = log(p2), family = poisson(link = "log"))
> z3 = glm(y ~ x1 + x2 + x5, offset = log(p3), family = poisson(link = "log"))
> z4 = glm(y ~ x1 + x2 + x5, offset = log(p1 + p2 + p3), family = poisson(link = "log"))
> z5 = glm(y ~ x1 + x2 +x5, offset = log(p1 + p3), family = poisson(link = "log"))

　結果は以下の通り．5通りのモデルを作成したが，AICが最小となったのは65歳以上人口のみをオフセット項として投入したモデルz3であった．

> summary(z1)
Call:
glm(formula = y ~ x1 + x2 + x5, family = poisson(link = "log"), 
    offset = log(p1))
Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-11.4213   -1.2344   -0.4067    0.9100   21.9089  
Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -1.933e+01  1.701e-02 -1136.4   <2e-16 ***
x1           2.843e-01  5.079e-04   559.8   <2e-16 ***
x2           7.917e-02  3.964e-04   199.7   <2e-16 ***
x5          -2.999e-01  1.654e-03  -181.3   <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: 982005  on 70107  degrees of freedom
Residual deviance: 276188  on 70104  degrees of freedom
AIC: 459136
Number of Fisher Scoring iterations: 5

> summary(z2)
Call:
glm(formula = y ~ x1 + x2 + x5, family = poisson(link = "log"), 
    offset = log(p2))
Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-12.8365   -1.1833   -0.3580    0.9729   21.3463  
Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -2.091e+01  1.699e-02 -1230.8   <2e-16 ***
x1           2.849e-01  5.096e-04   559.1   <2e-16 ***
x2           8.062e-02  4.008e-04   201.1   <2e-16 ***
x5          -3.062e-01  1.665e-03  -183.9   <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: 997200  on 70107  degrees of freedom
Residual deviance: 278900  on 70104  degrees of freedom
AIC: 461848
Number of Fisher Scoring iterations: 5

> summary(z3)
Call:
glm(formula = y ~ x1 + x2 + x5, family = poisson(link = "log"), 
    offset = log(p3))
Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-10.3979   -1.2493   -0.4561    0.7605   19.5694  
Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -2.007e+01  1.679e-02 -1195.6   <2e-16 ***
x1           2.821e-01  5.052e-04   558.4   <2e-16 ***
x2           8.262e-02  3.984e-04   207.4   <2e-16 ***
x5          -2.948e-01  1.667e-03  -176.8   <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: 971373  on 70107  degrees of freedom
Residual deviance: 235828  on 70104  degrees of freedom
AIC: 418776
Number of Fisher Scoring iterations: 5

> summary(z4)
Call:
glm(formula = y ~ x1 + x2 + x5, family = poisson(link = "log"), 
    offset = log(p1 + p2 + p3))
Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-12.0744   -1.2083   -0.3938    0.9011   20.9182  
Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -2.140e+01  1.694e-02 -1263.6   <2e-16 ***
x1           2.841e-01  5.083e-04   559.0   <2e-16 ***
x2           8.096e-02  3.996e-04   202.6   <2e-16 ***
x5          -3.024e-01  1.664e-03  -181.7   <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: 984470  on 70107  degrees of freedom
Residual deviance: 263358  on 70104  degrees of freedom
AIC: 446306
Number of Fisher Scoring iterations: 5

> summary(z5)

Call:
glm(formula = y ~ x1 + x2 + x5, family = poisson(link = "log"), 
    offset = log(p1 + p3))

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-10.7514   -1.2475   -0.4445    0.7999   20.2746  

Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept) -2.046e+01  1.686e-02 -1213.7   <2e-16 ***
x1           2.829e-01  5.061e-04   558.9   <2e-16 ***
x2           8.147e-02  3.977e-04   204.9   <2e-16 ***
x5          -2.964e-01  1.663e-03  -178.3   <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: 970727  on 70107  degrees of freedom
Residual deviance: 245088  on 70104  degrees of freedom
AIC: 428036

Number of Fisher Scoring iterations: 5

まとめ

　オフセット項として65歳以上人口の対数を投入すると，総人口を投入するよりも予測性能が改善する可能性がある．