ねこすたっと

ねこの気持ちと統計について悩む筆者の備忘録的ページ。

「割合」や「率」に対してオフセット項付きポアソン回帰モデルを当てはめる [R]

アウトカムが「割合」、というか「試行数と成功数」として与えられているときに二項回帰モデルを当てはめる方法は以前まとめました。

今回は「オフセット項」を使った回帰モデルに当てはめる方法をまとめてみようと思います。

オフセット項(offset)とは

オフセット項とは回帰係数か1に固定される説明変数のことです。 ポアソン回帰モデルはアウトカムが発生数などのカウント変数のときに使われますが、観察された総人数や期間が違えば発生数も違って当然ですよね。 オフセット項を使えば、この「分母」をモデルに組み入れることができます。

通常の対数-ポアソン回帰モデルは、次の式で表されます。

 
\begin{aligned}
Y_i &\sim Poisson(\lambda_i) \\
\log(\lambda_i) &= X_i \beta
\end{aligned}

ポアソン分布は強度(intensity)λというパラメータで決まり、分布の平均・分散はともにλです。 対数-ポアソンモデルでは、λの対数が説明変数の線形和(=Xβ)で決まる、として係数βが推定されます。 つまり、発生数の期待値が説明変数で説明されています。

これに対して、オフセット項が入った対数-ポアソン回帰モデルは、次のようになります。

 
\begin{aligned}
Y_i &\sim Poisson(\lambda_i) \\
\log(\lambda_i) &= X_i \beta + \log(offset_i)
\end{aligned}

ここで、オフセット変数(offset)は、発生数の分母にくる観察された総人数や期間で、対数変換して加えます。 第2式を変形すると、

 
\begin{aligned}
\log \left( \frac{\lambda_i}{offset_i} \right) &= X_i \beta  \\
\end{aligned}

となり、このモデルでは「発生数/オフセット」の期待値が説明変数で説明されていることになります。

オフセット変数が観察総数のときは「発生数/オフセット = 発生割合」になりますし、オフセット変数が期間のときは「発生数/オフセット = 単位期間あたりの発生数」、つまり発生率になります。

使用するデータ

今回はsegmentedパッケージのdownデータを使います。

library(segmented)
data(down)

このデータは母体年齢と出生児のダウン症(21トリソミー)の発生の関連を示したデータで、以下の3つの変数を含みます。

  • age:母体の平均年齢
  • births:出生数
  • cases:ダウン症を有する児の数

実際のデータはこんな感じ。

> head(down)
   age births cases
1 17.0  13555    16
2 18.5  13675    15
3 19.5  18752    16
4 20.5  22005    22
5 21.5  23896    16
6 22.5  24667    12

このデータのように、「試行数」と「成功数」が与えられている場合に二項回帰モデルを当てはめる方法は以前まとめました。

necostat.hatenablog.jp

glm( )でオフセット項付きポアソン回帰モデルを当てはめる

分母に来てほしいbirthsをオフセット項に指定します(logを忘れないように)。 書き方は2通り。

(1) offset( )を使ってモデル式の中に書く

fit.pois <- glm(cases ~ age + offset(log(births)),
                family=poisson(link=log), data=down)

(2) 引数offsetとして指定する

fit.pois <- glm(cases ~ age, offset=log(births),
                family=poisson(link=log), data=down)

結果はいずれも以下のとおり。

> summary(fit.pois)

Call:
glm(formula = cases ~ age + offset(log(births)), family = poisson(link = log), 
    data = down)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-3.4070  -1.9389   0.5429   2.1276   4.7578  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept) -10.55277    0.21385  -49.35   <2e-16 ***
age           0.13713    0.00645   21.26   <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: 622.56  on 29  degrees of freedom
Residual deviance: 182.19  on 28  degrees of freedom
AIC: 325.26

Number of Fisher Scoring iterations: 5

おわりに

  • segmentedパッケージの勉強をしたついでに書いた記事なので、モデルの当てはまりは適切ではないかもしれません。
  • 飼い猫の予防接種に初めて付き添ってきました。「ツメ切り希望」と書きたかったのですが、「瓜切り」と書いていたことに後で気づきました。

参考資料

  • p.688にオフセットの記載があります。