Chapter 21

ロジスティック回帰と一般化線形モデル

この章がなぜ必要なのか——0 か 1 かを予測する.
「この患者は病気か」「この顧客は解約するか」「このメールはスパムか」—— 目的変数が 0/1 の二値である問題は、実務で圧倒的に多い。
ここに普通の線形回帰を当てると壊れる。予測値が 1.3 や −0.2 になり、確率として意味をなさなくなるからだ。
ロジスティック回帰は、線形モデルを確率の世界に翻訳することでこれを解決する。そしてその仕組みを一般化したのが GLM（一般化線形モデル）で、カウントデータも比率も同じ枠組みで扱えるようになる。
この章で使う既出の用語（定義は各リンク先）. 目的変数・説明変数（01 章 4 節）、オッズ・オッズ比・リスク比（03 章 6 節）、有病率・感度・特異度・陽性的中率・基準率（05 章 3 節）、期待値・分散（06 章 3〜4 節）、ベルヌーイ・二項・ポアソン・負の二項分布（07 章）、正規分布と標準正規分布の分布関数 $\Phi$（08 章 4 節）、χ² 分布（08 章 7 節）、尤度・対数尤度・最尤法・漸近正規性（10 章 3 節, 5 節）、信頼区間・標準誤差（11 章）、検定・p 値・自由度（12 章、18 章 4 節）、マン・ホイットニーの $U$（16 章）、計画行列 $X$・$\mathbf x_i$・ベクトル微分・正規方程式（19 章 1〜2 節）、多重共線性（19 章 8 節）、Ridge 正則化（20 章 7 節）、ベイズ推定と弱情報事前分布（23 章）

1. なぜ線形回帰ではだめなのか

$y\in\{0,1\}$ に $y=\beta_0+\beta_1x+\varepsilon$ を当てはめると:

問題	内容
範囲を外れる	予測値が 0 未満・1 超になりうる
不均一分散	$V[y]=p(1-p)$ は $p$ に依存する（07 章）
正規性が成り立たない	誤差は 2 値しか取らない
線形性が不自然	確率が 0.5 付近と 0.95 付近で同じ傾きなのは不自然

4 番目が本質的である. 確率は 0 と 1 で頭打ちになる。「勉強時間を 1 時間増やすと合格率が 10 ポイント上がる」が合格率 95% の人にも成り立つはずがない。 端に近づくほど効きが鈍るという性質を、モデルに組み込む必要がある。

2. ロジスティック回帰のモデル

オッズとロジット

変換	式	範囲
確率 $p$	—	$[0,1]$
オッズ	$\dfrac{p}{1-p}$	$[0,\infty)$
ロジット	$\log\dfrac{p}{1-p}$	$(-\infty,\infty)$

ロジット変換で、確率が実数全体に引き伸ばされる。 ここに線形モデルを当てればよい。

\[ \boxed{\log\frac{p_i}{1-p_i}=\beta_0+\beta_1x_{i1}+\cdots+\beta_px_{ip}=\eta_i} \]

右辺の線形式に $\eta_i$（イータ）という名前を付けた。線形予測子と呼び、以下ずっと使う。

逆に解く（ロジスティック関数）

\[ \frac{p}{1-p}=e^{\eta} \quad\Longrightarrow\quad p = e^\eta(1-p) \quad\Longrightarrow\quad p(1+e^\eta)=e^\eta \]

\[ \boxed{p = \frac{e^\eta}{1+e^\eta}=\frac{1}{1+e^{-\eta}} = \sigma(\eta)} \]

（2 つ目の等号は分母・分子を $e^\eta$ で割った。）これがシグモイド関数（ロジスティック関数）である。名前が紛らわしいので整理する: ロジットは確率 $p$ を実数 $\eta$ に変える関数（$\log\frac{p}{1-p}$）、ロジスティック関数はその逆で $\eta$ を $p$ に戻す関数。

性質	内容
値域	$(0,1)$ — 確率として必ず正しい
$\eta=0$	$p=0.5$
単調増加	$\eta$ が増えれば $p$ も増える
導関数	$\sigma'(\eta)=\sigma(\eta)(1-\sigma(\eta))$ — 中央で最大、端で 0。導出: $\sigma = (1+e^{-\eta})^{-1}$ を合成関数の微分で $\sigma' = e^{-\eta}/(1+e^{-\eta})^2 = \dfrac{1}{1+e^{-\eta}}\cdot\dfrac{e^{-\eta}}{1+e^{-\eta}} = \sigma(1-\sigma)$（$\frac{e^{-\eta}}{1+e^{-\eta}} = 1 - \sigma$）

導関数の形が「端で効きが鈍る」を表現している. $p=0.5$ で傾き最大（0.25）、$p=0.95$ では 0.0475 と 5 分の 1 以下。先ほど「不自然だ」と述べた性質が、自然にモデルに入っている。

3. 係数の解釈 — オッズ比

$x_1$ を 1 増やすと:

\[ \log\frac{p'}{1-p'}-\log\frac{p}{1-p}=\beta_1 \quad\Longrightarrow\quad \frac{\text{odds}'}{\text{odds}}=e^{\beta_1} \]

\[ \boxed{e^{\beta_1}= \text{オッズ比 (odds ratio)}} \]

$\beta_1$	$e^{\beta_1}$	意味
0	1	効果なし
0.69	2.0	オッズが 2 倍
−0.69	0.5	オッズが半分
0.1	1.105	オッズが約 10% 増（$\beta$ が小さいなら、$e^x$ の $x = 0$ での接線近似 $e^\beta\approx 1+\beta$）

信頼区間: $\beta$ の区間を作ってから指数変換する。

\[ \left[e^{\hat\beta-1.96\text{SE}},\ e^{\hat\beta+1.96\text{SE}}\right] \]

（$\hat\beta$ は最尤推定量なので大標本で正規分布に近づく（10 章の漸近正規性）が、それを指数変換した $e^{\hat\beta}$ は右に歪んだ分布になり正規近似が効かない。区間は $e^{\hat\beta}$ を中心にしない。）

オッズ比をリスク比と読み替えてはいけない（03 章）。結果がまれ（$p<0.1$ 程度）なら $\text{OR}\approx\text{RR}$ だが、 $p$ が大きいと OR は RR を大幅に誇張する。
例: 対照群のリスク 40%、処置群 60%。RR = 1.5 だが OR = $\frac{0.6/0.4}{0.4/0.6}$ = 2.25。「2.25 倍起きやすい」と報告するのは誤りである。

4. パラメータ推定 — 最尤法

尤度関数

$y_i\in\{0,1\}$ がベルヌーイに従うので（07 章）

\[ L(\boldsymbol\beta)=\prod_{i=1}^n p_i^{y_i}(1-p_i)^{1-y_i} \]

\[ \boxed{\ell(\boldsymbol\beta)=\sum_{i=1}^n\left[y_i\log p_i + (1-y_i)\log(1-p_i)\right]} \]

これは機械学習で言う「交差エントロピー損失」の符号を反転したものである. ニューラルネットワークの二値分類で使う損失関数は、 ロジスティック回帰の対数尤度そのものである。深層学習の出力層に sigmoid + binary cross-entropy を使うのは、最後の層でロジスティック回帰をしているのと同じ意味になる。

尤度方程式

$p_i=\sigma(\mathbf x_i^\top\boldsymbol\beta)$、$\sigma'=\sigma(1-\sigma)$ を使って微分すると、驚くほどきれいな形になる。

\[ \frac{\partial\ell}{\partial\boldsymbol\beta}=\sum_{i=1}^n (y_i-p_i)\mathbf x_i = X^\top(\mathbf y-\mathbf p)=\mathbf 0 \]

導出: $\ell$ を $\eta_i$ で微分すると

\[ \frac{\partial\ell}{\partial\eta_i}=\left(\frac{y_i}{p_i}-\frac{1-y_i}{1-p_i}\right)p_i(1-p_i) = y_i(1-p_i)-(1-y_i)p_i = y_i-p_i \]

連鎖律で $\partial\eta_i/\partial\boldsymbol\beta=\mathbf x_i$ を掛ければよい。∎

線形回帰の正規方程式 $X^\top(\mathbf y - X\boldsymbol\beta)=\mathbf 0$ とそっくりである. 違いは $X\boldsymbol\beta$ が $\mathbf p = \sigma(X\boldsymbol\beta)$ になっただけ。 「残差 × 説明変数の和が 0」という構造は共通している。

解析解は存在しない

$\mathbf p$ が $\boldsymbol\beta$ の非線形関数なので、閉じた解はない。ニュートン・ラフソン法で反復的に解く。これは「現在の点 $\boldsymbol\beta^{(t)}$（上付き $(t)$ は反復回数）で目的関数を 2 次式で近似し、その 2 次式の頂点へ移る」を繰り返す方法で、1 変数なら $x^{(t+1)} = x^{(t)} - f'(x^{(t)})/f''(x^{(t)})$、多変数なら 1 階微分（勾配 $\mathbf g$）と 2 階微分を並べた行列（ヘッセ行列 $\mathcal H$）で $\boldsymbol\beta^{(t+1)} = \boldsymbol\beta^{(t)} - \mathcal H^{-1}\mathbf g$ とする。勾配は上で求めた $\mathbf g = X^\top(\mathbf y-\mathbf p)$。ヘッセ行列は $\mathbf g$ をもう一度 $\boldsymbol\beta$ で微分して、$\partial p_i/\partial\boldsymbol\beta = p_i(1-p_i)\mathbf x_i$ より

\[ \mathcal H = -\sum_i p_i(1-p_i)\,\mathbf x_i\mathbf x_i^\top = -X^\top WX, \qquad W=\text{diag}(p_i(1-p_i)) \]

よって更新式は

\[ \boldsymbol\beta^{(t+1)}=\boldsymbol\beta^{(t)}+(X^\top WX)^{-1}X^\top(\mathbf y-\mathbf p) \]

これは反復再重み付け最小二乗法 (IRLS) と呼ばれる。理由: $\boldsymbol\beta^{(t)} = (X^\top WX)^{-1}X^\top W X\boldsymbol\beta^{(t)}$ と書き直して右辺をまとめると

\[ \boldsymbol\beta^{(t+1)} = (X^\top WX)^{-1}X^\top W\mathbf z, \qquad \mathbf z = X\boldsymbol\beta^{(t)} + W^{-1}(\mathbf y-\mathbf p) \]

となり、これは「作業応答 $\mathbf z$ を $X$ に、重み $W$ で重み付き最小二乗回帰する」式そのものである。重み $W$ と $\mathbf z$ を毎回作り直しながら繰り返すので「反復再重み付け」。統計ソフトの内部実装はほぼこれである。

完全分離 (complete separation)

$x$ の値で 0/1 が完全に分かれてしまうと、尤度が最大値を持たず $\hat\beta\to\infty$ に発散する。

症状: 係数が 20 とか 40 という異常な値になり、標準誤差が桁違いに大きくなる。「完璧に予測できた」ように見えるが、実際には推定が破綻している。
対策: 変数を減らす、Firth のペナルティ付き尤度（対数尤度に係数の発散を抑える罰則項 $\frac{1}{2}\log\det(X^\top WX)$ を足して最大化する方法）、Ridge 正則化、ベイズ推定（弱情報事前分布）。小標本や希少事象で起きやすいので、係数の異常値には常に注意する。

5. モデルの評価

逸脱度 (deviance)

\[ D = -2\ell(\hat{\boldsymbol\beta}) \]

線形回帰の残差平方和に対応する量。小さいほど当てはまりが良い。

種類	意味
帰無逸脱度	切片だけのモデル
残差逸脱度	当てはめたモデル

尤度比検定

入れ子のモデルの比較:

\[ \Lambda = D_{\text{reduced}}-D_{\text{full}} \sim \chi^2_q \]

（$q$ = 増やしたパラメータ数。$\chi^2_q$ は自由度 $q$ の χ² 分布（08 章）で、$H_0$（増やしたパラメータはすべて 0）のもとで大標本近似として成り立つ。）12 章のネイマン・ピアソンの補題に対応する検定である。

擬似 $R^2$

以下、$\ell$ は対数尤度、$L = e^{\ell}$ は対数を取る前の尤度。添字 null は切片だけのモデル、full は当てはめたモデル。

名前	式
マクファデン	$1-\dfrac{\ell_{\text{full}}}{\ell_{\text{null}}}$
コックス・スネル	$1-(L_{\text{null}}/L_{\text{full}})^{2/n}$
ナゲルケルケ	コックス・スネルを最大 1 に正規化

線形回帰の $R^2$ と同じ感覚で読んではいけない. マクファデンの擬似 $R^2$ は 0.2〜0.4 でも「良好」とされる。 実務では判別性能（次項）の方が有用である。

判別性能

	予測 1	予測 0
実際 1	TP	FN
実際 0	FP	TN

指標	式	意味
感度（再現率）	$\frac{TP}{TP+FN}$	陽性を取りこぼさない率
特異度	$\frac{TN}{TN+FP}$	陰性を正しく陰性と言う率
適合率	$\frac{TP}{TP+FP}$	陽性と言ったうち本物の率
正解率	$\frac{TP+TN}{n}$	全体の正答率

不均衡データでは正解率が無意味になる. 陽性が 1% しかないなら、全部「陰性」と答えるだけで正解率 99% になる。 05 章の基準率の話がここでも効いている。 適合率も有病率に依存する（05 章の陽性的中率——陽性と判定されたうち本当に陽性である割合——と定義が同じ量である）。

ROC 曲線と AUC

閾値を 0 から 1 まで動かしながら、横軸に $1-$ 特異度、縦軸に感度をプロットする。

\[ \text{AUC} = \text{ROC 曲線の下側面積} \]

AUC	判断
0.5	ランダムと同じ
0.7〜0.8	許容範囲
0.8〜0.9	良好
0.9 以上	優秀

AUC の美しい解釈.
$$\text{AUC} = P(\text{陽性例のスコア} > \text{陰性例のスコア})$$ ランダムに陽性 1 つ、陰性 1 つを選んだとき、陽性の方が高いスコアを得る確率。
これは 16 章のマン・ホイットニー $U$ 統計量とまったく同じ量である（$\text{AUC}=U/(n_1n_2)$、$n_1, n_2$ は陽性例と陰性例の個数）。分野が違えば名前も違うが、中身は同じである。
AUC は閾値に依存しないのが利点。ただし不均衡データでは楽観的に見えるので、 PR 曲線（適合率-再現率曲線）も併用するとよい。

6. 一般化線形モデル (GLM)

ロジスティック回帰を含む、より広い枠組み。3 つの部品でできている。

部品	内容
確率分布	指数型分布族——確率（密度）が $\exp\{(y\theta - b(\theta))/\phi + c(y,\phi)\}$ の形に書ける分布の総称（正規・二項・ポアソン・ガンマ…がすべてこの形になる）
線形予測子	$\eta = X\boldsymbol\beta$
リンク関数	$g(\mu)=\eta$（$\mu = E[y]$ は目的変数の期待値）

主な GLM

モデル	分布	リンク	$y$ の型
線形回帰	正規	恒等 $\mu=\eta$	連続
ロジスティック回帰	二項	ロジット	0/1
プロビット回帰	二項	$\Phi^{-1}(\mu)$（$\Phi$ は標準正規分布の分布関数、08 章）	0/1
ポアソン回帰	ポアソン	$\log\mu=\eta$	カウント
負の二項回帰	負の二項	$\log$	過分散カウント
ガンマ回帰	ガンマ	$\log$ or 逆数	正の連続（保険金額など）

リンク関数の役割. $\mu$ の取りうる範囲（$[0,1]$ や $[0,\infty)$）を、 実数全体に引き伸ばすこと。そこに線形モデルを置く。ロジット関数がやっていたのはまさにこれで、GLM はそれを一般化しただけである。

ポアソン回帰

\[ \log \mu_i=\mathbf x_i^\top\boldsymbol\beta \quad\Longrightarrow\quad \mu_i = E[y_i]=e^{\mathbf x_i^\top\boldsymbol\beta} \]

係数の解釈: $e^{\beta_j}$ は「$x_j$ が 1 増えると期待件数が何倍になるか」（率比）。

オフセット項: 観測期間や人口が違う場合、

\[ \log \mu_i=\log(t_i)+\mathbf x_i^\top\boldsymbol\beta \]

とすると、実質的に「単位時間あたりの率」をモデル化できる。

過分散への対処

ポアソンは平均 = 分散（07 章）。実データでは分散の方が大きいことが多い。

対処	内容
準ポアソン	分散を $\phi\mu$ とし、標準誤差だけ調整
負の二項回帰	分散 $\mu+\mu^2/r$（$r$ は負の二項分布の形状パラメータ、07 章。小さいほど過分散が大きい）を明示的にモデル化。推奨
ゼロ過剰モデル	0 が異常に多い場合（ZIP, ZINB）

過分散を無視するとどうなるか. 標準誤差を過小評価し、 偽の有意を量産する。カウントデータを扱うときは、残差逸脱度 / 自由度（自由度 = データ数 − パラメータ数）が 1 を大きく超えていないか必ず確認する（モデルが正しければ残差逸脱度は近似的に自由度を期待値とする χ² 分布に従うので、比は 1 付近になるはず）。

7. 多クラス分類

$k$ 個のカテゴリを予測する場合。

多項ロジスティック回帰（ソフトマックス回帰）

\[ P(y=j) = \frac{e^{\eta_j}}{\sum_{m=1}^k e^{\eta_m}} \]

基準カテゴリを 1 つ決め、その係数を 0 に固定する。

これがニューラルネットワークの出力層の softmax である. 多クラス分類の深層学習モデルは、最終層で多項ロジスティック回帰をしている。

順序ロジスティック回帰

カテゴリに順序がある場合（満足度 1〜5 など）。比例オッズモデル:

\[ \log\frac{P(y\le j)}{P(y>j)}=\alpha_j - \mathbf x^\top\boldsymbol\beta \]

左辺は「$j$ 以下」のオッズなので、$\boldsymbol\beta$ の前をマイナスにしておくと「$\beta > 0$ なら $x$ が大きいほど高いカテゴリになりやすい」と、通常のロジスティック回帰と同じ向きで読める（符号の約束）。カテゴリの境目を表す定数 $\alpha_j$（5 節の ROC の閾値とは別物。ここでは潜在的な連続量をカテゴリに切る区切り）だけがカテゴリごとに違い、傾き $\boldsymbol\beta$ は共通。この仮定（比例オッズ仮定）は、境目 $j$ ごとに別々の 2 値ロジスティック回帰を当てはめて傾きが等しいかを検定するブラント検定などで確認する。

8. 実務での注意

注意点	内容
必要サンプル数	少ない方のクラスで、変数 1 個あたり 10〜20 例（EPV ルール）
不均衡データ	重み付け、閾値調整。アンダーサンプリングは慎重に
完全分離	係数が異常値。Firth 法や正則化
多重共線性	線形回帰と同じ問題が起きる
キャリブレーション	予測確率が実際の頻度と合っているか（信頼性曲線。下のキャリブレーションプロットと同じもの）

判別できることと確率が正しいことは別である. AUC が高くても、予測確率が実際より系統的に高い/低いことはよくある。「この患者の発症確率は 30%」という数字を意思決定に使うなら、 キャリブレーションプロット（予測確率をビンに分け、実際の発生率と比較）で確認する。

9. まとめ

項目	内容
モデル	$\log\frac{p}{1-p}=\mathbf x^\top\boldsymbol\beta$
確率	$p=\sigma(\eta)=1/(1+e^{-\eta})$
係数	$e^{\beta_j}$ = オッズ比
推定	最尤法（IRLS で数値解）
尤度	$\sum[y\log p+(1-y)\log(1-p)]$ = 交差エントロピー
評価	逸脱度、AUC、感度・特異度

確率をロジット変換して実数全体に伸ばし、そこに線形モデルを当てる。
係数の指数がオッズ比。リスク比と混同しない。
対数尤度は機械学習の交差エントロピー損失と同じもの。
GLM は「分布 + リンク関数 + 線形予測子」の 3 点セット。
AUC = $P(\text{陽性のスコア} > \text{陰性のスコア})$ = マン・ホイットニー $U$ の正規化版。
不均衡データでは正解率を使わない。カウントデータでは過分散を必ず確認する。

変換	式	範囲
確率 \(p\)	—	\([0,1]\)
オッズ	\(\dfrac{p}{1-p}\)	\([0,\infty)\)
ロジット	\(\log\dfrac{p}{1-p}\)	\((-\infty,\infty)\)

性質	内容
値域	\((0,1)\) — 確率として必ず正しい
\(\eta=0\)	\(p=0.5\)
単調増加	\(\eta\) が増えれば \(p\) も増える
導関数	\(\sigma'(\eta)=\sigma(\eta)(1-\sigma(\eta))\) — 中央で最大、端で 0。導出: \(\sigma = (1+e^{-\eta})^{-1}\) を合成関数の微分で \(\sigma' = e^{-\eta}/(1+e^{-\eta})^2 = \dfrac{1}{1+e^{-\eta}}\cdot\dfrac{e^{-\eta}}{1+e^{-\eta}} = \sigma(1-\sigma)\)（\(\frac{e^{-\eta}}{1+e^{-\eta}} = 1 - \sigma\)）

名前	式
マクファデン	\(1-\dfrac{\ell_{\text{full}}}{\ell_{\text{null}}}\)
コックス・スネル	\(1-(L_{\text{null}}/L_{\text{full}})^{2/n}\)
ナゲルケルケ	コックス・スネルを最大 1 に正規化

指標	式	意味
感度（再現率）	\(\frac{TP}{TP+FN}\)	陽性を取りこぼさない率
特異度	\(\frac{TN}{TN+FP}\)	陰性を正しく陰性と言う率
適合率	\(\frac{TP}{TP+FP}\)	陽性と言ったうち本物の率
正解率	\(\frac{TP+TN}{n}\)	全体の正答率

部品	内容
確率分布	指数型分布族——確率（密度）が \(\exp\{(y\theta - b(\theta))/\phi + c(y,\phi)\}\) の形に書ける分布の総称（正規・二項・ポアソン・ガンマ…がすべてこの形になる）
線形予測子	\(\eta = X\boldsymbol\beta\)
リンク関数	\(g(\mu)=\eta\)（\(\mu = E[y]\) は目的変数の期待値）

モデル	分布	リンク	\(y\) の型
線形回帰	正規	恒等 \(\mu=\eta\)	連続
ロジスティック回帰	二項	ロジット	0/1
プロビット回帰	二項	\(\Phi^{-1}(\mu)\)（\(\Phi\) は標準正規分布の分布関数、08 章）	0/1
ポアソン回帰	ポアソン	\(\log\mu=\eta\)	カウント
負の二項回帰	負の二項	\(\log\)	過分散カウント
ガンマ回帰	ガンマ	\(\log\) or 逆数	正の連続（保険金額など）

対処	内容
準ポアソン	分散を \(\phi\mu\) とし、標準誤差だけ調整
負の二項回帰	分散 \(\mu+\mu^2/r\)（\(r\) は負の二項分布の形状パラメータ、07 章。小さいほど過分散が大きい）を明示的にモデル化。推奨
ゼロ過剰モデル	0 が異常に多い場合（ZIP, ZINB）

項目	内容
モデル	\(\log\frac{p}{1-p}=\mathbf x^\top\boldsymbol\beta\)
確率	\(p=\sigma(\eta)=1/(1+e^{-\eta})\)
係数	\(e^{\beta_j}\) = オッズ比
推定	最尤法（IRLS で数値解）
尤度	\(\sum[y\log p+(1-y)\log(1-p)]\) = 交差エントロピー
評価	逸脱度、AUC、感度・特異度