「3つの薬の効果に差があるか」「5つの工場の品質に違いがあるか」— このような「複数グループの平均の差の検定」は、t検定を繰り返すのではなく分散分析(ANOVA)で行います。ANOVAの検定統計量が従う分布がF分布(F-distribution)です。
F分布はロナルド・フィッシャー(Ronald Fisher)にちなんで名付けられ、統計的推測の中核をなす分布です。
- 分散分析(ANOVA): 複数グループの平均の比較
- 等分散検定: 2つの母集団の分散の比の検定
- 回帰分析: モデル全体の有意性のF検定
- 実験計画法: 要因効果の検定
本記事の内容
- F分布の定義(カイ二乗分布の比)
- 密度関数と期待値・分散
- t分布との関係
- 分散分析(一元配置ANOVA)への応用
- Pythonでの実装と可視化
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
F分布の定義
カイ二乗分布の比としての導出
$U \sim \chi^2(d_1)$ と $V \sim \chi^2(d_2)$ が独立のとき
$$ \begin{equation} F = \frac{U/d_1}{V/d_2} \sim F(d_1, d_2) \end{equation} $$
$d_1$ は分子の自由度、$d_2$ は分母の自由度です。
直感的には、F分布は「2つの分散の比」の分布です。分子も分母もカイ二乗分布(正規分布の二乗和)をその自由度で割ったものであり、もし2つの母分散が等しければ $F \approx 1$ になるはずです。
密度関数
$$ f(x) = \frac{1}{B(d_1/2, d_2/2)} \left(\frac{d_1}{d_2}\right)^{d_1/2} \frac{x^{d_1/2 – 1}}{(1 + d_1 x/d_2)^{(d_1+d_2)/2}}, \quad x > 0 $$
期待値と分散
$$ E[F] = \frac{d_2}{d_2 – 2} \quad (d_2 > 2) $$
$$ \text{Var}(F) = \frac{2d_2^2(d_1 + d_2 – 2)}{d_1(d_2-2)^2(d_2-4)} \quad (d_2 > 4) $$
$d_2 > 2$ のとき $E[F] > 1$ であり、$d_2 \to \infty$ で $E[F] \to 1$ です。
t分布との関係
$T \sim t(\nu)$ のとき $T^2 \sim F(1, \nu)$ です。つまり、2つのグループの平均の比較(t検定)はF検定($d_1 = 1$ のANOVA)と等価です。
分散分析(一元配置ANOVA)
問題設定
$k$ 個のグループがあり、グループ $i$ から $n_i$ 個のデータ $X_{i1}, \ldots, X_{in_i}$ が得られています。
$$ X_{ij} = \mu_i + \varepsilon_{ij}, \quad \varepsilon_{ij} \sim N(0, \sigma^2) $$
帰無仮説 $H_0: \mu_1 = \mu_2 = \cdots = \mu_k$ を検定します。
F統計量
グループ間変動(SSB)とグループ内変動(SSW)を計算します。
$$ \text{SSB} = \sum_{i=1}^{k} n_i (\bar{X}_{i\cdot} – \bar{X}_{\cdot\cdot})^2, \quad \text{SSW} = \sum_{i=1}^{k} \sum_{j=1}^{n_i} (X_{ij} – \bar{X}_{i\cdot})^2 $$
F統計量は
$$ F = \frac{\text{SSB}/(k-1)}{\text{SSW}/(N-k)} = \frac{\text{MSB}}{\text{MSW}} $$
$H_0$ のもとで $F \sim F(k-1, N-k)$ に従います。
Pythonでの実装と可視化
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# (a) 異なる自由度でのPDF
ax = axes[0, 0]
x = np.linspace(0.01, 5, 500)
params_f = [(1, 1), (2, 5), (5, 5), (10, 10), (10, 50), (50, 50)]
for d1, d2 in params_f:
pdf = stats.f.pdf(x, d1, d2)
ax.plot(x, pdf, linewidth=2, label=f'F({d1},{d2})')
ax.set_xlabel('x', fontsize=12)
ax.set_ylabel('f(x)', fontsize=12)
ax.set_title('F-Distribution PDF', fontsize=13)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
ax.set_ylim(0, 1.2)
# (b) T^2 = F(1,ν) の確認
ax = axes[0, 1]
np.random.seed(42)
nu = 10
n_sim = 100000
t_samples = np.random.standard_t(nu, n_sim)
f_from_t = t_samples**2
x_range = np.linspace(0, 8, 300)
ax.hist(f_from_t, bins=100, density=True, range=(0, 8),
color='steelblue', alpha=0.6, label='$T^2$ (simulation)')
ax.plot(x_range, stats.f.pdf(x_range, 1, nu), 'r-', linewidth=2.5,
label=f'F(1, {nu}) theory')
ax.set_xlabel('x', fontsize=12)
ax.set_ylabel('Density', fontsize=12)
ax.set_title(f'$T^2 \\sim F(1, \\nu)$ verification ($\\nu$={nu})', fontsize=12)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
# (c) ANOVAのシミュレーション
ax = axes[1, 0]
np.random.seed(42)
n_anova_sim = 10000
# H0のもとでのF統計量
f_null = []
k = 3
n_per_group = 20
for _ in range(n_anova_sim):
groups = [np.random.normal(0, 1, n_per_group) for _ in range(k)]
all_data = np.concatenate(groups)
grand_mean = np.mean(all_data)
ssb = sum(n_per_group * (np.mean(g) - grand_mean)**2 for g in groups)
ssw = sum(np.sum((g - np.mean(g))**2) for g in groups)
f_stat = (ssb/(k-1)) / (ssw/(k*n_per_group-k))
f_null.append(f_stat)
x_range = np.linspace(0, 6, 300)
ax.hist(f_null, bins=80, density=True, color='steelblue', alpha=0.6,
label='Simulation (H0 true)')
ax.plot(x_range, stats.f.pdf(x_range, k-1, k*n_per_group-k), 'r-',
linewidth=2.5, label=f'F({k-1}, {k*n_per_group-k}) theory')
# 棄却域
f_crit = stats.f.ppf(0.95, k-1, k*n_per_group-k)
ax.axvline(f_crit, color='green', linestyle='--', linewidth=2,
label=f'Critical value (α=0.05): {f_crit:.2f}')
ax.fill_between(x_range[x_range >= f_crit],
stats.f.pdf(x_range[x_range >= f_crit], k-1, k*n_per_group-k),
alpha=0.3, color='red', label='Rejection region')
ax.set_xlabel('F statistic', fontsize=12)
ax.set_ylabel('Density', fontsize=12)
ax.set_title(f'One-way ANOVA (k={k}, n={n_per_group})', fontsize=12)
ax.legend(fontsize=8)
ax.grid(True, alpha=0.3)
# (d) 検出力の比較
ax = axes[1, 1]
effect_sizes = np.linspace(0, 2, 50)
powers = []
for delta in effect_sizes:
rejections = 0
for _ in range(2000):
groups = [np.random.normal(i*delta, 1, n_per_group) for i in range(k)]
all_data = np.concatenate(groups)
grand_mean = np.mean(all_data)
ssb = sum(n_per_group * (np.mean(g) - grand_mean)**2 for g in groups)
ssw = sum(np.sum((g - np.mean(g))**2) for g in groups)
f_stat = (ssb/(k-1)) / (ssw/(k*n_per_group-k))
if f_stat > f_crit:
rejections += 1
powers.append(rejections / 2000)
ax.plot(effect_sizes, powers, 'b-', linewidth=2.5)
ax.axhline(0.05, color='red', linestyle='--', linewidth=1.5,
label='$\\alpha$ = 0.05')
ax.axhline(0.80, color='green', linestyle=':', linewidth=1.5,
label='Power = 0.80')
ax.set_xlabel('Effect size ($\\delta$)', fontsize=12)
ax.set_ylabel('Power', fontsize=12)
ax.set_title('ANOVA Power Curve', fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('f_distribution.png', dpi=150, bbox_inches='tight')
plt.show()
この可視化から、F分布とANOVAの関係が確認できます。
-
左上(PDF): 自由度が大きくなるにつれてピークが1に近づき、分布が集中します。$F(1,1)$ は右に長い裾を持ちます。
-
右上($T^2 = F(1,\nu)$): t分布の二乗がF分布に一致することが数値的に確認されています。
-
左下(ANOVA): $H_0$ のもとでのF統計量のヒストグラムが理論的なF分布(赤い実線)と一致しています。臨界値(緑の破線)を超えた領域が棄却域です。
-
右下(検出力): 効果量が大きくなるにつれて検出力が上昇します。一般的な目安である検出力80%を達成するのに必要な効果量が読み取れます。
まとめ
本記事では、F分布の導出から分散分析への応用まで解説しました。
- F分布は2つの独立なカイ二乗分布の比(各自由度で割ったもの)として定義される
- t分布との関係: $T^2 \sim F(1, \nu)$
- 分散分析(ANOVA)のF統計量は「グループ間変動/グループ内変動」の比であり、帰無仮説のもとでF分布に従う
- 自由度が大きくなると分布は1の周りに集中する
次のステップとして、以下の記事も参考にしてください。