電球はいつ切れるでしょうか。窓口にいる客は次に何分後に来るでしょうか。原子核はいつ崩壊するでしょうか。これらに共通するのは、「いつ起こるかは決まっていないが、平均的なペースは分かっている」という状況です。ある事象がランダムに、しかし一定の平均ペースで起こるとき、「次にそれが起こるまでの待ち時間」を記述する分布が指数分布(exponential distribution)です。
指数分布は、単に「待ち時間の分布」というだけでなく、確率論・統計学・工学の至るところに顔を出す基礎部品です。具体的には次のような場面で活躍します。
- 信頼性工学: 電子部品の寿命や故障までの時間をモデル化する。「故障率が一定」という仮定のもとでの寿命分布はまさに指数分布です。
- 待ち行列理論(M/M/1など): コールセンターやサーバへのリクエスト到着間隔・サービス時間を指数分布で表し、待ち時間や混雑を解析する。
- 生存時間解析: 患者の生存時間や機械の稼働時間を扱う。指数分布は最も単純な生存モデルとして出発点になります。
ちなみに、いわゆる「指数型分布族(exponential family)」と「指数分布」は名前は似ていますが概念が異なります。指数分布は指数型分布族の一員ではありますが、本記事では指数型分布族の一般論には立ち入りません。あくまで「待ち時間の分布」としての指数分布を、定義から導出、応用まで丁寧に追いかけます。
本記事の内容
- 指数分布の直感的な理解と確率密度関数の定義
- 平均 $1/\lambda$・分散 $1/\lambda^2$ の導出（省略なし）
- 無記憶性、ポアソン過程・ガンマ分布との関係、一定ハザード
- 最尤推定とPythonでの実装・可視化
- 信頼性・待ち行列・生存解析への応用
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
特にポアソン分布は指数分布と表裏一体の関係にあります。「単位時間あたり平均 $\lambda$ 回起こる事象が、ある時間内に何回起こるか」を数えるのがポアソン分布、「次の1回が起こるまで何秒かかるか」を測るのが指数分布です。両者は同じ現象を「回数」と「時間」という別の側面から見ているにすぎません。
指数分布の直感 — なぜ「指数」なのか
まずは数式の前に、なぜ待ち時間が指数関数で表されるのかをイメージしておきましょう。
放射性原子核を例にとります。原子核は「年をとる」ことがありません。1秒間に崩壊する確率はいつでも同じで、過去にどれだけ存在し続けたかには一切依存しません。いま手元に大量の原子核があり、ごく短い時間 $\Delta t$ の間に、各核が確率 $\lambda \Delta t$ で崩壊するとします。
すると、時刻 $t$ までに崩壊していない確率を $S(t)$ と書くと、次の $\Delta t$ の間にも生き残る確率は $1 – \lambda \Delta t$ ですから、
$$ \begin{equation} S(t + \Delta t) = S(t)\,(1 – \lambda \Delta t) \end{equation} $$
が成り立ちます。これを整理すると $S(t+\Delta t) – S(t) = -\lambda \Delta t\, S(t)$ となり、$\Delta t \to 0$ の極限で微分方程式
$$ \begin{equation} \frac{dS}{dt} = -\lambda S(t) \end{equation} $$
になります。「減り方が現在の量に比例する」というこの式の解は、$S(0)=1$ のもとで $S(t) = e^{-\lambda t}$ です。指数関数が現れるのは、減衰率が一定という性質から自然に導かれた結果なのです。
この $S(t) = e^{-\lambda t}$ は「時刻 $t$ までイベントが起きていない確率」、すなわち後で出てくる生存関数です。ここから待ち時間そのものの密度を取り出すには、$S(t)$ を時間で微分すればよく、その符号を反転した $\lambda e^{-\lambda t}$ が確率密度になります。これがまさに指数分布の正体です。次の節では、この密度関数を正式な定義として書き下します。
指数分布の定義
直感で見たとおり、指数分布の確率密度関数(PDF)は減衰率 $\lambda$ をもつ指数関数です。
単位時間あたりにイベントが平均 $\lambda$ 回起こるとき、次にイベントが起こるまでの待ち時間 $x \ge 0$ は、次の密度に従います。
$$ \begin{equation} f(x) = \lambda e^{- \lambda x} \qquad (x \ge 0) \end{equation} $$
ここで $\lambda > 0$ はレート(率)パラメータと呼ばれ、「単位時間あたりのイベント発生回数」を表します。$\lambda$ が大きいほどイベントは頻繁に起こり、待ち時間は短くなります。$x < 0$ では $f(x) = 0$ です（待ち時間は負になりえません）。
一方で、同じ分布を次のように書くこともあります。
$$ \begin{equation} f(x) = \frac{1}{\theta}e^{ – \frac{x}{\theta}} \end{equation} $$
この表現の $\theta$ は「イベントが発生するまでの平均時間」を表します。式 (3) と式 (4) は $\theta = 1/\lambda$ の関係で結ばれた、まったく同じ分布です。レートで見るか、平均待ち時間（スケール)で見るかの違いにすぎません。SciPyの scipy.stats.expon は後者の scale パラメータ $\theta = 1/\lambda$ を引数にとるので、実装ではこの対応を常に意識しましょう。
念のため、これが確率密度として正しい（全区間で積分すると1になる）ことを確認しておきます。
$$ \begin{equation} \int_0^{\infty} \lambda e^{-\lambda x}\, dx = \left[ -e^{-\lambda x} \right]_0^{\infty} = 0 – (-1) = 1 \end{equation} $$
たしかに全確率が1になりました。今考えている対象が「単位時間あたりの平均発生回数」ならレート $\lambda$ を使い、「発生までの平均時間」を直接扱いたいなら $\theta$ を使う、と覚えておけば混乱しません。
この辺りは、実際にグラフを見ると一気にイメージできます。次の節で $\lambda$ を変えながら形を眺めてみましょう。
指数分布を可視化する
$\lambda$ を変えると密度の形がどう変わるかを見てみます。
import numpy as np
from scipy.stats import expon
import matplotlib.pyplot as plt
x = np.linspace(0, 6, 500)
fig, ax = plt.subplots(figsize=(10, 5))
for lam in [0.5, 1.0, 2.0, 3.0]:
# expon は scale = 1/lam を引数にとる
y = expon.pdf(x, scale=1 / lam)
ax.plot(x, y, lw=2, label=f"λ = {lam}")
ax.set_xlabel("待ち時間 x")
ax.set_ylabel("密度 f(x)")
ax.legend()
ax.grid(alpha=0.3)
plt.show()
このグラフから、指数分布の本質的な特徴が読み取れます。すべての曲線は $x=0$ で最大値 $\lambda$ をとり、そこから単調に減少します。$\lambda$ が大きいほど $x=0$ での値は高くなり（発生がすぐ起こりやすい）、しかも急速に減衰します。つまり $\lambda$ が大きいほど「短い待ち時間」に確率が集中するのです。逆に $\lambda$ が小さいと裾が長く伸び、長い待ち時間も起こりやすくなります。
注意したいのは、密度の値 $\lambda$ は $\lambda > 1$ のとき1を超えることです。確率密度は確率そのものではなく「単位幅あたりの確率の濃さ」なので、1を超えても問題ありません（面積が1になっていればよい）。
形が分かったところで、次は「では平均的にどれくらい待つのか」という統計量を計算してみましょう。
指数分布の期待値と分散
指数分布の期待値（平均）と分散は、レート $\lambda$ だけでシンプルに書けます。
期待値
$$ \begin{equation} E[X] = \frac{1}{\lambda} \end{equation} $$
分散
$$ \begin{equation} \mathrm{Var}[X] = \frac{1}{\lambda^2} \end{equation} $$
期待値が $1/\lambda$ なのは直感的に納得できます。1秒あたり平均 $\lambda$ 回イベントが起こるなら、1回あたりの平均間隔は $1/\lambda$ 秒です。たとえば1分間に平均3人の客が来るなら、客が来る平均間隔は $1/3$ 分 $= 20$ 秒、という具合です。
興味深いのは標準偏差 $\sqrt{\mathrm{Var}[X]} = 1/\lambda$ が平均とぴったり等しいことです。つまり指数分布は「平均と同じだけのばらつき」をもつ、非常に散らばりの大きい分布です。変動係数（標準偏差÷平均）は常に1になります。
上の図は $\lambda=1$ の場合に、平均 $1/\lambda$ の位置（赤破線）と、平均 $\pm$ 標準偏差の帯（緑）を描いたものです。平均がちょうど密度の重心になっており、ばらつきの幅（標準偏差）が平均と同じスケールであることが視覚的に分かります。平均より左側は密度が高く、右側に長い裾を引いているため、平均は中央値より右に位置します。
これらの値がどこから来るのかを、次の節で省略せずに導出します。
期待値 1/λ の導出
期待値の定義は $E[X] = \int_0^{\infty} x\, f(x)\, dx$ です。指数分布を代入すると、
$$ \begin{equation} E[X] = \int_0^{\infty} x\, \lambda e^{-\lambda x}\, dx \end{equation} $$
を計算すればよいことになります。この積分は部分積分で処理します。$u = x$, $dv = \lambda e^{-\lambda x}\, dx$ と置くと、$du = dx$, $v = -e^{-\lambda x}$ です。部分積分の公式 $\int u\, dv = uv – \int v\, du$ を使うと、
$$ \begin{equation} E[X] = \left[ -x e^{-\lambda x} \right]_0^{\infty} + \int_0^{\infty} e^{-\lambda x}\, dx \end{equation} $$
となります。第1項を評価します。$x \to \infty$ では指数関数 $e^{-\lambda x}$ が $x$ より速く0に向かうので $-x e^{-\lambda x} \to 0$、$x=0$ では $-x e^{-\lambda x} = 0$ です。したがって第1項は0になります。残る第2項を計算すると、
$$ \begin{equation} \int_0^{\infty} e^{-\lambda x}\, dx = \left[ -\frac{1}{\lambda} e^{-\lambda x} \right]_0^{\infty} = \frac{1}{\lambda} \end{equation} $$
です。以上より $E[X] = 1/\lambda$ が得られました。「1単位時間に平均 $\lambda$ 回起こるなら、1回あたりの間隔は $1/\lambda$」という直感が、積分でもきちんと確認できました。
分散を求めるには2次のモーメント $E[X^2]$ が必要なので、続けて計算します。
分散 1/λ² の導出
分散は $\mathrm{Var}[X] = E[X^2] – (E[X])^2$ で計算できます。すでに $E[X] = 1/\lambda$ は分かっているので、$E[X^2]$ を求めます。定義に代入すると、
$$ \begin{equation} E[X^2] = \int_0^{\infty} x^2\, \lambda e^{-\lambda x}\, dx \end{equation} $$
です。これも部分積分します。$u = x^2$, $dv = \lambda e^{-\lambda x}\, dx$ と置くと $du = 2x\, dx$, $v = -e^{-\lambda x}$ なので、
$$ \begin{equation} E[X^2] = \left[ -x^2 e^{-\lambda x} \right]_0^{\infty} + \int_0^{\infty} 2x\, e^{-\lambda x}\, dx \end{equation} $$
となります。第1項は先ほどと同じ理由で0になります（$x \to \infty$ で $x^2 e^{-\lambda x} \to 0$、$x=0$ で0）。第2項に注目すると、$\int_0^{\infty} 2x\, e^{-\lambda x}\, dx = \frac{2}{\lambda}\int_0^{\infty} x\, \lambda e^{-\lambda x}\, dx$ と書き直せます。ここで右側の積分は、まさに先ほど求めた期待値 $E[X] = 1/\lambda$ そのものです。これを代入すると、
$$ \begin{equation} E[X^2] = \frac{2}{\lambda} \cdot \frac{1}{\lambda} = \frac{2}{\lambda^2} \end{equation} $$
が得られます。あとは分散の公式に代入するだけです。
$$ \begin{equation} \mathrm{Var}[X] = E[X^2] – (E[X])^2 = \frac{2}{\lambda^2} – \left(\frac{1}{\lambda}\right)^2 = \frac{2}{\lambda^2} – \frac{1}{\lambda^2} = \frac{1}{\lambda^2} \end{equation} $$
これで分散 $1/\lambda^2$ が導けました。標準偏差は $1/\lambda$ なので平均と一致し、変動係数が1であることも確認できます。
ここまでで指数分布の「位置」と「ばらつき」が分かりました。次は指数分布を指数分布たらしめている、最も特徴的な性質「無記憶性」を見ていきます。
無記憶性(memoryless property)
指数分布の最大の特徴は無記憶性(memoryless property)です。これは「すでにどれだけ待ったかに関係なく、これから先の待ち時間の分布が変わらない」という性質です。
数式で書くと、任意の $s, t \ge 0$ に対して次が成り立ちます。
$$ \begin{equation} P(X > s + t \mid X > s) = P(X > t) \end{equation} $$
左辺は「すでに時刻 $s$ までイベントが起きていない、という条件のもとで、さらに $t$ だけ起きずに待つ確率」です。これが、最初から $t$ だけ待つ確率 $P(X > t)$ に等しい、と主張しています。「過去を忘れて、いつもまっさらな状態からやり直す」というイメージです。
証明してみましょう。条件付き確率の定義から、
$$ \begin{equation} P(X > s + t \mid X > s) = \frac{P(X > s + t,\ X > s)}{P(X > s)} \end{equation} $$
です。ここで $s+t > s$ なので、事象 $\{X > s+t\}$ は $\{X > s\}$ に含まれます。よって同時確率は $P(X > s+t)$ に一致します。
$$ \begin{equation} P(X > s + t \mid X > s) = \frac{P(X > s + t)}{P(X > s)} \end{equation} $$
指数分布では $P(X > a) = \int_a^{\infty} \lambda e^{-\lambda x}\, dx = e^{-\lambda a}$ です（これが生存関数 $S(a) = e^{-\lambda a}$）。これを代入すると、
$$ \begin{equation} \frac{P(X > s + t)}{P(X > s)} = \frac{e^{-\lambda (s+t)}}{e^{-\lambda s}} = e^{-\lambda t} = P(X > t) \end{equation} $$
となり、無記憶性が示されました。指数の割り算で $s$ がきれいに消えるのがポイントです。連続分布で無記憶性をもつのは指数分布だけ（離散では幾何分布だけ）であることも知られています。
上の図は無記憶性を視覚化したものです。灰色が元の密度、赤色は「すでに $s=2$ まで待った」という条件のもとでの残り待ち時間の密度です。赤い曲線は元の密度を $s$ だけ右に平行移動して引き伸ばしただけで、形そのものは最初の指数分布とまったく同じです。10分待ったバス停でも、待ち始めたばかりのバス停でも、「あと何分待つか」の分布は変わらない — これが無記憶性の意味です。
無記憶性は便利な反面、「待てば待つほど来やすくなる（あるいは来にくくなる）」という現実の感覚とは合わないこともあります。なぜ指数分布だけがこの性質をもつのか。その答えは、指数分布がポアソン過程と深く結びついていることにあります。次の節で見ていきましょう。
ポアソン過程との関係 — 到着間隔としての指数分布
指数分布が自然界・社会のあちこちに現れる理由は、それがポアソン過程の到着間隔だからです。
ポアソン過程とは、「イベントが互いに独立に、単位時間あたり平均 $\lambda$ 回のレートでランダムに起こる」確率過程です。電話の着信、放射性崩壊、Webサーバへのアクセスなどがその典型例です。
このとき、時間 $[0, t]$ の間に起こるイベント数 $N(t)$ はポアソン分布に従います。
$$ \begin{equation} P(N(t) = k) = \frac{(\lambda t)^k e^{-\lambda t}}{k!} \end{equation} $$
ここで「最初のイベントが起こるまでの待ち時間 $X$」を考えます。$X > t$ であることは、「時間 $[0, t]$ の間に1回もイベントが起きていない」ことと同じです。すなわち $N(t) = 0$ です。上のポアソン分布で $k = 0$ とすると、
$$ \begin{equation} P(X > t) = P(N(t) = 0) = \frac{(\lambda t)^0 e^{-\lambda t}}{0!} = e^{-\lambda t} \end{equation} $$
が得られます。これはまさに指数分布の生存関数です。両辺を1から引けば累積分布関数 $F(t) = 1 – e^{-\lambda t}$、さらに微分すれば密度 $f(t) = \lambda e^{-\lambda t}$ が出てきます。つまりポアソン過程の到着間隔は指数分布に従うのです。
しかもポアソン過程の独立増分性から、最初の間隔だけでなく、2番目の間隔、3番目の間隔……どれも独立に同じ指数分布に従います。
import numpy as np
from scipy.stats import expon
import matplotlib.pyplot as plt
rng = np.random.default_rng(3)
lam = 1.5
T = 30.0
# 到着間隔を指数分布から生成し、累積して到着時刻を作る
gaps, t = [], 0.0
while True:
g = rng.exponential(1 / lam)
t += g
if t > T:
break
gaps.append(g)
gaps = np.array(gaps)
# 到着間隔のヒストグラムと理論密度
xx = np.linspace(0, gaps.max(), 200)
plt.figure(figsize=(9, 4))
plt.hist(gaps, bins=15, density=True, alpha=0.5, label="到着間隔のヒストグラム")
plt.plot(xx, expon.pdf(xx, scale=1 / lam), lw=2, label=f"指数分布 λ={lam}")
plt.xlabel("到着間隔")
plt.ylabel("密度")
plt.legend()
plt.show()
上の図は、レート $\lambda=1.5$ のポアソン過程を1本シミュレーションした結果です。上段はタイムライン上のイベント発生（赤）と到着間隔（青矢印）を示し、間隔がてんでバラバラであることが分かります。下段では、それら到着間隔のヒストグラムが理論曲線（赤）の指数分布によく一致しています。「回数を数えるとポアソン、間隔を測ると指数」という双対関係が、シミュレーションでも確認できました。
では、1回目ではなく $k$ 回目のイベントまでの時間はどんな分布になるのでしょうか。これは指数分布の「和」を考えることになり、ガンマ分布へとつながります。
ガンマ分布との関係 — 指数分布の和
ポアソン過程で「$k$ 回目のイベントが起こるまでの時間」を $T_k$ とすると、これは独立な $k$ 個の指数分布の間隔を足し合わせたものです。
$$ \begin{equation} T_k = X_1 + X_2 + \cdots + X_k, \qquad X_i \overset{\text{iid}}{\sim} \mathrm{Exp}(\lambda) \end{equation} $$
独立な指数分布の和はガンマ分布に従うことが知られています。
$$ \begin{equation} T_k \sim \mathrm{Gamma}(k,\ \theta = 1/\lambda) \end{equation} $$
ここで第1引数 $k$ は形状(shape)パラメータ、$\theta = 1/\lambda$ はスケールパラメータです。形状パラメータが $k=1$ のとき、ガンマ分布はちょうど指数分布に一致します。つまり指数分布はガンマ分布の特別な場合($k=1$)なのです。
ガンマ分布の密度は次のとおりです。
$$ \begin{equation} f_{T_k}(x) = \frac{\lambda^k}{(k-1)!}\, x^{k-1} e^{-\lambda x} \qquad (x \ge 0) \end{equation} $$
$k=1$ を代入すると $(k-1)! = 0! = 1$, $x^{k-1} = x^0 = 1$ なので、$f(x) = \lambda e^{-\lambda x}$ となり、確かに指数分布に戻ります。
import numpy as np
from scipy.stats import gamma
import matplotlib.pyplot as plt
x = np.linspace(0, 12, 400)
lam = 1.0
plt.figure(figsize=(9, 4))
for k in [1, 2, 3, 5]:
# k 個の指数(scale=1/lam)の和 = Gamma(a=k, scale=1/lam)
y = gamma.pdf(x, a=k, scale=1 / lam)
label = "k=1 : 指数分布" if k == 1 else f"k={k} : Gamma"
plt.plot(x, y, lw=2, label=label)
plt.xlabel("x （k回目のイベントまでの時間）")
plt.ylabel("密度")
plt.legend()
plt.show()
上の図から、$k$ が増えるにつれて分布の山が右へ移動し、ピークがはっきりした釣鐘型に近づくことが分かります。$k=1$ では原点で最大の指数分布、$k$ が大きくなると複数の待ち時間が積み重なるため、平均 $k/\lambda$ のあたりに山ができます。中心極限定理を予感させる挙動です。指数分布が「1つの待ち時間」、ガンマ分布が「複数の待ち時間の合計」を表す、という関係が視覚的にも納得できます。
ここで「故障率」という観点から指数分布を見ると、もう一つ重要な性質が見えてきます。次はハザード関数を扱います。
ハザード関数 — 故障率が一定であること
信頼性工学や生存解析ではハザード関数(hazard function) $h(t)$ が重要です。これは「時刻 $t$ まで生き残った（壊れていない）という条件のもとで、ちょうど今この瞬間に壊れる率」を表します。直感的には「瞬間故障率」です。定義は次のとおりです。
$$ \begin{equation} h(t) = \frac{f(t)}{S(t)} \end{equation} $$
ここで $f(t)$ は密度、$S(t) = P(X > t)$ は生存関数です。指数分布では $f(t) = \lambda e^{-\lambda t}$, $S(t) = e^{-\lambda t}$ なので、これを代入すると、
$$ \begin{equation} h(t) = \frac{\lambda e^{-\lambda t}}{e^{-\lambda t}} = \lambda \end{equation} $$
となります。指数の項が分子・分母で打ち消し合い、ハザードが時刻によらない定数 $\lambda$ になるのです。これは無記憶性のもう一つの言い換えにほかなりません。「使い古した部品も新品も、壊れる率は同じ」 — これが指数分布が表す世界です。
import numpy as np
from scipy.stats import expon, weibull_min
import matplotlib.pyplot as plt
t = np.linspace(0.01, 3, 300)
plt.figure(figsize=(9, 4))
# 指数分布(ワイブルの形状 k=1)・ワイブル k=0.6, 2 のハザードを比較
for k, name in [(0.6, "k=0.6 故障率が減少"),
(1.0, "k=1 一定（指数分布）"),
(2.0, "k=2 故障率が増加")]:
d = weibull_min(c=k, scale=1.0)
h = d.pdf(t) / d.sf(t)
plt.plot(t, h, lw=2, label=name)
plt.xlabel("時刻 t")
plt.ylabel("ハザード h(t)")
plt.ylim(0, 4)
plt.legend()
plt.show()
左側の図は、指数分布のハザードが時刻によらず水平な直線（一定値 $\lambda$）になることを示しています。右側はワイブル分布との対比です。ワイブルは形状パラメータ $k$ によってハザードが増加（摩耗故障)・一定・減少（初期故障)のいずれにもなり、$k=1$ のときちょうど指数分布に一致します。つまり指数分布は「故障率が一定」という理想化された、最も単純な寿命モデルだと位置づけられます。現実の部品は初期故障や摩耗で故障率が変動するため、より一般的なワイブルが使われることも多いですが、指数分布はその出発点です。
理論的な性質がそろったので、次はデータから $\lambda$ を推定する方法、すなわち最尤推定を扱います。
最尤推定 — データから λ を推定する
観測データ $x_1, x_2, \dots, x_n$ が指数分布に従うとき、レート $\lambda$ をデータから推定したいとします。最尤推定(maximum likelihood estimation)を使いましょう。
まず尤度（データが得られる確率密度の積）は、各データが独立なので、
$$ \begin{equation} L(\lambda) = \prod_{i=1}^{n} \lambda e^{-\lambda x_i} = \lambda^n \exp\!\left( -\lambda \sum_{i=1}^{n} x_i \right) \end{equation} $$
です。積のままでは扱いにくいので、対数をとって対数尤度にします。対数は単調増加なので、最大化問題の答えは変わりません。
$$ \begin{equation} \ell(\lambda) = \log L(\lambda) = n \log \lambda – \lambda \sum_{i=1}^{n} x_i \end{equation} $$
これを $\lambda$ で微分して0と置きます。最大値を探すためです。
$$ \begin{equation} \frac{d\ell}{d\lambda} = \frac{n}{\lambda} – \sum_{i=1}^{n} x_i = 0 \end{equation} $$
この式を $\lambda$ について解きます。$\frac{n}{\lambda} = \sum x_i$ より、
$$ \begin{equation} \hat{\lambda} = \frac{n}{\sum_{i=1}^{n} x_i} = \frac{1}{\bar{x}} \end{equation} $$
が得られます。ここで $\bar{x}$ は標本平均です。$\ell$ の2階微分は $-n/\lambda^2 < 0$ なので、この点は確かに最大値です。結果は非常にきれいで、「レートの最尤推定値は標本平均の逆数」になります。期待値が $1/\lambda$ だったことを思い出せば、これも自然な結果です（標本平均で $1/\lambda$ を推定し、その逆数をとる)。
import numpy as np
from scipy.stats import expon
rng = np.random.default_rng(7)
lam_true = 1.5
data = rng.exponential(1 / lam_true, size=60)
# 最尤推定: lambda_hat = 1 / 標本平均
lam_hat = 1 / data.mean()
print(f"真値 lambda = {lam_true}")
print(f"最尤推定 lambda_hat = {lam_hat:.3f}")
# SciPy の fit でも確認（loc=0 に固定）
loc, scale = expon.fit(data, floc=0)
print(f"SciPy fit scale = {scale:.3f} -> lambda = {1/scale:.3f}")
このコードを実行すると、lambda_hat がおよそ 1.4〜1.6 程度になり、真値 1.5 に近い値が得られます。SciPyの expon.fit で求めた scale の逆数も同じ値になり、手計算の公式 $\hat{\lambda} = 1/\bar{x}$ と一致することが確認できます。サンプル数が少ないと推定値は真値からずれますが、$n$ を増やせば収束していきます。
上の図は対数尤度 $\ell(\lambda)$ を $\lambda$ の関数としてプロットしたものです。曲線が最大になる点（赤破線）がちょうど $\hat{\lambda} = 1/\bar{x}$ であり、真値（緑点線）の近くに位置しています。最尤推定が「データを最もよく説明する $\lambda$ を選ぶ」操作であることが、山のピークとして直感的に見てとれます。
推定の道具がそろったので、ここから先は指数分布が実際にどう使われるのかを見ていきます。まず累積分布関数と生存関数を整理しておきましょう。
累積分布関数・生存関数・分位点
これまで何度か登場した累積分布関数(CDF)と生存関数を、ここで一度まとめておきます。CDFは「時刻 $x$ までにイベントが起こる確率」です。
$$ \begin{equation} F(x) = P(X \le x) = 1 – e^{-\lambda x} \end{equation} $$
生存関数(SF)はその補集合で、「時刻 $x$ まで起きない確率」です。
$$ \begin{equation} S(x) = 1 – F(x) = e^{-\lambda x} \end{equation} $$
中央値(50%点)は $F(x) = 0.5$ を解いて求めます。$1 – e^{-\lambda x} = 0.5$ より $e^{-\lambda x} = 0.5$、両辺の対数をとって、
$$ \begin{equation} x_{\text{med}} = \frac{\ln 2}{\lambda} \approx \frac{0.693}{\lambda} \end{equation} $$
となります。中央値 $0.693/\lambda$ は平均 $1/\lambda$ より小さいことに注目してください。これは分布が右に長い裾を引いているためで、「半数のイベントは平均より早く起こる」ことを意味します。
上の図はCDF $F(x) = 1 – e^{-\lambda x}$ を描き、中央値（赤破線、$F=0.5$ の交点)と平均（緑破線)の位置を示したものです。中央値が平均より明確に左にあることが見てとれ、指数分布の右に偏った（右裾の長い)性質を反映しています。CDFは $x=0$ で0から立ち上がり、$x \to \infty$ で1に漸近します。
生存関数は信頼性・生存解析の主役なので、$\lambda$ を変えて改めて描いてみましょう。
import numpy as np
from scipy.stats import expon
import matplotlib.pyplot as plt
x = np.linspace(0, 6, 400)
plt.figure(figsize=(9, 4))
for lam in [0.5, 1.0, 2.0]:
S = expon.sf(x, scale=1 / lam) # 生存関数 = 1 - CDF
plt.plot(x, S, lw=2, label=f"λ={lam} （平均寿命 {1/lam:.2f}）")
plt.axhline(0.5, color="gray", ls=":")
plt.xlabel("時刻 t")
plt.ylabel("生存確率 S(t)")
plt.legend()
plt.show()
上の図は生存曲線 $S(t) = e^{-\lambda t}$ を3つの $\lambda$ で比較したものです。$\lambda$ が小さい（平均寿命が長い)ほど曲線はゆっくり減衰し、長く生き残ることが分かります。どの曲線も同じ「指数的減衰」の形をしていますが、減り方の速さだけが $\lambda$ で決まっています。生存曲線が0.5を横切る時刻が、先ほど計算した中央値 $\ln 2/\lambda$ です。
最後に、これらの道具を使った代表的な応用と、実データへのフィットを見ていきます。
応用 — 信頼性・待ち行列・生存解析
指数分布は理論的な美しさだけでなく、実務で広く使われています。代表的な3つの応用を整理します。
信頼性工学では、故障率が一定とみなせる電子部品の寿命モデルに指数分布を使います。平均故障時間(MTBF, mean time between failures)はちょうど $1/\lambda$ です。たとえばMTBFが10000時間の部品なら $\lambda = 1/10000$ で、5000時間後に生き残っている確率は $S(5000) = e^{-5000/10000} = e^{-0.5} \approx 0.607$ と即座に計算できます。
待ち行列理論では、客の到着間隔とサービス時間をともに指数分布でモデル化したものがM/M/1待ち行列です。到着レート $\lambda$、サービスレート $\mu$ のとき、システムが安定する条件 $\lambda < \mu$ のもとで、平均待ち時間や平均系内人数を閉じた式で求められます。指数分布の無記憶性が、こうした解析を劇的に簡単にしています。
生存時間解析では、患者の生存時間や機械の稼働時間を扱います。指数分布は最も単純な生存モデルで、ハザードが一定という仮定が現実に合うかをまず確認し、合わなければワイブル分布やCox比例ハザードモデルなどに進みます。
実データが指数分布にどれだけ従うかを確認するには、ヒストグラムへのフィットとQ-Qプロットが定番です。
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
rng = np.random.default_rng(11)
lam_true = 0.8
data = rng.exponential(1 / lam_true, size=400)
lam_hat = 1 / data.mean()
fig, axes = plt.subplots(1, 2, figsize=(12, 4))
# 左: ヒストグラムと最尤フィット
xx = np.linspace(0, data.max(), 300)
axes[0].hist(data, bins=30, density=True, alpha=0.5, label="観測データ")
axes[0].plot(xx, stats.expon.pdf(xx, scale=1 / lam_hat), lw=2,
label=f"フィット λ={lam_hat:.3f}")
axes[0].set_xlabel("到着間隔"); axes[0].set_ylabel("密度"); axes[0].legend()
# 右: Q-Qプロット（点が直線に乗れば適合良好）
stats.probplot(data, dist=stats.expon(scale=1 / lam_hat),
fit=False, plot=axes[1])
axes[1].set_title("Q-Qプロット")
plt.show()
左の図では、観測データのヒストグラムに最尤フィットした指数分布（赤線)がよく重なっています。右のQ-Qプロットでは、標本分位点と理論分位点を対応させた点が直線 $y=x$ にほぼ乗っており、データが指数分布によく従っていることが確認できます。実務では、まずこのようにフィットの良さを目で確認してから、平均故障時間や待ち時間といった量を推定します。点が直線から大きく外れる場合は、指数分布の仮定（一定ハザード）が成り立っていないサインで、別の分布を検討すべきです。
まとめ
本記事では、指数分布について直感から応用まで解説しました。
- 指数分布は「ランダムに一定ペースで起こる事象の待ち時間」を表し、密度は $f(x) = \lambda e^{-\lambda x}$
- 期待値は $1/\lambda$、分散は $1/\lambda^2$ で、いずれも部分積分で導出できる。標準偏差が平均と等しく変動係数は常に1
- 無記憶性 $P(X>s+t \mid X>s) = P(X>t)$ をもち、連続分布ではこの性質をもつのは指数分布だけ
- ポアソン過程の到着間隔であり、その和はガンマ分布（$k=1$ が指数分布)になる
- ハザードが一定 $h(t)=\lambda$ で、ワイブル分布の $k=1$ の特別な場合
- 最尤推定は $\hat{\lambda} = 1/\bar{x}$（標本平均の逆数）と非常に単純
- 信頼性・待ち行列・生存解析で出発点となる基礎的な分布
指数分布は単独でも重要ですが、ポアソン分布・ポアソン過程・ガンマ分布・ワイブル分布をつなぐハブのような存在です。ここを足場にすると、確率過程や信頼性工学の世界が一気に見通しよくなります。
次のステップとして、以下の記事も参考にしてください。