ロジスティック回帰モデルを学び直す

分かったつもりで使っているけど、よくよく勉強し直したらあんまり分かってなかったことがよくあります。今回はロジスティック回帰モデル（ロジットモデル）について学び直してみました。

読んだ文献は、以下のものです（学び直しと言いながら1本だけです。すいません。）。

Norton, E.C. and Dowd, B.E. (2018), Log Odds and the Interpretation of Logit Models. Health Serv Res, 53: 859-878.

登場する関数いろいろ
2値アウトカムのモデル化
ロジスティック回帰モデルの解釈に関する注意点
おわりに
参考資料

登場する関数いろいろ

標準ロジスティック分布（standard logistic distribution）の確率密度関数（probability density function, PDF）は以下の数式で表されます。

$\begin{aligned} f(x) = \frac{e^x}{(e^x + 1)^2} \end{aligned}$

グラフを描くとこんな感じ。

curve(exp(x)/(exp(x)+1)^2,
      from=-10,to=10, 
      main="PDF of standard logistic distribution",
      ylab="density")

平均が0で左右対称な分布であるという点で標準正規分布（standard normal distribution）と同じ特徴を有しています。

標準ロジスティック分布の累積分布関数（cumulative distribution function, CDF）は、

$\begin{aligned} F(x) &= \int_{- \infty}^{x} \frac{e^{\alpha}}{(e^{\alpha} + 1)^2} d\alpha \\ &= \frac{e^x}{e^x +1} \\ &= \frac{1}{1 +e^{-x}} \end{aligned}$

となります。この関数は、標準ロジスティック関数（standard logistic function）と呼ばれます。

グラフを描くと、

curve(exp(x)/(exp(x)+1),
      from=-10,to=10, 
      main="CDF of standard logistic distribution",
      ylab="density")

となります。 y切片が0.5の単調増加関数で、標準正規分布のCDFと似たような形になります。ちなみに、標準正規分布のCDFは標準プロビット関数（standard probit function）と呼ばれます。

2値アウトカムのモデル化

潜在的な連続変数を介してモデル化する

線形モデルでは、連続変数 $y^{*}$ を説明変数の線型結合と誤差 $\varepsilon_i$ でモデル化します。

$\begin{aligned} y^*_i = \beta_0 + \beta_1 x_i + \varepsilon_i \end{aligned}$

アウトカム $y$ が2値変数のときはこの方法が使えないので、連続変数 $y^{*}$ を潜在変数として、次のように変換するとしましょう。

${ \begin{eqnarray} y_{i} = \left\{ \begin{array}{ll} 1 & (y^*_i > 0) \\ 0 & (y^*_i \leq 0) \end{array}\right. \end{eqnarray} }$

このときアウトカムが発生する確率は、

$\begin{aligned} Pr(y_i = 1 | x_i) &= Pr(y^*_i > 0 | x_i) \\ &= Pr(\beta_0 + \beta_1 x_i + \varepsilon_i > 0 | x_i) \\ &= Pr(\varepsilon_i > -(\beta_0 + \beta_1 x_i) | x_i) \end{aligned}$

となります。ここで、誤差 $\varepsilon$ の確率分布が0を中心に左右対称であれば、

$\begin{aligned} Pr(\varepsilon_i > -(\beta_0 + \beta_1 x_i) | x_i) = Pr(\varepsilon_i < (\beta_0 + \beta_1 x_i) | x_i) \end{aligned}$

となるので、以下のようになります。

$\begin{aligned} Pr(y_i = 1 | x_i) = Pr(\varepsilon_i < (\beta_0 + \beta_1 x_i) | x_i) \end{aligned}$

右辺は $\varepsilon$ のCDFを使えば計算できます。

標準化誤差の分布を仮定する

ここで問題になるのが、誤差 $\varepsilon$ の分布が分からないということです。仮に誤差 $\varepsilon$ が正規分布に従うと仮定したとしても、分散が分からなければ $\beta_0 + \beta_1 x_i$ における累積確率がどれくらいなのか分かりません。

そこで、誤差 $\varepsilon$ の標準誤差 $\sigma$ で割って標準化するステップが必要になります。

$\begin{aligned} Pr(y_i = 1 | x_i) = Pr \left( \frac{\varepsilon_i}{\sigma} < \frac{\beta_0 + \beta_1 x_i}{\sigma} | x_i \right) \end{aligned}$

こうしておいて、標準化した誤差 $\varepsilon/\sigma$ が標準ロジスティック分布に従うと考えれば、

$\begin{aligned} Pr(y_i = 1 | x_i) = \frac{1}{1+\exp \left( -\frac{\beta_0 + \beta_1 x_i}{\sigma} \right)} \end{aligned}$

となりますし、標準正規分布に従うと仮定すれば、

$\begin{aligned} Pr(y_i = 1 | x_i) = \Phi \left( \frac{\beta_0 + \beta_1 x_i}{\sigma} \right) \end{aligned}$

となります（ $\Phi(x)$ は標準正規分布のCDF）。

オッズ比の導出

ロジスティック回帰モデルを使うとオッズは次のようになります。

$\begin{aligned} \frac{p_i}{1-p_i} &= \frac{ \frac{1}{1+\exp \left( -\frac{\beta_0 + \beta_1 x_i}{\sigma} \right)}}{ \frac{\exp \left( -\frac{\beta_0 + \beta_1 x_i}{\sigma} \right)}{1+\exp \left( -\frac{\beta_0 + \beta_1 x_i}{\sigma} \right)}} \\ &= \frac{1}{\exp \left( -\frac{\beta_0 + \beta_1 x_i}{\sigma} \right)} \\ &= \exp \left( \frac{\beta_0 + \beta_1 x_i}{\sigma} \right) \end{aligned}$

さらにx=0のときに対するx=1のときのオッズ比ORは、

$\begin{aligned} OR &= \frac{\exp \left( \frac{\beta_0 + \beta_1}{\sigma} \right)}{\exp \left( \frac{\beta_0}{\sigma} \right)} \\ &= \exp \left( \frac{\beta_1}{\sigma} \right) \end{aligned}$