新薬の臨床試験で「治療群と対照群の平均値に差があるか」を検定したいとき、私たちは通常 $t$ 検定を使います。しかし、$t$ 検定には「データが正規分布に従う」という仮定が必要です。もしデータの分布が正規分布から大きくずれていたら、$t$ 検定の結果はどこまで信頼できるのでしょうか。
このような「分布の仮定に頼りたくない」場面で力を発揮するのが並べ替え検定(permutation test)です。置換検定(permutation test)やランダム化検定(randomization test)とも呼ばれます。この検定の核心的なアイデアは驚くほどシンプルです。「もし帰無仮説が正しければ、グループのラベル(治療群/対照群)を入れ替えても統計量の値は大きく変わらないはずだ」という発想です。
並べ替え検定を理解すると、以下のような幅広い応用が開けます。
- 臨床試験の解析: 正規性の仮定が疑わしい小標本での治療効果の検定に使われます
- ゲノミクス・バイオインフォマティクス: 遺伝子発現データなど、分布の仮定が困難な高次元データの解析で標準的な手法です
- 機械学習のモデル評価: 予測精度の差が統計的に有意かどうかを、分布の仮定なしに検定できます
- 生態学・環境科学: 種の分布パターンや生態系の多様性指標の比較など、正規性を仮定しにくい場面で広く使われています
本記事の内容
- 並べ替え検定の基本アイデアと直感的な理解
- 帰無分布の正確な構成と漸近的性質
- 正確検定と近似検定(モンテカルロ並べ替え検定)
- 条件付き検定としての理論的位置づけ
- Pythonでの実装と$t$検定との比較シミュレーション
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
- 仮説検定の基礎 — p値、有意水準、第1種・第2種の過誤の定義
- t検定の理論 — 二標本検定の基本的な枠組み
- 大数の法則と中心極限定理 — 漸近理論の基礎
並べ替え検定の基本アイデア
カードゲームで理解する並べ替え検定
並べ替え検定のアイデアを、身近な例で理解しましょう。
10人の患者がいて、5人が新薬を、5人がプラセボを投与されたとします。10人の血圧の変化量が次のように記録されました。
| 新薬群 | -8, -5, -12, -3, -7 |
|---|---|
| プラセボ群 | -2, +1, -4, +2, -1 |
新薬群の平均変化量は $-7.0$、プラセボ群の平均変化量は $-0.8$ で、差は $-6.2$ です。この差は「本当に新薬の効果」なのでしょうか、それとも「偶然の割り当てで生じた差」にすぎないのでしょうか。
ここで、こう考えてみます。もし新薬にまったく効果がないなら、10人の血圧変化量は「誰が新薬を飲んだか」に関係なく同じはずです。つまり、10枚の数値カード(-8, -5, -12, -3, -7, -2, +1, -4, +2, -1)をランダムに5枚ずつ2つの山に分けても、平均の差は今回観測された $-6.2$ と同程度になりうるはずです。
そこで、10枚のカードのすべての分け方($\binom{10}{5} = 252$ 通り)について平均の差を計算し、観測された差 $-6.2$ 以上に極端な差がどれだけの割合で現れるかを調べます。この割合がp値です。
もしこのp値が非常に小さければ(例えば $5\%$ 以下)、「帰無仮説(薬の効果なし)のもとでは今回のような大きな差はめったに起きない」と結論し、帰無仮説を棄却します。
正式な定式化
上の直感を数学的に定式化しましょう。
データの構造: $n$ 個の観測値 $z_1, z_2, \dots, z_n$ が2つの群に分かれているとします。群1に $n_1$ 個、群2に $n_2$ 個($n_1 + n_2 = n$)の観測値があり、群の割り当ては $\bm{g} = (g_1, g_2, \dots, g_n)$ で表されます($g_i = 1$ または $2$)。
帰無仮説: $H_0$: 群の割り当ては観測値に影響しない。つまり、$z_i$ の同時分布は $g_i$ に依存しない。
検定統計量: 群間の差を測る統計量 $T = T(z_1, \dots, z_n, g_1, \dots, g_n)$。例えば、群平均の差 $T = \bar{z}_1 – \bar{z}_2$ です。
帰無分布の構成: 帰無仮説のもとでは、群の割り当て $\bm{g}$ はデータの値に影響しないため、すべての可能な割り当て(置換)が等確率で生じます。$n$ 個の観測値を $n_1$ 個と $n_2$ 個に分けるすべての方法 $\binom{n}{n_1}$ 通りについて検定統計量 $T$ を計算し、その分布(並べ替え分布、permutation distribution)を帰無分布として使います。
p値の計算:
$$ \begin{equation} p = \frac{\#\{T_\pi : |T_\pi| \geq |T_{\text{obs}}|\}}{\binom{n}{n_1}} \end{equation} $$
ここで $T_{\text{obs}}$ は観測された検定統計量の値、$T_\pi$ は各置換 $\pi$ に対する検定統計量の値です。
この方法の美しさは、確率分布に関する一切の仮定を必要としないことです。正規分布である必要も、等分散である必要もありません。帰無仮説のもとでデータの値が群のラベルに依存しないという交換可能性(exchangeability)だけが必要です。
しかし、すべての置換を列挙する正確検定は、$n$ が大きいと計算量が爆発します。この問題にどう対処するかを次のセクションで見ていきましょう。
正確検定と近似検定
正確並べ替え検定
すべての可能な置換を列挙して帰無分布を構成する方法を正確並べ替え検定(exact permutation test)と呼びます。
$\binom{n}{n_1}$ 通りのすべての置換について検定統計量を計算するため、p値は離散的な値のみを取ります。例えば、$n = 10$, $n_1 = 5$ の場合、可能な値は $\frac{k}{252}$($k = 0, 1, \dots, 252$)に限られます。
正確検定の利点は次の通りです。
- 有限標本で正確なサイズを持つ: 帰無仮説のもとで、棄却確率が名目水準以下であることが保証されます
- 分布の仮定が不要: パラメトリックな仮定を一切置きません
- 小標本でも適用可能: 漸近近似に頼らないため、$n$ が小さくても使えます
一方、計算量の問題があります。$n = 20$, $n_1 = 10$ でも $\binom{20}{10} = 184{,}756$ 通りの計算が必要で、$n = 30$, $n_1 = 15$ では $\binom{30}{15} \approx 1.55 \times 10^8$ 通りになります。$n$ が40を超えると、すべての置換を列挙するのは現実的ではありません。
モンテカルロ並べ替え検定
計算量の問題を解決するのがモンテカルロ並べ替え検定(Monte Carlo permutation test)です。すべての置換を列挙する代わりに、$B$ 回のランダムな置換をサンプリングして帰無分布を近似します。
具体的な手順は以下の通りです。
- 観測データから検定統計量 $T_{\text{obs}}$ を計算する
- $b = 1, 2, \dots, B$ について: – データの群ラベルをランダムにシャッフルする(ランダムな置換 $\pi_b$ を生成する) – シャッフルしたデータで検定統計量 $T_{\pi_b}$ を計算する
- p値を次のように計算する:
$$ \begin{equation} \hat{p} = \frac{\#\{b : |T_{\pi_b}| \geq |T_{\text{obs}}|\} + 1}{B + 1} \end{equation} $$
分子と分母に $+1$ を加えているのは、観測されたデータ自体も置換の一つとして数えるためです。この補正により、$\hat{p}$ の値が0にならないことが保証されます。
モンテカルロ近似の精度
モンテカルロ並べ替え検定のp値は、正確なp値の推定量です。$B$ 回のランダム置換に基づくp値の推定精度は、二項分布の性質から次のように評価できます。
真のp値を $p^*$ とすると、$\hat{p}$ の標準誤差は近似的に次のようになります。
$$ \text{SE}(\hat{p}) \approx \sqrt{\frac{p^*(1 – p^*)}{B}} $$
$p^* = 0.05$ のとき、$B = 9{,}999$ とすれば $\text{SE} \approx 0.0022$ であり、$p^*$ の $95\%$ 信頼区間の幅は約 $0.009$ です。つまり、$B = 10{,}000$ 程度のランダム置換で、$p = 0.05$ 付近の判定には十分な精度が得られます。
$B = 999$ でも多くの場合は実用上十分ですが、有意水準の境界付近での判定がクリティカルな場合は $B = 9{,}999$ や $B = 99{,}999$ に増やすことが推奨されます。
正確検定とモンテカルロ近似の関係を整理しました。次に、並べ替え検定の理論的な性質を詳しく見ていきましょう。
並べ替え検定の理論的性質
条件付き検定としての位置づけ
並べ替え検定は、データの順序統計量(観測値の値そのもの)を条件として固定した条件付き検定として理解できます。帰無仮説のもとでは、群のラベル $\bm{g}$ は観測値 $\bm{z} = (z_1, \dots, z_n)$ に対して等確率で割り当てられるため、条件付き分布は離散一様分布になります。
$$ P(\bm{g} = \bm{g}_\pi | \bm{z}) = \frac{1}{\binom{n}{n_1}} \quad \text{for all permutations } \pi $$
この条件付き分布を帰無分布として使うため、並べ替え検定のサイズ(第1種の過誤率)は条件付きで正確に制御されます。
$$ \begin{equation} P(T_\pi \geq c_\alpha | \bm{z}) \leq \alpha \quad \text{for all } \bm{z} \end{equation} $$
条件付きでサイズが制御されることは、無条件でもサイズが制御されることを意味します(全確率の法則により)。つまり、並べ替え検定は正確検定(exact test)です。
交換可能性の条件
並べ替え検定が正しく機能するための理論的な条件は、帰無仮説のもとでデータが交換可能(exchangeable)であることです。交換可能性とは、データの同時分布が添え字の入れ替えに対して不変であること、すなわち任意の置換 $\pi$ に対して次が成り立つことです。
$$ (Z_1, Z_2, \dots, Z_n) \overset{d}{=} (Z_{\pi(1)}, Z_{\pi(2)}, \dots, Z_{\pi(n)}) $$
独立同一分布(i.i.d.)のデータは交換可能ですが、交換可能性はi.i.d.より弱い条件です。例えば、同一の分布からの標本であれば、たとえ独立でなくても交換可能であり得ます。
しかし、交換可能性が成り立たない場合には注意が必要です。例えば、2つの群の分散が帰無仮説のもとでも異なる場合、群平均の差を検定統計量とする並べ替え検定は正しいサイズを持ちません。この問題については後で詳しく議論します。
検出力に関する性質
並べ替え検定は、適切に構成された場合、無条件の最も強力な検定に近い性能を持つことが知られています。特に以下の結果が知られています。
命題(Hoeffding, 1952): 帰無仮説のもとでデータがi.i.d.であり、検定統計量として十分統計量を用いた場合、並べ替え検定は漸近的に最も強力な非条件付き検定と同等の検出力を持ちます。
つまり、並べ替え検定は「分布の仮定を置かないから検出力が劣る」ということは必ずしもなく、むしろ分布の仮定が正しくない場合にはパラメトリック検定よりも高い検出力を示すことがあります。
漸近的な等価性
大標本のもとでは、並べ替え検定と対応するパラメトリック検定は漸近的に同等になることが多いです。例えば、二標本の群平均の差に基づく並べ替え検定は、$n_1, n_2 \to \infty$ のもとで二標本 $t$ 検定と漸近的に同じ棄却領域を持ちます。
これは、並べ替え分布が漸近的に正規分布に近づくためです。中心極限定理の類似により、置換の数が十分に大きければ、並べ替え分布はガウス分布で良く近似されます。
並べ替え検定の理論的基盤を理解したところで、次にさまざまな検定統計量の選択と応用例を見ていきましょう。
検定統計量の選択
二標本問題での選択肢
並べ替え検定の枠組みでは、検定統計量の選択は自由です。帰無仮説のもとでの分布が並べ替えによって正確に求まるため、どのような統計量を使っても正確なサイズが保証されます。しかし、検定統計量の選択は検出力に影響します。
よく使われる統計量には以下のものがあります。
平均の差:
$$ T_1 = \bar{X}_1 – \bar{X}_2 $$
最もシンプルで、位置の差(location shift)に対して自然な統計量です。
t統計量:
$$ T_2 = \frac{\bar{X}_1 – \bar{X}_2}{s_p\sqrt{1/n_1 + 1/n_2}} $$
ここで $s_p$ はプールされた標準偏差です。$t$ 統計量はスケールの違いを標準化するため、分散が群間で異なる場合にも頑健です。並べ替え検定では $t$ 統計量の分布に関する仮定($t$ 分布に従うこと)は不要で、帰無分布は並べ替えから直接構成されます。
ウィルコクソンの順位和統計量:
$$ T_3 = \sum_{i \in \text{群1}} R_i $$
ここで $R_i$ は $i$ 番目の観測値の全体における順位です。順位を使うことで外れ値の影響を軽減できます。実は、ウィルコクソンの順位和検定自体が並べ替え検定の特殊ケースです。
コルモゴロフ・スミルノフ統計量:
$$ T_4 = \sup_x |F_{1,n_1}(x) – F_{2,n_2}(x)| $$
2つの群の経験分布関数の最大の乖離を測ります。分布の形状の違い(位置だけでなくスケールや形状も含む)を検出できますが、位置の差のみに焦点を当てる場合は検出力が低下します。
検定統計量の選択指針
検定統計量の選択は、対立仮説で検出したい差の種類に依存します。
| 検出したい差 | 推奨する統計量 |
|---|---|
| 平均の差(正規に近い分布) | 平均の差、$t$ 統計量 |
| 位置の差(歪んだ分布) | ウィルコクソン順位和 |
| 分布全体の差 | コルモゴロフ・スミルノフ |
| 分散の差 | ルビーン検定統計量 |
一般的には、$t$ 統計量が最も汎用的な選択肢です。並べ替え検定では分布の仮定が不要なので、$t$ 統計量を使うことの不都合はありません。むしろ、$t$ 統計量は分散を考慮するため、純粋な平均の差よりも多くの場合で高い検出力を示します。
具体的な実装に移る前に、並べ替え検定を一般的な対応のある(ペアの)データにどう適用するかを見ましょう。
対応のあるデータへの拡張
ペアの並べ替え検定
対応のある(paired)データでは、各ペアの差を計算し、その差の符号をランダムに入れ替えることで帰無分布を構成します。
$n$ 組の対応するデータ $(X_{i1}, X_{i2})$ から差 $D_i = X_{i1} – X_{i2}$ を計算します。帰無仮説 $H_0: E[D_i] = 0$ のもとでは、$D_i$ と $-D_i$ は同じ確率で生じるはずです。
したがって、各 $D_i$ の符号をランダムに $+1$ または $-1$ にする $2^n$ 通りの符号反転について検定統計量を計算し、帰無分布を構成します。
$$ T_\pi = \frac{1}{n}\sum_{i=1}^{n} s_i D_i, \quad s_i \in \{-1, +1\} $$
これは符号並べ替え検定(sign permutation test)と呼ばれ、$2^n$ 通りの帰無分布を持ちます。$n = 20$ では $2^{20} \approx 10^6$ 通りとなり、正確計算が可能な境界にあります。
多群への拡張
$k$ 群の比較(一元配置の設定)への拡張も自然です。$n$ 個のデータを $k$ つの群に割り当てるすべての方法 $\frac{n!}{n_1! n_2! \cdots n_k!}$ 通りについて検定統計量を計算します。
検定統計量としては、通常の $F$ 統計量(分散分析の統計量)がよく使われます。
$$ T = \frac{\text{群間変動} / (k-1)}{\text{群内変動} / (n-k)} $$
この検定は並べ替えに基づく分散分析(permutation ANOVA)と呼ばれ、正規性の仮定なしに群間差を検定できます。
次に、Pythonで並べ替え検定を実装し、$t$ 検定との比較を行いましょう。
Pythonによる実装
正確並べ替え検定の実装
まず、小さなデータで正確並べ替え検定(すべての置換を列挙する方法)を実装します。
import numpy as np
from itertools import combinations
import matplotlib.pyplot as plt
def exact_permutation_test(x, y, alternative="two-sided"):
"""
正確並べ替え検定(平均の差を検定統計量として使用)
Parameters
----------
x : array-like 群1のデータ
y : array-like 群2のデータ
alternative : str 'two-sided', 'greater', 'less'
Returns
-------
obs_stat : float 観測された統計量
p_value : float p値
perm_stats : array 帰無分布の統計量
"""
x = np.asarray(x)
y = np.asarray(y)
n1, n2 = len(x), len(y)
n = n1 + n2
# 観測された検定統計量
obs_stat = np.mean(x) - np.mean(y)
# すべてのデータを結合
combined = np.concatenate([x, y])
# すべての置換について統計量を計算
perm_stats = []
for idx in combinations(range(n), n1):
idx = list(idx)
rest = [i for i in range(n) if i not in idx]
stat = np.mean(combined[idx]) - np.mean(combined[rest])
perm_stats.append(stat)
perm_stats = np.array(perm_stats)
# p値の計算
if alternative == "two-sided":
p_value = np.mean(np.abs(perm_stats) >= np.abs(obs_stat))
elif alternative == "greater":
p_value = np.mean(perm_stats >= obs_stat)
else:
p_value = np.mean(perm_stats <= obs_stat)
return obs_stat, p_value, perm_stats
# --- 臨床試験データの例 ---
drug = np.array([-8, -5, -12, -3, -7]) # 新薬群
placebo = np.array([-2, 1, -4, 2, -1]) # プラセボ群
obs_stat, p_value, perm_stats = exact_permutation_test(drug, placebo)
print("=" * 50)
print("正確並べ替え検定の結果")
print("=" * 50)
print(f"新薬群の平均: {np.mean(drug):.1f}")
print(f"プラセボ群の平均: {np.mean(placebo):.1f}")
print(f"観測された差: {obs_stat:.1f}")
print(f"並べ替えの総数: {len(perm_stats)}")
print(f"p値 (両側): {p_value:.4f}")
print(f"有意水準5%で: {'棄却' if p_value < 0.05 else '棄却しない'}")
# --- 帰無分布の可視化 ---
fig, ax = plt.subplots(figsize=(9, 5))
ax.hist(perm_stats, bins=30, density=True, alpha=0.6, color="steelblue",
edgecolor="white", label="Permutation distribution")
ax.axvline(obs_stat, color="red", linewidth=2, linestyle="--",
label=f"Observed = {obs_stat:.1f}")
ax.axvline(-obs_stat, color="red", linewidth=2, linestyle="--", alpha=0.5)
# 棄却域の着色
mask = np.abs(perm_stats) >= np.abs(obs_stat)
extreme_stats = perm_stats[mask]
ax.hist(extreme_stats, bins=30, density=True, alpha=0.4, color="red",
edgecolor="white", label=f"Extreme values (p={p_value:.4f})")
ax.set_xlabel("Difference in means", fontsize=12)
ax.set_ylabel("Density", fontsize=12)
ax.set_title("Exact permutation test: null distribution", fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig("permutation_test_exact.png", dpi=150, bbox_inches="tight")
plt.show()
上のグラフから、以下のことが読み取れます。
-
帰無分布は離散的でほぼ対称な形状をしている。$\binom{10}{5} = 252$ 通りの置換から得られた平均の差の分布は、ゼロを中心に対称です。これは帰無仮説のもとで2群が区別されないという対称性を反映しています
-
観測された差 $-7.0$(赤い破線)は帰無分布の裾にある。この位置が帰無分布の中でどれだけ「極端」かがp値に対応します。赤く着色された領域が帰無分布に占める割合がp値です
-
正確なp値は252通りの列挙から得られている。漸近近似を使っていないため、小標本でも正確な結論を得ることができます
モンテカルロ並べ替え検定の実装
次に、サンプルサイズが大きい場合に使えるモンテカルロ並べ替え検定を実装します。
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def monte_carlo_permutation_test(x, y, n_perm=9999, stat_func=None,
alternative="two-sided"):
"""
モンテカルロ並べ替え検定
Parameters
----------
x, y : array-like 2つの群のデータ
n_perm : int ランダム置換の回数
stat_func : callable 検定統計量の関数
alternative : str 'two-sided', 'greater', 'less'
Returns
-------
obs_stat, p_value, perm_stats
"""
x = np.asarray(x, dtype=float)
y = np.asarray(y, dtype=float)
n1 = len(x)
combined = np.concatenate([x, y])
n = len(combined)
# デフォルトの統計量: 平均の差
if stat_func is None:
stat_func = lambda a, b: np.mean(a) - np.mean(b)
obs_stat = stat_func(x, y)
# モンテカルロ置換
perm_stats = np.zeros(n_perm)
for i in range(n_perm):
perm = np.random.permutation(combined)
perm_stats[i] = stat_func(perm[:n1], perm[n1:])
# p値
if alternative == "two-sided":
p_value = (np.sum(np.abs(perm_stats) >= np.abs(obs_stat)) + 1) / (n_perm + 1)
elif alternative == "greater":
p_value = (np.sum(perm_stats >= obs_stat) + 1) / (n_perm + 1)
else:
p_value = (np.sum(perm_stats <= obs_stat) + 1) / (n_perm + 1)
return obs_stat, p_value, perm_stats
# --- 大きめのデータで検証 ---
np.random.seed(42)
n1, n2 = 30, 35
x = np.random.exponential(scale=3.0, size=n1) + 1.5 # 指数分布(歪んだ分布)
y = np.random.exponential(scale=3.0, size=n2)
# モンテカルロ並べ替え検定
obs_stat, perm_p, perm_stats = monte_carlo_permutation_test(x, y, n_perm=9999)
# 比較: t検定
t_stat, t_p = stats.ttest_ind(x, y, equal_var=False)
# 比較: マン・ホイットニーU検定
u_stat, u_p = stats.mannwhitneyu(x, y, alternative="two-sided")
print("=" * 50)
print("検定結果の比較(指数分布データ)")
print("=" * 50)
print(f"群1: n={n1}, 平均={np.mean(x):.3f}, 中央値={np.median(x):.3f}")
print(f"群2: n={n2}, 平均={np.mean(y):.3f}, 中央値={np.median(y):.3f}")
print(f"\n並べ替え検定: p = {perm_p:.4f}")
print(f"Welch t検定: p = {t_p:.4f}")
print(f"Mann-Whitney: p = {u_p:.4f}")
# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(13, 5))
# 帰無分布
ax = axes[0]
ax.hist(perm_stats, bins=60, density=True, alpha=0.6, color="steelblue",
edgecolor="white", label="Permutation distribution")
ax.axvline(obs_stat, color="red", linewidth=2, linestyle="--",
label=f"Observed = {obs_stat:.2f}")
# 正規近似
mu_perm, sigma_perm = np.mean(perm_stats), np.std(perm_stats)
x_norm = np.linspace(mu_perm - 4*sigma_perm, mu_perm + 4*sigma_perm, 200)
ax.plot(x_norm, stats.norm.pdf(x_norm, mu_perm, sigma_perm), "g-",
linewidth=1.5, label="Normal approx.")
ax.set_xlabel("Difference in means", fontsize=11)
ax.set_ylabel("Density", fontsize=11)
ax.set_title("Monte Carlo permutation distribution", fontsize=12)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
# データの分布
ax = axes[1]
ax.hist(x, bins=15, density=True, alpha=0.5, color="blue",
edgecolor="white", label=f"Group 1 (n={n1})")
ax.hist(y, bins=15, density=True, alpha=0.5, color="orange",
edgecolor="white", label=f"Group 2 (n={n2})")
ax.set_xlabel("Value", fontsize=11)
ax.set_ylabel("Density", fontsize=11)
ax.set_title("Data distribution (exponential)", fontsize=12)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig("permutation_test_mc.png", dpi=150, bbox_inches="tight")
plt.show()
この結果から、以下の点が読み取れます。
-
並べ替え検定のp値とWelch t検定のp値は近い値を示している。大標本では並べ替え検定とt検定は漸近的に同等になるという理論的予測と整合しています
-
並べ替え分布は正規近似(緑線)によく近似されている。中心極限定理の効果により、$n_1 = 30$, $n_2 = 35$ の標本サイズでは、並べ替え分布はほぼ正規分布です。これが、大標本でt検定と並べ替え検定が同等になる理由です
-
データの分布は明らかに歪んでいる(右図)。指数分布は右に裾の長い歪んだ分布であり、正規分布とは大きく異なります。にもかかわらず、並べ替え検定は分布の仮定なしに正確なp値を提供しています
サイズと検出力のシミュレーション比較
並べ替え検定とt検定のサイズ(第1種の過誤率)と検出力(第2種の過誤率の補数)を、正規分布と歪んだ分布の両方で比較します。
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
np.random.seed(42)
def quick_perm_test(x, y, n_perm=999):
"""高速なモンテカルロ並べ替え検定"""
combined = np.concatenate([x, y])
n1 = len(x)
obs = np.mean(x) - np.mean(y)
count = 1 # 観測値自身を含む
for _ in range(n_perm):
perm = np.random.permutation(combined)
if abs(np.mean(perm[:n1]) - np.mean(perm[n1:])) >= abs(obs):
count += 1
return count / (n_perm + 1)
n1, n2 = 15, 15
n_sim = 2000
alpha = 0.05
effect_sizes = np.linspace(0, 1.5, 10)
results = {
"normal": {"perm": [], "t": []},
"lognormal": {"perm": [], "t": []}
}
for dist_name in ["normal", "lognormal"]:
for delta in effect_sizes:
reject_perm = 0
reject_t = 0
for _ in range(n_sim):
if dist_name == "normal":
x = np.random.randn(n1) + delta
y = np.random.randn(n2)
else:
x = np.random.lognormal(mean=delta * 0.5, sigma=1.0, size=n1)
y = np.random.lognormal(mean=0, sigma=1.0, size=n2)
# 並べ替え検定
p_perm = quick_perm_test(x, y, n_perm=499)
if p_perm < alpha:
reject_perm += 1
# t検定
_, p_t = stats.ttest_ind(x, y, equal_var=False)
if p_t < alpha:
reject_t += 1
results[dist_name]["perm"].append(reject_perm / n_sim)
results[dist_name]["t"].append(reject_t / n_sim)
fig, axes = plt.subplots(1, 2, figsize=(13, 5.5))
for ax, dist_name, title in zip(axes,
["normal", "lognormal"],
["Normal distribution", "Log-normal distribution"]):
ax.plot(effect_sizes, results[dist_name]["perm"], "b-o", markersize=5,
linewidth=1.5, label="Permutation test")
ax.plot(effect_sizes, results[dist_name]["t"], "r-s", markersize=5,
linewidth=1.5, label="Welch t-test")
ax.axhline(alpha, color="gray", linestyle="--", linewidth=0.8,
label=f"Nominal level = {alpha}")
ax.set_xlabel("Effect size", fontsize=12)
ax.set_ylabel("Rejection rate", fontsize=12)
ax.set_title(title, fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
ax.set_ylim(0, 1.05)
plt.tight_layout()
plt.savefig("permutation_vs_t_power.png", dpi=150, bbox_inches="tight")
plt.show()
このシミュレーション結果から、以下の重要な知見が得られます。
-
正規分布の場合(左図): 並べ替え検定とWelch t検定はほぼ同じ性能を示しています。効果サイズ $= 0$(帰無仮説が正しい)では両方とも棄却率が名目水準 $0.05$ に近く、効果サイズが増加するにつれて検出力が同じように上昇しています。正規分布の場合、t検定の仮定が正しいため、並べ替え検定を使う利点はありませんが、損失もありません
-
対数正規分布の場合(右図): 歪んだ分布では、2つの検定の間にいくつかの違いが見られます。効果サイズ $= 0$ でのサイズ(第1種の過誤率)は両方とも $0.05$ に近いですが、検出力の上がり方にわずかな差があります。対数正規分布は強い右歪みを持つため、平均の差を検定統計量とする場合、外れ値の影響を受けやすくなります
-
並べ替え検定は分布によらず安定したサイズを持つ。これが並べ替え検定の最大の利点です。データの分布が何であっても、帰無仮説のもとでの棄却率は名目水準を超えません
並べ替え検定の注意点と限界
交換可能性の仮定の重要性
並べ替え検定は「分布フリー」とよく言われますが、これは「一切の仮定なし」という意味ではありません。帰無仮説のもとでのデータの交換可能性は本質的な仮定です。
以下のような場合、交換可能性が成り立たず、標準的な並べ替え検定が不正確になる可能性があります。
-
帰無仮説のもとでも分散が異なる場合: 例えば、$H_0: \mu_1 = \mu_2$ を検定したいが、$\sigma_1^2 \neq \sigma_2^2$ の場合。平均の差の並べ替え検定はリベラルに(棄却しすぎ)なることがあります。この場合は、スチューデント化した統計量(t統計量)を使うことで問題を軽減できます
-
時系列データ: 時間的な相関があるデータは交換可能ではありません。このような場合は、ブロック並べ替えや循環並べ替え(circular permutation)などの修正が必要です
-
クラスター構造: データがクラスター(学校内の生徒、病院内の患者など)に分かれている場合、クラスター内の相関を考慮した並べ替えが必要です
計算コスト
モンテカルロ並べ替え検定は、$B$ 回のランダム置換について検定統計量を計算する必要があるため、計算コストは $O(Bn)$ です。大規模データで複雑な統計量を使う場合、計算時間が問題になることがあります。
高速化のテクニックとしては、以下のものがあります。
- 十分統計量の利用: 統計量の更新式を利用して、各置換での計算をインクリメンタルに行う
- 早期打ち切り: p値が明らかに大きい(または小さい)場合に、置換の回数を途中で打ち切る
- 並列計算: 各置換は独立に計算できるため、並列処理が容易
ブートストラップ検定との違い
並べ替え検定はしばしばブートストラップ検定と混同されますが、両者は異なる方法です。
| 特徴 | 並べ替え検定 | ブートストラップ検定 |
|---|---|---|
| 再標本の方法 | ラベルの置換(非復元) | 復元抽出 |
| 帰無仮説の構成 | ラベルの交換可能性 | 帰無仮説に合わせたリサンプリング |
| 正確性 | 正確検定(有限標本) | 漸近的に正確 |
| 主な用途 | 仮説検定 | 信頼区間、分散推定 |
並べ替え検定は「帰無仮説のもとでラベルを入れ替えても統計量の分布は変わらない」という性質を利用するのに対し、ブートストラップ検定は「データの再標本化により統計量の分布を推定する」という異なるアプローチを取ります。
まとめ
本記事では、並べ替え検定(置換検定)の理論を直感的な理解から実装まで解説しました。
- 並べ替え検定の核心的アイデアは、帰無仮説のもとでは群のラベルを入れ替えても検定統計量の分布が変わらないという交換可能性を利用することです
- 正確並べ替え検定はすべての可能な置換を列挙して帰無分布を構成するため、分布の仮定なしに有限標本で正確なサイズを保証します
- モンテカルロ並べ替え検定はランダムな置換で帰無分布を近似し、$B = 999 \sim 9{,}999$ 回の置換で実用上十分な精度が得られます
- 並べ替え検定はパラメトリック検定と漸近的に同等であり、分布の仮定が正しい場合でも検出力の損失はほとんどありません。分布の仮定が正しくない場合には、並べ替え検定の方が頑健です
- 注意点として、帰無仮説のもとでのデータの交換可能性は必要な仮定であり、分散が異なる場合や時系列データでは修正が必要です
次のステップとして、以下の記事も参考にしてください。
- マン・ホイットニーU検定の理論 — 順位に基づくノンパラメトリック検定
- ブートストラップ法 — 再標本化法のもう一つの柱
- コルモゴロフ・スミルノフ検定 — 分布全体の比較
- 偽発見率の制御 — 多重検定問題への対処法