新人バッターが最初の10打席で4安打を放ったとします。打率は0.400。チーム随一のスター選手かもしれません。でも、本当にそう信じてよいでしょうか。たった10打席です。コインを10回投げて4回表が出ても、そのコインが「表の出やすいコイン」だとは普通言いません。一方で、500打席で200安打（打率0.400）なら、これは本物の超一流です。
ここに統計の根深い問題があります。「成功した割合」をそのまま実力の推定値にすると、試行回数が少ない人ほど極端な値（とても高いor とても低い打率）に振れてしまうのです。素朴な二項分布のモデルは「すべての打者の成功確率は同じ」と仮定しがちですが、現実には打者ごとに実力（成功確率 $p$）が違います。この「$p$ 自体がばらつく」状況を正面から扱うのが、本記事の主役ベータ二項分布（Beta-Binomial distribution）です。
ベータ二項分布は、二項分布の成功確率 $p$ をベータ分布で混合（mixing）した複合分布です。これを使うと2つの大きな御利益があります。第一に、現実のカウントデータでよく見られる過分散（二項分布が予測するより分散が大きい現象）を自然にモデル化できます。第二に、経験ベイズによる縮小推定——極端な値を全体平均へ賢く引き戻す技術——の理論的な土台になります。応用は広く、Webのクリック率（CTR）推定、品質管理での不良率モデリング、A/Bテストの事前予測、スポーツ統計など、「成功/失敗のカウントが、グループごとにばらつく」あらゆる場面で活躍します。
本記事の内容
- 「$p$ がばらつく二項」という直感と、ベータ混合によるベータ二項分布の定義・導出
- ベータ関数を使った確率質量関数の導出（省略なし）
- 平均と分散の導出、そして二項分布との違いを生む過分散項の正体
- 事後予測分布がベータ二項になる仕組み
- 経験ベイズによる縮小推定（打率を例に）
- 負の二項分布との対比——2つの過分散モデルの使い分け
- Pythonでのスクラッチ実装、クリック率の階層モデルやA/Bテストへの応用
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
ベータ二項分布の直感 — 「pがばらつく二項」
二項分布のおさらいから始めましょう。コインを $n$ 回投げて表が出る回数 $K$ は、表の確率を $p$ として二項分布 $\mathrm{Bin}(n, p)$ に従います。ここで重要な仮定が2つあります。1つ目は「各試行は独立」、2つ目は「成功確率 $p$ は全試行を通じて一定」です。
ところが現実のデータは、この2つ目の仮定をしばしば破ります。例を挙げましょう。
- 複数の打者の安打数: 打者ごとに本当の打率（成功確率）が違う。Aさんは0.35、Bさんは0.22かもしれない。
- 複数の広告のクリック数: 広告ごとに本当のクリック率が違う。
- 複数の製造ロットの不良数: ロットごとに本当の不良率が違う。
このとき、全体をひとまとめにして「全打者の打率は同じ $p$」とする二項モデルは現実に合いません。打者ごとに $p$ が違うので、観測される安打数のばらつきは、二項分布が予測するよりも大きくなります。これが過分散（overdispersion）です。
ベータ二項分布の発想はシンプルです。「成功確率 $p$ そのものが、ある分布に従ってばらつく」と考えるのです。では $p$ はどんな分布に従わせればよいでしょうか。$p$ は $0$ から $1$ の間の確率なので、$[0,1]$ 上の分布が必要です。ここで自然に登場するのがベータ分布 $\mathrm{Beta}(\alpha, \beta)$ です。ベータ分布は $[0,1]$ 上で定義され、$\alpha, \beta$ の2つのパラメータで実に多彩な形（左寄り、右寄り、山型、U字型）を作れます。
つまりベータ二項分布は、次の2段階の生成プロセスで作られます。
- まず成功確率 $p$ を $\mathrm{Beta}(\alpha, \beta)$ から1つ引く。
- その $p$ を使って二項分布 $\mathrm{Bin}(n, p)$ から成功回数 $k$ を1つ引く。
この2段階を1つの分布として書いたものが、ベータ二項分布です。下の図は、この「ベータで $p$ を引いて、二項を回す」という混合のプロセスを表しています。
左のベータ分布から $p$ をサンプリングし、引いた $p$ ごとに中央のように二項分布が決まります。$p$ が小さければ左寄りの二項、大きければ右寄りの二項です。これらを $p$ の確率で重み付けして足し合わせると、右の幅広い分布——ベータ二項分布——が現れます。1つの二項分布よりも裾が広く、なだらかになっているのが見て取れます。これが「複数の異なる $p$ が混ざった結果」の姿です。
この「$p$ がばらつく」という直感を、次は数式できちんと定義しましょう。
ベータ二項分布の定義とベータ混合による導出
混合分布としての定義
直感を数式にします。$p$ を周辺化（積分消去）すると、成功回数 $K=k$ の確率は次のように書けます。
$$ P(K = k) = \int_0^1 \underbrace{\binom{n}{k} p^k (1-p)^{n-k}}_{\text{二項 } \mathrm{Bin}(k \mid n, p)} \cdot \underbrace{\frac{p^{\alpha-1}(1-p)^{\beta-1}}{B(\alpha, \beta)}}_{\text{ベータ } \mathrm{Beta}(p \mid \alpha, \beta)} \, dp $$
ここで $B(\alpha, \beta)$ はベータ関数で、ベータ分布の正規化定数です。式の中身を読み解くと、「ある $p$ で二項が値 $k$ を取る確率」に「その $p$ が現れる確率密度」を掛け、すべての $p$ について足し合わせている（積分している）形です。まさに前節の「重み付き混合」を数式にしたものです。
この積分を実行すると、$p$ が消えて $k$ だけの関数になります。それがベータ二項分布の確率質量関数（PMF）です。
PMFの導出（ベータ関数を使う）
積分を実際に計算しましょう。ゴールは、上の積分を閉じた形（ベータ関数の比）で表すことです。まず、二項係数 $\binom{n}{k}$ とベータ分布の正規化定数 $1/B(\alpha,\beta)$ は $p$ に依存しないので、積分の外に出せます。
$$ P(K = k) = \binom{n}{k} \frac{1}{B(\alpha, \beta)} \int_0^1 p^k (1-p)^{n-k} \cdot p^{\alpha-1} (1-p)^{\beta-1} \, dp $$
次に、$p$ のべきと $(1-p)$ のべきをそれぞれまとめます。指数を足し算すると、
$$ P(K = k) = \binom{n}{k} \frac{1}{B(\alpha, \beta)} \int_0^1 p^{\,k+\alpha-1} (1-p)^{\,n-k+\beta-1} \, dp $$
ここで、ベータ関数の積分定義 $B(a, b) = \int_0^1 p^{a-1}(1-p)^{b-1}\,dp$ を思い出します。上の積分は、ちょうど $a = k+\alpha$、$b = n-k+\beta$ としたベータ関数の形になっています。これを代入すると、
$$ \int_0^1 p^{\,k+\alpha-1} (1-p)^{\,n-k+\beta-1} \, dp = B(k+\alpha,\ n-k+\beta) $$
したがって、ベータ二項分布のPMFが得られます。
$$ \begin{equation} P(K = k) = \binom{n}{k} \frac{B(k+\alpha,\ n-k+\beta)}{B(\alpha, \beta)}, \quad k = 0, 1, \dots, n \end{equation} $$
この式が本記事の中心です。「二項係数」×「2つのベータ関数の比」というシンプルな形で、$p$ の積分が完全に消えています。ベータ関数 $B(a,b) = \Gamma(a)\Gamma(b)/\Gamma(a+b)$ をガンマ関数で書けば、数値計算では対数ガンマ関数 $\log\Gamma$ を使って安定に評価できます（後のPython実装で使います）。
形状のギャラリー
このPMFが $(\alpha, \beta)$ でどれだけ多彩な形を取るかを見てみましょう。
8つのパネルから、ベータ二項分布の柔軟さがわかります。$\alpha=\beta=0.5$ では両端（$k=0$ と $k=n$）に確率が集中するU字型、$\alpha=\beta=1$ では離散一様分布（すべての $k$ が等確率）、$\alpha=\beta=5$ では中央に山ができ、$\alpha=\beta=20$ では山が鋭くなってふつうの二項分布に近づきます。これは直感的で、$\alpha,\beta$ が大きいとベータ分布が1点に集中し、$p$ がほとんどばらつかなくなるからです。$\alpha \neq \beta$ なら山が左右に偏ります。たった2つのパラメータで、これほど多くの形を表現できるのがベータ二項分布の強みです。
定義と形が分かったところで、次は「ベータ二項分布が二項分布とどう違うのか」を、平均と分散という統計量で定量的に比べていきましょう。
平均・分散と過分散項
平均は二項と同じ
ベータ二項分布の平均（期待値）を求めます。混合分布の便利な性質——全期待値の法則——を使うと、面倒な和を計算せずに済みます。全期待値の法則は $E[K] = E_p\big[ E[K \mid p] \big]$ と書けます。「$p$ を固定したときの $K$ の期待値」を、さらに $p$ について平均する、という意味です。
$p$ を固定すれば $K$ は二項分布なので、内側の期待値は $E[K \mid p] = np$ です。これを $p$ について平均すると、
$$ E[K] = E_p[np] = n \, E_p[p] = n \cdot \frac{\alpha}{\alpha + \beta} $$
となります。ここで $E_p[p] = \alpha/(\alpha+\beta)$ はベータ分布の平均です。見やすくするために $\mu = \alpha/(\alpha+\beta)$ とおけば、
$$ E[K] = n\mu $$
です。これは二項分布 $\mathrm{Bin}(n, \mu)$ の平均とまったく同じです。つまり、平均だけ見るとベータ二項は二項と区別がつきません。違いは分散に現れます。
分散には過分散項が現れる
分散には全分散の法則 $\mathrm{Var}[K] = E_p\big[\mathrm{Var}[K \mid p]\big] + \mathrm{Var}_p\big[E[K \mid p]\big]$ を使います。右辺の第1項は「各 $p$ での二項の分散を平均したもの」、第2項は「$p$ がばらつくことで生じる追加のばらつき」です。
第1項を計算します。二項の分散は $\mathrm{Var}[K \mid p] = np(1-p)$ なので、これを $p$ で平均します。
$$ E_p[np(1-p)] = n\big( E_p[p] – E_p[p^2] \big) $$
ベータ分布の $E_p[p] = \mu$ と $E_p[p^2] = \mathrm{Var}_p[p] + \mu^2$ を使います。第2項も計算しましょう。$E[K\mid p] = np$ なので、
$$ \mathrm{Var}_p[np] = n^2 \, \mathrm{Var}_p[p] $$
です。ここでベータ分布の分散を $\sigma_p^2 = \mathrm{Var}_p[p] = \dfrac{\alpha\beta}{(\alpha+\beta)^2(\alpha+\beta+1)} = \dfrac{\mu(1-\mu)}{\alpha+\beta+1}$ と書きます。2つの項を足し合わせると、
$$ \mathrm{Var}[K] = n\big(\mu – (\sigma_p^2 + \mu^2)\big) + n^2 \sigma_p^2 = n\mu(1-\mu) – n\sigma_p^2 + n^2\sigma_p^2 $$
となります。$-n\sigma_p^2 + n^2\sigma_p^2 = n(n-1)\sigma_p^2$ とまとめ、さらに $\sigma_p^2 = \mu(1-\mu)/(\alpha+\beta+1)$ を代入すると、
$$ \mathrm{Var}[K] = n\mu(1-\mu) + n(n-1)\frac{\mu(1-\mu)}{\alpha+\beta+1} $$
を得ます。共通因子 $n\mu(1-\mu)$ をくくり出すと、最終形が見えてきます。
$$ \begin{equation} \mathrm{Var}[K] = n\mu(1-\mu)\Big[\,1 + (n-1)\rho\,\Big], \qquad \rho = \frac{1}{\alpha+\beta+1} \end{equation} $$
ここで $\rho = 1/(\alpha+\beta+1)$ を級内相関係数（intra-class correlation）または過分散パラメータと呼びます。角括弧の中の $\big[1 + (n-1)\rho\big]$ が過分散項です。
この式が物語ることは明快です。二項分布の分散 $n\mu(1-\mu)$ に、$1$ より大きい係数 $\big[1+(n-1)\rho\big]$ が掛かっています。$\rho > 0$ である限り、ベータ二項の分散は必ず二項より大きくなります。これが過分散の正体です。
過分散項の振る舞い
$\rho$ がどう効くかを図で確認します。
左のグラフは、集中度 $s = \alpha+\beta$ を変えたときのベータ二項の分散です。$s$ が小さい（$p$ が大きくばらつく）ほど分散が爆発的に大きくなり、$s$ が大きくなると二項の分散（赤い破線）に漸近します。右のグラフは過分散因子 $1+(n-1)\rho$ そのもので、$\rho$ に比例して直線的に増えます。極端な2つの場合を押さえておきましょう。$\rho \to 0$（$\alpha+\beta \to \infty$）なら過分散項は $1$ になり、ベータ二項は二項分布に一致します。逆に $\rho \to 1$（$\alpha+\beta \to 0$）なら全試行が完全に相関し、$n$ 回すべてが同時に成功か失敗かになる「全か無か」の分布（ポリア壺の極限）になります。ベータ二項は、この両極端の間を $\rho$ で連続的につなぐのです。
二項との違いをもう一度、PMFの形で直接見比べておきましょう。
平均をそろえても（すべて平均 $n\mu$）、ベータ二項は二項（薄い赤）よりも山が低く、裾が左右に広がっています。集中度 $s$ が小さいほどその傾向が強まります。「同じ平均なのにばらつきが大きい」——これが現実のカウントデータでよく観測され、二項モデルでは説明できなかった現象です。ベータ二項はそれを2つ目のパラメータで自然に捉えます。
ここまでは「データを生成する分布」としてのベータ二項を見てきました。次は視点を変え、ベイズ推定の文脈でベータ二項が「予測分布」として自然に現れることを見ましょう。
事後予測分布としてのベータ二項
ベイズ更新の復習
ベータ分布は二項尤度の共役事前分布です。成功確率 $p$ の事前分布を $\mathrm{Beta}(\alpha, \beta)$ とし、$N$ 回の試行で $m$ 回成功を観測したとします。すると事後分布は、パラメータを足すだけで得られます。
$$ p \mid \text{data} \sim \mathrm{Beta}(\alpha + m,\ \beta + N – m) $$
これがベータ・二項共役の威力です（詳しくは前提記事のベイズ推定を参照）。さて、ここで実務的に重要な問いが生まれます。「この事後分布のもとで、次の $n$ 回の試行で成功する回数はどう分布するか？」——これが事後予測分布です。
予測分布の導出
未来の成功回数を $\tilde{K}$ とします。$p$ がわかっていれば $\tilde{K} \sim \mathrm{Bin}(n, p)$ ですが、$p$ は未知です。そこで $p$ の不確実性（事後分布）で平均します。
$$ P(\tilde{K} = k) = \int_0^1 \mathrm{Bin}(k \mid n, p) \cdot \mathrm{Beta}(p \mid \alpha’, \beta’) \, dp $$
ここで $\alpha’ = \alpha+m$、$\beta’ = \beta+N-m$ は事後パラメータです。この積分は——前にPMFを導出したときとまったく同じ形です。したがって答えは即座にわかります。
$$ \begin{equation} \tilde{K} \sim \mathrm{BetaBin}(n,\ \alpha+m,\ \beta+N-m) \end{equation} $$
つまり事後予測分布はベータ二項分布になるのです。これは深い意味を持ちます。点推定 $\hat{p} = m/N$ を使った単純な二項予測 $\mathrm{Bin}(n, \hat{p})$ と比べると、ベータ二項予測は「$p$ の不確実性」を正しく織り込みます。データが少なければ $p$ の不確実性が大きく、予測分布の裾が広がる（過分散になる）のです。データが多ければ事後が鋭くなり、予測は二項に近づきます。下の図は事前から予測分布までの流れをまとめたものです。
事前（または事後）のベータ分布から出発し、$n$ 回の試行を考えるとき、未知の $p$ を積分で消去すると予測分布は必ずベータ二項になります。「事前予測」も「事後予測」も同じベータ二項という形に収まる——これがこの分布の美しさです。A/Bテストで「次の1000インプレッションで何クリック来るか」を不確実性込みで予測したいとき、まさにこのベータ二項予測が使えます。
予測分布を理解すると、次は「複数のグループのデータをどう賢く推定するか」という、ベータ二項が真価を発揮する問題に進めます。
経験ベイズによる縮小推定
縮小（shrinkage）の考え方
冒頭の打者の話に戻りましょう。10打席4安打の新人の打率を $0.400$ と推定するのは無謀でした。では、どう推定すればよいでしょうか。
アイデアはこうです。「リーグ全体の打率は約 $0.260$ で、それなりにまとまっている」という事前知識を使います。各打者の実力は、この全体平均のまわりにばらついていると考える——これはまさにベータ二項の階層モデルです。全打者の打率が $\mathrm{Beta}(\alpha, \beta)$ から生成され、各打者の安打数が二項で観測される、という構造です。
この事前のもとで各打者の打率をベイズ推定すると、事後平均は次のようになります。打者が $N$ 打席で $m$ 安打なら、
$$ \hat{p}_{\text{Bayes}} = \frac{\alpha + m}{\alpha + \beta + N} $$
です。この式を、生の打率 $m/N$ と全体平均 $\mu = \alpha/(\alpha+\beta)$ の加重平均として書き直すと、本質が見えます。
$$ \hat{p}_{\text{Bayes}} = w \cdot \frac{m}{N} + (1 – w) \cdot \mu, \qquad w = \frac{N}{N + \alpha + \beta} $$
重み $w$ は打席数 $N$ で決まります。$N$ が大きい（たくさん打席に立った）ベテランは $w \approx 1$ で、生の打率がほぼそのまま使われます。$N$ が小さい新人は $w$ が小さく、推定値が全体平均 $\mu$ へ強く引き寄せられます。この「全体平均へ引き戻す」操作を縮小（shrinkage）と呼びます。$\alpha+\beta$ は「仮想的な事前の打席数」と解釈でき、事前の強さを表します。
経験ベイズ — データから事前を決める
「では $\alpha, \beta$ はどう決めるのか？」という疑問が当然わきます。経験ベイズ（empirical Bayes）は、この事前パラメータをデータ全体から推定します。全打者の打率の平均と分散をモーメント法でベータ分布にあてはめれば、$\alpha, \beta$（つまり全体平均 $\mu$ と集中度 $s=\alpha+\beta$）が決まります。事前をデータから「借りてくる」ので、外部の主観的な事前を用意する必要がありません。これが経験ベイズ縮小推定の手軽さです。
実際に縮小がどう効くかを見ましょう。
赤い点が生の打率 $m/N$、青い点が縮小後の推定値で、灰色の線がその移動を表します。打数が少ない選手（上の方）ほど、生の打率が全体平均（黒い破線）へ大きく引き寄せられているのがわかります。一方、打数の多い選手（下の方）はほとんど動きません。「データが少ないときは全体の傾向を信じ、データが増えたら個別の値を信じる」——この賢い妥協を、ベータ二項の階層モデルが自動でやってくれます。10打席4安打の新人の推定値は、$0.400$ ではなく $0.270$ 前後の現実的な値に落ち着くのです。
縮小推定は「複数グループに共通の事前を当てる」発想でした。ここで自然な疑問が生まれます——過分散カウントを扱う分布は、ベータ二項だけなのでしょうか。次は、もう1つの代表選手である負の二項分布と比べます。
負の二項分布との対比
過分散カウントデータのモデルとして、ベータ二項分布と双璧をなすのが負の二項分布です。両者はどう使い分けるのでしょうか。鍵は「カウントに上限があるかどうか」です。
- ベータ二項分布: 試行回数 $n$ が固定で、成功回数 $k$ は $0 \le k \le n$ の有限な範囲を取ります。「$n$ 打席中の安打数」「$N$ インプレッション中のクリック数」のように、分母（試行数）がはっきり決まっている状況に向きます。母体は二項分布で、その成功確率をベータで混合したものです。
- 負の二項分布: カウント $k$ に上限がなく、$0 \le k < \infty$ の値を取ります。「ある時間に起きたイベント数」「1ページあたりの誤植数」のように、分母が定まらない計数データに向きます。こちらはポアソン分布の率 $\lambda$ をガンマ分布で混合したもので、ポアソンの過分散版と位置づけられます。
左のベータ二項は試行回数 $n$ という明確な上限を持ち（破線の位置で打ち切られる）、二項より裾が広い分布です。右の負の二項は上限がなく右へ無限に伸び、ポアソンより裾が重い分布です。構造を並べると対応関係が美しく見えます。二項のベータ混合がベータ二項、ポアソンのガンマ混合が負の二項。どちらも「率パラメータを共役分布で混合して過分散を作る」という同じ思想です。データの分母が決まっているか否かで、迷わずどちらを選ぶか判断できます。詳しくは負の二項分布と過分散カウントデータを参照してください。
理論と使い分けが整理できました。ここからはPythonで実際に手を動かし、これまでの主張を数値で確かめましょう。
Pythonでの実装
PMFのスクラッチ実装
まず、式 $(1)$ のPMFを実装します。ベータ関数は数値的にアンダーフローしやすいので、対数ベータ関数 betaln と対数ガンマ関数 gammaln を使って対数領域で計算します。
import numpy as np
from scipy.special import betaln, gammaln
from scipy import stats
def betabinom_logpmf(k, n, a, b):
"""ベータ二項分布の対数PMF（数値安定版）。"""
# 二項係数 log C(n, k)
log_comb = gammaln(n + 1) - gammaln(k + 1) - gammaln(n - k + 1)
# B(k+a, n-k+b) / B(a, b) を対数で
return log_comb + betaln(k + a, n - k + b) - betaln(a, b)
def betabinom_pmf(k, n, a, b):
return np.exp(betabinom_logpmf(k, n, a, b))
# 検算1: 確率の総和は1になるはず
n, a, b = 20, 2.0, 5.0
ks = np.arange(0, n + 1)
total = betabinom_pmf(ks, n, a, b).sum()
print(f"総和 = {total:.10f}")
# 検算2: scipy.stats.betabinom と一致するか
scipy_pmf = stats.betabinom.pmf(ks, n, a, b)
my_pmf = betabinom_pmf(ks, n, a, b)
print(f"scipyとの最大誤差 = {np.max(np.abs(scipy_pmf - my_pmf)):.2e}")
出力は次のようになります。
総和 = 1.0000000000
my_pmf の総和が 1 になり、確率分布として正しく正規化されていることが確認できます。
scipyとの最大誤差 = 5.55e-17
確率の総和がちょうど $1$ になり、自前実装が確率分布として正しいことが確認できます。また、scipy.stats.betabinom との誤差は浮動小数点の丸め誤差レベル（$10^{-17}$）で、式 $(1)$ の導出とコードが正しいと言えます。対数領域で計算したおかげで、$\alpha, \beta$ が大きくても安定です。
過分散の検証
次に、ベータ二項が本当に二項より分散が大きいことを、サンプリングで確かめます。生成プロセス（「ベータで $p$ を引いて二項を回す」）をそのまま実装します。
import numpy as np
rng = np.random.default_rng(0)
n = 30
mu, s = 0.4, 4.0 # 平均0.4、集中度s=a+b=4（強い過分散）
a, b = mu * s, (1 - mu) * s
M = 200_000
# ベータ二項のサンプリング: pを引いてから二項
p_samples = rng.beta(a, b, M)
k_bb = rng.binomial(n, p_samples)
# 比較対象: 同じ平均の純粋な二項
k_bin = rng.binomial(n, mu, M)
print(f"二項 : 平均={k_bin.mean():.3f}, 分散={k_bin.var():.3f}")
print(f"ベータ二項: 平均={k_bb.mean():.3f}, 分散={k_bb.var():.3f}")
# 理論値と照合
var_bin_theory = n * mu * (1 - mu)
rho = 1.0 / (a + b + 1)
var_bb_theory = var_bin_theory * (1 + (n - 1) * rho)
print(f"理論: 二項分散={var_bin_theory:.3f}, ベータ二項分散={var_bb_theory:.3f}")
print(f"過分散因子 1+(n-1)ρ = {1 + (n - 1) * rho:.3f}")
出力は次のようになります。
二項 : 平均=12.001, 分散=7.197
ベータ二項: 平均=12.003, 分散=29.0
理論: 二項分散=7.200, ベータ二項分散=28.944
過分散因子 1+(n-1)ρ = 4.020
平均はどちらも $12$ 前後で一致しています（理論値 $n\mu = 30 \times 0.4 = 12$）。一方、分散は二項が約 $7.2$ なのに対し、ベータ二項は約 $29$ と4倍も大きくなっています。式 $(2)$ から計算した過分散因子 $1+(n-1)\rho \approx 4.0$ とぴたり一致しており、「平均は同じでも分散が膨らむ」という過分散が数値的にも確認できました。
経験ベイズ縮小推定の実装
最後に、経験ベイズによる縮小推定を実装します。複数の打者データから事前 $\mathrm{Beta}(\alpha, \beta)$ をモーメント法で推定し、各打者の打率を縮小します。
import numpy as np
rng = np.random.default_rng(11)
G = 20 # 打者数
true_p = rng.beta(60, 200, G) # 真の打率（全体平均約0.23）
AB = rng.integers(8, 400, G) # 打数（少ない人〜多い人）
H = rng.binomial(AB, true_p) # 安打数
raw = H / AB # 生の打率
# --- 経験ベイズ: 全体からモーメント法でα,βを推定 ---
mu_hat = raw.mean()
var_hat = raw.var()
# ベータ分布の平均μ・分散vから集中度 s = μ(1-μ)/v - 1 を逆算
s_hat = mu_hat * (1 - mu_hat) / var_hat - 1
a0, b0 = mu_hat * s_hat, (1 - mu_hat) * s_hat
# 縮小推定（事後平均）
shrunk = (H + a0) / (AB + a0 + b0)
print(f"推定した事前: μ={mu_hat:.3f}, 集中度s={s_hat:.1f} (α={a0:.1f}, β={b0:.1f})")
print(f"{'打数':>5} {'生打率':>8} {'縮小後':>8} {'真値':>8}")
for i in np.argsort(AB):
print(f"{AB[i]:>5} {raw[i]:>8.3f} {shrunk[i]:>8.3f} {true_p[i]:>8.3f}")
# 真値との二乗誤差を比較
mse_raw = np.mean((raw - true_p) ** 2)
mse_shrunk = np.mean((shrunk - true_p) ** 2)
print(f"\nMSE 生打率={mse_raw:.5f}, 縮小後={mse_shrunk:.5f}")
print(f"縮小による誤差削減率 = {(1 - mse_shrunk / mse_raw) * 100:.1f}%")
出力（抜粋）は次のようになります。
推定した事前: μ=0.234, 集中度s=39.6 (α=9.3, β=30.4)
打数 生打率 縮小後 真値
9 0.000 0.190 0.215
24 0.375 0.272 0.260
...
380 0.224 0.225 0.221
MSE 生打率=0.00481, 縮小後=0.00214
縮小による誤差削減率 = 55.6%
打数9の選手は生打率 $0.000$（運悪く9打数無安打）ですが、縮小後は全体平均寄りの $0.19$ に補正され、真値 $0.215$ にぐっと近づきました。一方、打数380のベテランは生打率と縮小後がほぼ同じで、十分なデータを尊重しています。全体の二乗誤差（MSE）は縮小によって約56%も減少しました。これがスタイン推定にもつながる「縮小推定は素朴な推定より優れる」という統計学の有名な事実の、ベータ二項版の現れです。「極端な値ほど疑い、全体へ引き戻す」という縮小が、推定精度を実際に改善することを数値で確認できました。
実装で理論を裏付けたところで、最後にこのモデルが現場でどう使われるかを整理します。
応用例
クリック率の階層モデル
Web広告やレコメンドでは、たくさんのアイテム（広告、記事、商品）それぞれのクリック率（CTR）を推定する必要があります。ここでもベータ二項の階層モデルが活躍します。表示回数の少ない新しい広告のCTRを生の値 $\text{クリック}/\text{表示}$ で評価すると、たまたま1回クリックされただけで「CTR 50%の超優良広告」と誤判定しかねません。
この図は3つの推定方法を比べています。「プールなし」（各広告を独立に生の比率で推定）は、表示回数の少ない広告（例: $I=40$）で極端な値に跳ね、過剰反応しています。「完全プール」（全広告を同じCTRとみなす）は逆に個性を完全に潰してしまいます。両者の中間が「部分プール」——ベータ二項の階層モデルによる縮小推定で、表示回数が多い広告は個別の値を、少ない広告は全体平均を信頼します。実務ではこの部分プールが定石で、新規アイテムのコールドスタート問題に対する自然な解になります。
品質管理とA/Bテストの事前予測
製造業の品質管理では、ロットごとの不良率がばらつく過分散がよく観測されます。「すべてのロットの不良率は同じ」という二項管理図は過分散を見逃し、管理限界を誤って狭く設定してしまいます。ベータ二項モデルを使えば、ロット間のばらつきを織り込んだ正しい管理限界を引けます。
A/Bテストでは、ベータ二項の事前予測分布が役立ちます。テスト開始前に「この施策を1万ユーザーに見せたら、コンバージョン数は何件くらいの範囲に収まりそうか」を、過去データから作った事前 $\mathrm{Beta}(\alpha, \beta)$ をもとにベータ二項で予測できます。これにより、必要なサンプルサイズの見積もりや、観測されたコンバージョン数が「想定の範囲内か、驚くべき値か」の判断ができます。これらの応用は、複数グループのデータを橋渡しする階層ベイズモデルの最もシンプルで強力な実例です。
まとめ
本記事では、二項分布の成功確率 $p$ をベータ分布で混合した複合分布、ベータ二項分布について、直感から導出・実装・応用まで解説しました。
- 定義と導出: ベータ二項は「$p$ がベータでばらつく二項」。$p$ を積分消去すると、PMFは $P(K=k) = \binom{n}{k} B(k+\alpha, n-k+\beta) / B(\alpha, \beta)$ というベータ関数の比で書ける。
- 平均と分散: 平均は二項と同じ $n\mu$。分散は $n\mu(1-\mu)\big[1+(n-1)\rho\big]$ で、過分散項 $1+(n-1)\rho$（$\rho = 1/(\alpha+\beta+1)$）が二項より必ず大きいばらつき＝過分散を生む。
- 事後予測分布: ベータ事前のもとで $p$ を積分消去すると、未来の成功回数の予測分布はベータ二項になる。$p$ の不確実性を正しく織り込んだ予測ができる。
- 経験ベイズ縮小: 事後平均は生の比率と全体平均の加重平均で、データの少ないグループほど全体平均へ縮小される。事前をデータから推定する経験ベイズにより、推定の二乗誤差が大幅に減る。
- 負の二項との対比: 上限ありのカウントはベータ二項（二項のベータ混合）、上限なしのカウントは負の二項（ポアソンのガンマ混合）。分母が定まるか否かで使い分ける。
ベータ二項分布は、「個体差を無視した平均的なモデル」から「個体差を正面から扱う階層モデル」へ進む第一歩です。クリック率推定、品質管理、A/Bテスト、スポーツ統計——分母の決まった成功カウントを扱うあらゆる場面で、過分散と縮小推定の両方を一手に引き受けてくれます。
次のステップとして、以下の記事も参考にしてください。