統計 23 · ベイズ統計 — 事前分布から事後分布へ

Chapter 23

ベイズ統計 — 事前分布から事後分布へ

この章がなぜ必要なのか——「知りたいこと」に直接答える枠組み.

11 章で見たとおり、頻度論の 95% 信頼区間は 「真値がこの範囲にある確率が 95%」とは言えなかった。 12 章の p 値も「\(H_0\) が正しい確率」ではなかった。

だが私たちが本当に知りたいのは、まさにその 「このパラメータがこの範囲にある確率」なのである。

ベイズ統計はパラメータ自体を確率変数と見ることで、その問いに直接答える。 代償として事前分布を置く必要があり、そこが 100 年にわたる論争の的でもあった。 現在では、計算機の発達によって実務でも広く使われる標準的な道具になっている。

この章で使う既出の用語(定義は各リンク先). 標本(01 章 2 節)、オッズ・オッズ比(03 章 6 節)、ベイズの定理・事前確率・事後確率・尤度(05 章 2 節)、確率変数(06 章 1 節)、期待値・分散(06 章 3〜4 節)、正規分布(08 章 4 節)、ガンマ分布(08 章 6 節)、ベータ分布(08 章 10 節)、標本平均の分散(09 章)、最尤推定・フィッシャー情報量 \(I(\theta)\)(10 章 2〜3 節)、信頼区間・予測区間(11 章)、仮説 \(H_0, H_1\)・p 値(12 章)、オプショナル・ストッピング(17 章 4 節)、Ridge・Lasso・正則化(20 章 7〜8 節)、ロジスティック回帰の係数 \(\beta\) とオッズ比(21 章)

1. 基本の枠組み

05 章のベイズの定理を、パラメータ \(\theta\) とデータ \(D\) に適用する。

\[ \boxed{p(\theta\mid D)=\frac{p(D\mid\theta)p(\theta)}{p(D)}\propto \underbrace{p(D\mid\theta)}_{\text{尤度}}\times\underbrace{p(\theta)}_{\text{事前分布}}} \]

(\(\propto\) は「比例する」——\(\theta\) によらない定数倍を除いて等しい、の意味。分母 \(p(D)\) は \(\theta\) を含まないので省ける。) 05 章の事後確率は 1 つの事象についての値だったが、ここでは \(\theta\) の値ごとの事後確率(密度)を並べたもの、すなわち分布として扱う。

記号名前意味
\(p(\theta)\)事前分布 (prior)データを見る前の \(\theta\) についての信念
\(p(D\mid\theta)\)尤度10 章と同じもの
\(p(\theta\mid D)\)事後分布 (posterior)データを見た後の信念。これが答え
\(p(D)=\int p(D\mid\theta)p(\theta)d\theta\)周辺尤度(証拠)正規化定数

頻度論との根本的な違い

頻度論ベイズ
\(\theta\)固定された未知の定数確率変数(06 章。値が確率的に決まる量)
データ確率変数(繰り返しを想像)観測されたもの(固定)
出力点推定 + 信頼区間事後分布そのもの
区間の意味手順の長期成功率その区間に \(\theta\) がある確率
事前情報使えない明示的に使う

ベイズの最大の利点は「答えが分布として出る」ことである. 点推定と区間を別々に計算するのではなく、 \(\theta\) についての完全な情報が事後分布 1 つに入っている。 そこから平均でも中央値でも区間でも、任意の確率でも取り出せる。

「この施策の効果が 5% を超えている確率は?」という問いに、 事後分布があれば \(P(\theta>0.05\mid D)\) として直接答えられる。 頻度論ではこの問い自体が定義できない。

2. 共役事前分布 — 計算が閉じる場合

事前分布と事後分布が同じ分布族になるとき、その事前分布を共役 (conjugate) と呼ぶ。

尤度共役事前事後
二項 \(\text{Bin}(n,p)\)\(\text{Beta}(\alpha,\beta)\)\(\text{Beta}(\alpha+s,\ \beta+f)\)
ポアソン \(\text{Po}(\lambda)\)\(\text{Gamma}(\alpha,\beta)\)\(\text{Gamma}(\alpha+\sum x_i,\ \beta+n)\)
正規(\(\sigma^2\) 既知)\(N(\mu_0,\tau_0^2)\)\(N(\mu_n,\tau_n^2)\)
正規(\(\mu\) 既知、\(\sigma^2\) を推定)逆ガンマ(\(X\) がガンマ分布に従うときの \(1/X\) の分布)逆ガンマ

ベータ・二項モデル(完全な導出)

事前分布を \(\theta\sim\text{Beta}(\alpha,\beta)\) とする(08 章)。

\[ p(\theta)\propto \theta^{\alpha-1}(1-\theta)^{\beta-1} \]

\(n\) 回中 \(s\) 回成功を観測。尤度は

\[ p(D\mid\theta)\propto\theta^s(1-\theta)^{n-s} \]

掛け合わせると

\[ p(\theta\mid D)\propto\theta^{\alpha+s-1}(1-\theta)^{\beta+n-s-1} \]
\[ \boxed{\theta\mid D \sim \text{Beta}(\alpha+s,\ \beta+n-s)} \]

∎

「ただ足すだけ」である. 成功回数を \(\alpha\) に、失敗回数を \(\beta\) に足す。 理由を形で見ると、事前分布の指数 \(\theta^{\alpha-1}(1-\theta)^{\beta-1}\) は「\(\alpha-1\) 回成功、\(\beta-1\) 回失敗」の尤度 \(\theta^{\alpha-1}(1-\theta)^{\beta-1}\) と同じ形をしている。 そこへ実際の \(s\) 回成功・\(n-s\) 回失敗を掛けると指数が \((\alpha-1)+s\)、\((\beta-1)+(n-s)\) になる——これが「足すだけ」の正体である。 だから \(\text{Beta}(\alpha,\beta)\) の事前分布は 「すでに \(\alpha-1\) 回成功、\(\beta-1\) 回失敗を見たのと同じ情報量」と解釈できる。 \(\text{Beta}(1,1)\) = 一様分布は「0 回成功、0 回失敗 = 何も知らない」に対応する。

事後平均は加重平均になる

\(\text{Beta}(a,b)\) の平均は \(a/(a+b)\)(08 章)なので、事後分布 \(\text{Beta}(\alpha+s,\beta+n-s)\) の平均は

\[ E[\theta\mid D]=\frac{\alpha+s}{\alpha+\beta+n} \]

変形すると

\[ E[\theta\mid D]=\underbrace{\frac{n}{\alpha+\beta+n}}_{w}\cdot\underbrace{\frac{s}{n}}_{\text{データ}}+\underbrace{\frac{\alpha+\beta}{\alpha+\beta+n}}_{1-w}\cdot\underbrace{\frac{\alpha}{\alpha+\beta}}_{\text{事前平均}} \]
\[ \boxed{\text{事後平均} = \text{データの推定値と事前平均の加重平均}} \]

\(n\) が大きくなるほどデータの重みが増す. \(n\to\infty\) で事後平均は最尤推定値 \(s/n\) に収束する。 データが十分にあれば、事前分布の影響は消える。 これが「事前分布は主観的だ」という批判への、実務的な回答の 1 つである。

逆に \(n\) が小さいときは事前分布が効く。それは欠点ではなく、 少ないデータで極端な結論を出さないための安全装置として働く。

具体例: 新しい薬を 10 人に投与して 8 人に効果があった。

最尤推定の 0.8 は「10 人中 8 人」だけを根拠にした値で、過信しやすい. ベイズは自動的に中央(0.5)寄りに引き戻す。これは 20 章の正則化と同じ働きであり、 18 章の「平均への回帰」とも通じる考え方である。

正規モデル(\(\sigma^2\) 既知)

事前 \(\mu\sim N(\mu_0,\tau_0^2)\)、データ \(x_1,\dots,x_n\) は独立に \(N(\mu,\sigma^2)\)(\(\sigma^2\) は既知)、標本平均を \(\bar x\) とする。

導出: 尤度は \(\prod_i \exp\{-(x_i-\mu)^2/2\sigma^2\}\) に比例する。\(\sum_i(x_i-\mu)^2 = \sum_i(x_i-\bar x)^2 + n(\bar x-\mu)^2\)(交差項は \(\sum(x_i-\bar x)=0\) で消える)で、第 1 項は \(\mu\) を含まないので

\[ p(D\mid\mu)\propto\exp\left\{-\frac{n(\bar x-\mu)^2}{2\sigma^2}\right\} \]

事前分布 \(\exp\{-(\mu-\mu_0)^2/2\tau_0^2\}\) を掛け、指数の中を \(\mu\) について整理する(\(\mu\) を含まない項は定数として落とす):

\[ -\frac{1}{2}\left[\mu^2\left(\frac{1}{\tau_0^2}+\frac{n}{\sigma^2}\right)-2\mu\left(\frac{\mu_0}{\tau_0^2}+\frac{n\bar x}{\sigma^2}\right)\right]+\text{const} \]

これは \(-\frac{(\mu-\mu_n)^2}{2\tau_n^2}\) の形(平方完成)で、\(\mu^2\) の係数と \(\mu\) の係数を見比べると

\[ \mu\mid D\sim N(\mu_n,\tau_n^2) \]
\[ \frac{1}{\tau_n^2}=\frac{1}{\tau_0^2}+\frac{n}{\sigma^2}, \qquad \mu_n = \tau_n^2\left(\frac{\mu_0}{\tau_0^2}+\frac{n\bar x}{\sigma^2}\right) \]

∎

精度(分散の逆数)が足し算になる. 情報は精度として蓄積される、というのがベイズの美しい表現である。 事後平均は「精度で重み付けした平均」になっている。

09 章の \(V[\bar X]=\sigma^2/n\) がここに現れている。 データの精度は \(n/\sigma^2\) ——標本を増やせば精度が線形に増える。

3. 事前分布の選び方

種類内容使いどころ
無情報事前一様分布など事前知識がない
ジェフリーズ事前\(p(\theta)\propto\sqrt{I(\theta)}\)(\(I(\theta)\) はフィッシャー情報量。10 章)変数変換に不変。理論的に自然
弱情報事前広めだが極端を排除(例: \(N(0,10^2)\))実務での既定。推奨
情報事前過去の研究や専門知識を反映情報がある場合

「無情報事前」は思ったより無情報ではない. \(\theta\) に一様分布を置くと、\(\log\theta\) や \(1/\theta\) には一様でなくなる。 どの尺度で「平ら」にするかという選択がすでに情報を含んでいる。 ジェフリーズ事前はこの問題を解決するために作られた。

現代の実務では、完全な無情報を目指すより 「常識的にありえない値を排除する弱情報事前」を置く方が良いとされる。 例えば、ロジスティック回帰(21 章)の係数 \(\beta\) はオッズ比 \(e^\beta\) に対応するが、「薬の効果のオッズ比が 1000 倍(\(\beta\approx 6.9\))」はありえないので、 \(\beta\sim N(0, 2.5^2)\) 程度(\(|\beta|>5\)、オッズ比 150 倍超が 2 標準偏差の外)の事前を置くのは合理的である。

事前分布の敏感性分析

必ず行うべき作法: 複数の事前分布で計算し、結論が変わらないことを確認する。 変わるなら、その結論はデータではなく事前分布が決めている—— それを明記した上で報告する必要がある。

4. 事後分布からの推論

点推定

事後分布から 1 つの値 \(a\) を答えとして選ぶとき、「外れたときの損失」\(L(\theta,a)\) の事後期待値が最小になる \(a\) を選ぶ、と考えると次の対応になる。

名前定義対応する損失関数
事後平均\(E[\theta\mid D]\)二乗損失 \((\theta-a)^2\)
事後中央値50% 点絶対損失 \(\lvert\theta-a\rvert\)
MAP事後最頻値(事後密度が最大の点)0-1 損失(\(a\) が \(\theta\) にぴったり当たれば 0、外れれば 1)

(二乗損失の場合: \(E[(\theta-a)^2\mid D]\) を \(a\) で微分して 0 とおくと \(-2E[\theta\mid D]+2a=0\)、よって \(a=E[\theta\mid D]\)。絶対損失では \(a\) の左右の確率が等しくなる点=中央値が最小にする。0-1 損失では「当たる確率」を最大にする点、すなわち密度最大の点が選ばれる。)

MAP 推定は正則化と等価である(20 章)。 事前分布を正規にすれば Ridge、ラプラス分布(密度が \(e^{-|\theta|/b}\) に比例する、0 で尖った分布)にすれば Lasso になる。 正則化とは、暗黙のうちにベイズをやっていたのである。

信用区間 (credible interval)

\[ P(\theta\in[a,b]\mid D)=0.95 \]

これは文字通り「\(\theta\) がこの区間にある確率が 95%」である。

種類定義
等裾区間2.5% 点と 97.5% 点
HDI(最高密度区間)密度が最も高い領域。同じ確率で最も狭い

信頼区間 vs 信用区間. 数値はしばしばほぼ同じになる(特に \(n\) が大きく事前が弱いとき)。 だが主張している内容がまったく違う。

11 章で「本当は言えないこと」を、ベイズでは正々堂々と言える。 ただしその代償として事前分布を置いた——そこが取引条件である。

事後予測分布

次の観測 \(\tilde y\) の分布は、\(\theta\) の不確かさを積分して消したもの。

\[ p(\tilde y\mid D)=\int p(\tilde y\mid\theta)p(\theta\mid D)d\theta \]

これがベイズの実務的な強みである. 頻度論の予測区間(11 章・18 章)は「\(\hat\theta\) が正しい」ことを暗に前提にしがちだが、 ベイズの事後予測分布はパラメータの不確かさも自動的に織り込む。

「来月の需要は何個か」という問いに対し、 パラメータ推定の誤差と個体のばらつきの両方を含んだ分布が得られる。

5. MCMC — 共役でないときどうするか

現実の問題では共役事前が使えないことがほとんどである。\(p(D)=\int p(D\mid\theta)p(\theta)d\theta\) が計算できない。

解決策: 積分を計算せず、事後分布からサンプルを取る。 サンプルが十分あれば、平均も分位点も何でも計算できる。

メトロポリス・ヘイスティングス法

  1. 現在の値 \(\theta^{(t)}\) から、提案分布 \(q(\theta^*\mid\theta^{(t)})\)(「今いる場所から次の候補をどう選ぶか」の分布。たとえば \(\theta^{(t)}\) を中心とする正規分布 \(N(\theta^{(t)},s^2)\))に従って候補 \(\theta^*\) を生成
  2. 受容確率を計算:
\[ r = \min\left(1,\ \frac{p(\theta^*\mid D)q(\theta^{(t)}\mid\theta^*)}{p(\theta^{(t)}\mid D)q(\theta^*\mid\theta^{(t)})}\right) \]
  1. 確率 \(r\) で \(\theta^{(t+1)}=\theta^*\)、そうでなければ \(\theta^{(t+1)}=\theta^{(t)}\)

こうして得られる列 \(\theta^{(1)},\theta^{(2)},\dots\) をチェーンと呼ぶ。 \(q\) が比の中に両方向で入っているのは、「行きやすいが戻りにくい」提案の偏りを打ち消して、長く続けたときの滞在分布がちょうど事後分布になるようにするためである。 正規分布のような左右対称な提案(\(q(a\mid b)=q(b\mid a)\))なら \(q\) は約分されて消え、\(r=\min(1, p(\theta^*\mid D)/p(\theta^{(t)}\mid D))\) になる。

正規化定数 \(p(D)\) が比を取ると消える——これが決定的である。 計算できない積分を計算せずに済む。

直感: 事後確率が高い場所へは積極的に移動し、低い場所へは確率的にしか移動しない。 長く歩き回れば、滞在時間の分布が事後分布に一致する。 山の高いところに長くいる、というだけの仕組みである。

主なアルゴリズム

手法特徴
メトロポリス・ヘイスティングス最も基本的。汎用
ギブスサンプリング条件付き分布から順に。共役性が部分的に使えるとき高速
ハミルトニアン MC (HMC)事後密度の対数の勾配(各パラメータで偏微分した傾き)を使って効率的に探索。高次元に強い
NUTSHMC の自動調整版。Stan / PyMC の既定

収束診断(必須)

診断内容
トレースプロット毛虫のようにランダムに見えれば良好。傾向があればダメ
\(\hat R\)(Gelman-Rubin)異なる初期値から走らせた複数チェーンについて、「チェーン間のばらつき」と「チェーン内のばらつき」の比。収束していれば同じ分布を見ているので 1 に近づく。1.01 未満が目安
有効サンプルサイズ隣り合うサンプルは似た値になりやすい(自己相関: 列の中で隣どうしの相関)ので、独立なサンプルに換算した実質的な個数。数百以上が目安
バーンイン最初の数百〜数千回を捨てる

収束していない MCMC の結果は、単なる乱数である. 必ずトレースプロットを見て、\(\hat R\) を確認する。 「計算が終わった = 正しい」ではない。

import pymc as pm
with pm.Model() as model:
    theta = pm.Beta("theta", alpha=2, beta=2)      # 弱情報事前
    y = pm.Binomial("y", n=10, p=theta, observed=8)
    idata = pm.sample(2000, chains=4)              # NUTS
    print(pm.summary(idata))                        # r_hat を確認

6. ベイズファクター — 仮説の比較

比較したい 2 つの仮説(またはモデル)を \(H_0\)、\(H_1\) とする。ベイズファクターは、それぞれの仮説のもとでデータが出る確率(周辺尤度)の比:

\[ \text{BF}_{10}=\frac{p(D\mid H_1)}{p(D\mid H_0)} \]
\[ \underbrace{\frac{P(H_1\mid D)}{P(H_0\mid D)}}_{\text{事後オッズ}}=\text{BF}_{10}\times\underbrace{\frac{P(H_1)}{P(H_0)}}_{\text{事前オッズ}} \]

(05 章のオッズ形式そのもの。)

BF証拠の強さ(Jeffreys、Kass と Raftery による慣例的な区切り。境界に理論的な根拠はない)
1〜3ほとんど意味なし
3〜10中程度
10〜30強い
30〜100非常に強い
100 以上決定的

p 値との決定的な違い: \(H_0\) を支持する証拠を表現できる. p 値は「\(H_0\) を棄却できない」としか言えず、 「\(H_0\) が正しそうだ」を主張できない(12 章)。 BF なら \(\text{BF}_{10}=0.1\) で「\(H_0\) の方が 10 倍もっともらしい」と言える。

また、逐次的に見ても問題が起きない。 17 章のオプショナル・ストッピング問題は、ベイズでは原理的に生じない (データを追加しても事後分布が更新されるだけで、誤り率の概念がない)。 A/B テストの逐次モニタリングでベイズが好まれる理由の 1 つである。

ただしBF は事前分布に敏感である。特に \(H_1\) のもとでのパラメータの事前分布の広さで 値が大きく変わるので、報告時には必ず事前分布を明示する。

7. 階層ベイズモデル

複数のグループがあるとき、グループごとのパラメータにさらに事前分布を置く。

\[ y_{ij}\sim N(\theta_j,\sigma^2), \qquad \theta_j\sim N(\mu,\tau^2), \qquad \mu,\tau\sim \text{(超事前分布)} \]

(\(j\) がグループ、\(i\) がグループ \(j\) の中の観測の番号。\(\theta_j\) がグループごとの平均、\(\mu,\tau\) はグループ平均たちの全体平均とばらつき。)

部分プーリング (partial pooling) という考え方.

例: 野球選手の打率を推定する。10 打数 5 安打の選手を「5 割打者」とは推定しない。 階層モデルは自動的に全体平均(リーグの全選手の打率の平均で、おおよそ 2 割 5 分)に引き戻し、 打数が増えるほど個人のデータを信用するようになる。

これは 20 章の縮小推定 (shrinkage) と同じ働きであり、 統計学における最も有用なアイデアの 1 つである。 小地域推定、教育の学校効果、医療機関の成績評価など、 「グループごとにデータ量が違う」あらゆる場面で使われる。

8. ベイズと頻度論の使い分け

場面推奨
事前情報が明確にあるベイズ
客観性が強く要求される(規制当局)頻度論(ただしベイズも認められつつある)
小標本ベイズ(事前分布が安定化する)
逐次モニタリングしたいベイズ
階層構造・欠測・複雑なモデルベイズ(柔軟)
大標本・単純なモデルどちらでもほぼ同じ結果

敵対する 2 つの流派ではなく、道具箱の 2 つの道具である. 実際、大標本では事前分布の影響が消えて両者はほぼ一致する(ベルンシュタイン・フォン・ミーゼスの定理)。 違いが問題になるのは、データが少ないときと、複雑なモデルを扱うときである。

そして 20 章で見たように、頻度論の正則化はベイズの MAP 推定と等価であった。 実務家は無意識のうちに両方を使っている。

9. まとめ

概念内容
ベイズの定理事後 \(\propto\) 尤度 × 事前
共役事前事後が同じ族。ベータ-二項、ガンマ-ポアソン
事後平均データと事前の加重平均。\(n\) 増で事前の影響が消える
信用区間「\(\theta\) がこの区間にある確率」— 文字通りの意味
MAP正則化と等価(Ridge = 正規事前、Lasso = ラプラス事前)
MCMC事後分布からサンプリング。\(\hat R < 1.01\) を確認
ベイズファクター仮説の比較。\(H_0\) 支持も表現できる
階層モデル部分プーリング。データの少ないグループを引き戻す