負の二項分布の導出・性質とポアソン分布との関係

コインを投げて初めて表が出るまでの回数は幾何分布に従います。では、$r$ 回目の表が出るまでの回数はどんな分布に従うでしょうか?

この疑問への答えが負の二項分布(negative binomial distribution)です。幾何分布が「1回成功するまで」の待ち時間なら、負の二項分布は「$r$ 回成功するまで」の待ち時間の一般化です。

しかし、負の二項分布の重要性は待ち時間の問題に留まりません。統計学において、カウントデータがポアソン分布より「ばらつきが大きい」(過分散)場合に、負の二項分布はポアソンの代替として広く使われています。

  • 生態学: 生物の個体数のカウント(集中分布を示すデータ)
  • 疫学: 感染症の発生件数(過分散を伴うカウントデータ)
  • 自然言語処理: 文書中の単語の出現回数
  • RNA-seq解析: 遺伝子発現量の統計モデリング
  • 保険数理: 保険金請求件数の分布

本記事の内容

  • 負の二項分布の定義と幾何分布からの一般化
  • PMFの導出と期待値・分散の計算
  • ポアソン-ガンマ混合としての解釈
  • 過分散データへの応用
  • Pythonでの実装と可視化

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

負の二項分布とは — 「$r$回成功するまでの試行」

幾何分布からの一般化

幾何分布は「1回の成功」を待つ分布でした。負の二項分布はこれを「$r$ 回の成功」に拡張します。

成功確率 $p$ のベルヌーイ試行を繰り返し、$r$ 回目の成功が起きるまでの総試行回数 $X$ が負の二項分布に従います。$X = k$ となるには

  1. 最初の $k – 1$ 回のうちちょうど $r – 1$ 回成功している
  2. $k$ 回目で $r$ 回目の成功が起きる

PMFの導出

$k$ 回目が $r$ 回目の成功であるためには、最初の $k-1$ 回で $r-1$ 回成功する必要があります。その組み合わせ数は $\binom{k-1}{r-1}$ です。$r-1$ 回の成功と $k-r$ 回の失敗の確率は $p^{r-1} q^{k-r}$ です。最後の $k$ 回目が成功する確率は $p$ です。

これらを掛け合わせると

$$ \begin{equation} P(X = k) = \binom{k-1}{r-1} p^r (1-p)^{k-r}, \quad k = r, r+1, r+2, \ldots \end{equation} $$

$r = 1$ のとき $\binom{k-1}{0} = 1$ なので $P(X = k) = p (1-p)^{k-1}$ となり、幾何分布に一致します。

失敗回数による定義

$Y = X – r$(成功前の失敗の総回数)で定義する慣習もよく使われます。

$$ P(Y = k) = \binom{k + r – 1}{k} p^r (1-p)^k, \quad k = 0, 1, 2, \ldots $$

この形の方が $k = 0$ から始まるため扱いやすく、統計ソフトウェアではこちらが標準的です。本記事では両方の定義を状況に応じて使い分けます。

「負の二項」という名前の由来

名前の「負の二項」は、PMFが負の二項係数を含むことに由来します。

$$ \binom{k+r-1}{k} = \binom{-r}{k}(-1)^k $$

つまり、負の二項展開 $(1-q)^{-r} = \sum_{k=0}^{\infty} \binom{-r}{k}(-q)^k = \sum_{k=0}^{\infty} \binom{k+r-1}{k} q^k$ に現れる係数です。

期待値と分散の導出

幾何分布の和としての期待値

$r$ 回目の成功までの試行回数 $X$ は、$r$ 個の独立な幾何分布の和です。$X_i$ を $(i-1)$ 回目の成功から $i$ 回目の成功までの追加試行回数とすると、$X_i \sim \text{Geom}(p)$ かつ $X = X_1 + X_2 + \cdots + X_r$ です。

$$ E[X] = r \cdot E[X_1] = \frac{r}{p} $$

$$ \text{Var}(X) = r \cdot \text{Var}(X_1) = \frac{r(1-p)}{p^2} $$

失敗回数 $Y = X – r$ で表すと

$$ \begin{equation} E[Y] = \frac{r(1-p)}{p}, \quad \text{Var}(Y) = \frac{r(1-p)}{p^2} \end{equation} $$

過分散の性質

負の二項分布の重要な特徴は過分散(overdispersion)です。

$$ \text{Var}(Y) = \frac{r(1-p)}{p^2} = E[Y] + \frac{(E[Y])^2}{r} > E[Y] $$

ポアソン分布では $\text{Var}(Y) = E[Y]$(平均 = 分散)ですが、負の二項分布では分散が平均より大きくなります。$r$ が大きいほどポアソン分布に近づき($r \to \infty$ でポアソンに収束)、$r$ が小さいほど過分散が強くなります。

この過分散の性質が、実データのモデリングで負の二項分布が重宝される理由です。

ポアソン-ガンマ混合としての解釈

動機

負の二項分布がポアソン分布より適切な場面を理解するには、ポアソン-ガンマ混合(Poisson-Gamma mixture)としての解釈が重要です。

ポアソン分布は「平均 $\lambda$ が全てのユニットで同じ」と仮定します。しかし現実には、各ユニット(個人、地域など)で $\lambda$ が異なることが多いです。この $\lambda$ のばらつきをガンマ分布でモデル化すると、結果として負の二項分布が得られます。

導出

各ユニット $i$ のカウント $Y_i$ がポアソン分布 $Y_i | \lambda_i \sim \text{Poisson}(\lambda_i)$ に従い、$\lambda_i$ がガンマ分布 $\lambda_i \sim \text{Gamma}(r, p/(1-p))$ に従うとします。

$Y_i$ の周辺分布は、$\lambda$ を積分消去して

$$ P(Y = k) = \int_0^{\infty} \frac{e^{-\lambda} \lambda^k}{k!} \cdot \frac{\beta^r}{\Gamma(r)} \lambda^{r-1} e^{-\beta\lambda} d\lambda $$

ここで $\beta = p/(1-p)$ です。この積分を実行すると(ガンマ関数の性質を使って)

$$ P(Y = k) = \binom{k + r – 1}{k} \left(\frac{\beta}{1+\beta}\right)^r \left(\frac{1}{1+\beta}\right)^k $$

$p = \beta/(1+\beta)$ と定義すれば、これは負の二項分布のPMFです。

つまり、ポアソン分布の平均パラメータ自体がランダム(ガンマ分布)である場合、周辺分布は負の二項分布になるのです。

Pythonでの実装と可視化

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats

fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# (a) 異なるr,pでのPMF
ax = axes[0, 0]
params = [(1, 0.3), (3, 0.3), (5, 0.3), (10, 0.3), (5, 0.5)]
k = np.arange(0, 40)

for r, p in params:
    pmf = stats.nbinom.pmf(k, r, p)
    ax.plot(k, pmf, 'o-', markersize=4, linewidth=1.5,
            label=f'r={r}, p={p}')

ax.set_xlabel('k (number of failures)', fontsize=12)
ax.set_ylabel('P(Y = k)', fontsize=12)
ax.set_title('Negative Binomial PMF', fontsize=13)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)

# (b) ポアソン分布との比較(過分散)
ax = axes[1, 0]
np.random.seed(42)
n_data = 1000

# ポアソン-ガンマ混合の生成
r_param = 3
p_param = 0.3
lambdas = np.random.gamma(r_param, (1-p_param)/p_param, n_data)
y_nbinom = np.array([np.random.poisson(l) for l in lambdas])

# 同じ平均のポアソン
mean_val = r_param * (1-p_param) / p_param
y_poisson = np.random.poisson(mean_val, n_data)

k_range = np.arange(0, 30)
ax.hist(y_nbinom, bins=k_range - 0.5, density=True, alpha=0.5,
        color='coral', label=f'NegBin (mean={np.mean(y_nbinom):.1f}, var={np.var(y_nbinom):.1f})')
ax.hist(y_poisson, bins=k_range - 0.5, density=True, alpha=0.5,
        color='steelblue', label=f'Poisson (mean={np.mean(y_poisson):.1f}, var={np.var(y_poisson):.1f})')

# 理論PMF
ax.plot(k_range, stats.nbinom.pmf(k_range, r_param, p_param), 'r-',
        linewidth=2, label='NegBin theory')
ax.plot(k_range, stats.poisson.pmf(k_range, mean_val), 'b--',
        linewidth=2, label='Poisson theory')

ax.set_xlabel('k', fontsize=12)
ax.set_ylabel('Density', fontsize=12)
ax.set_title('Negative Binomial vs Poisson\n(Same mean, different variance)',
             fontsize=12)
ax.legend(fontsize=8)
ax.grid(True, alpha=0.3)

# (c) rが大きくなるとポアソンに収束
ax = axes[0, 1]
mean_fixed = 5  # 平均を5に固定
r_values = [1, 3, 10, 50, 200]
k = np.arange(0, 20)

for r in r_values:
    p = r / (r + mean_fixed)
    pmf = stats.nbinom.pmf(k, r, p)
    label = f'NegBin(r={r})'
    if r == 200:
        label += ' ≈ Poisson'
    ax.plot(k, pmf, 'o-', markersize=4, linewidth=1.5, label=label)

ax.plot(k, stats.poisson.pmf(k, mean_fixed), 'k--', linewidth=2.5,
        label=f'Poisson({mean_fixed})', zorder=0)
ax.set_xlabel('k', fontsize=12)
ax.set_ylabel('P(Y = k)', fontsize=12)
ax.set_title('NegBin → Poisson as r → ∞\n(mean = 5 fixed)', fontsize=12)
ax.legend(fontsize=8)
ax.grid(True, alpha=0.3)

# (d) ポアソン-ガンマ混合の可視化
ax = axes[1, 1]
lambda_range = np.linspace(0, 25, 200)
r_param = 3
p_param = 0.3
beta_param = p_param / (1 - p_param)

# ガンマ分布の密度(λの事前分布)
gamma_pdf = stats.gamma.pdf(lambda_range, r_param, scale=1/beta_param)
ax.plot(lambda_range, gamma_pdf, 'b-', linewidth=2.5,
        label=f'Gamma(r={r_param}, β={beta_param:.2f})')
ax.fill_between(lambda_range, gamma_pdf, alpha=0.15, color='blue')

# いくつかのλでのポアソン分布を重ねる
for lam_val in [2, 5, 10]:
    k_small = np.arange(0, 25)
    poisson_pmf = stats.poisson.pmf(k_small, lam_val)
    ax.bar(k_small, poisson_pmf * 0.3, alpha=0.3, width=0.5,
           label=f'Poisson($\\lambda$={lam_val})')

ax.set_xlabel('$\\lambda$ / k', fontsize=12)
ax.set_ylabel('Density / Probability', fontsize=12)
ax.set_title('Poisson-Gamma Mixture\n$\\lambda \\sim$ Gamma, Y|$\\lambda \\sim$ Poisson',
             fontsize=12)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('negative_binomial.png', dpi=150, bbox_inches='tight')
plt.show()

この可視化から、負の二項分布の重要な性質が読み取れます。

  1. 左上(PMF): $r$ が大きくなるほど分布のピークが右に移動し、形が対称的になります。$p$ が小さいと裾が長くなり、過分散が顕著になります。

  2. 右上(ポアソンへの収束): 平均を5に固定して $r$ を増やすと、負の二項分布が徐々にポアソン分布(黒い破線)に近づきます。$r = 200$ ではほぼ完全に一致しています。

  3. 左下(過分散の比較): 同じ平均のポアソン分布と負の二項分布を比較すると、負の二項分布の方が裾が長く(大きな値が出やすい)、分散が大きいことが明確です。

  4. 右下(ポアソン-ガンマ混合): ガンマ分布で $\lambda$ のばらつきをモデル化し、各 $\lambda$ でのポアソン分布を重ねると、全体として負の二項分布が得られることが視覚的に理解できます。

まとめ

本記事では、負の二項分布の定義・導出・性質・応用を解説しました。

  • 負の二項分布は幾何分布の一般化であり、$r$ 回成功するまでの試行回数の分布
  • 期待値 $E[Y] = r(1-p)/p$、分散 $\text{Var}(Y) = r(1-p)/p^2$
  • 過分散: $\text{Var}(Y) > E[Y]$ であり、ポアソン分布の「平均 = 分散」の制約を緩和する
  • ポアソン-ガンマ混合として解釈でき、ユニット間で平均が異なる場合のカウントデータモデルとして自然
  • $r \to \infty$ でポアソン分布に収束する

次のステップとして、以下の記事も参考にしてください。