ホテリングの $T^2$ 統計量は、多変量データに対する異常検知手法の中で最も基本的かつ重要なものの1つです。正常データの分布を多変量正規分布で仮定し、その分布の中心からの「距離」が大きいデータ点を異常として検出します。
本手法はシンプルでありながら強力で、製造業の品質管理や宇宙機のヘルスモニタリングなど、幅広い分野で活用されています。
本記事の内容
- ホテリングの $T^2$ 統計量の定義と導出
- マハラノビス距離との関係
- 異常判定の閾値設定（F分布による方法）
- Pythonでのスクラッチ実装と可視化
前提知識
この記事を読む前に、以下の概念を押さえておくと理解が深まります。
- 正規分布とは？定義やグラフの見方を徹底解説
- 多変量正規分布の確率密度関数
- 共分散行列の意味
ホテリングの $T^2$ 統計量とは
ホテリングの $T^2$ 統計量は、多変量データの平均ベクトルからの「統計的な距離」を測る指標です。
大雑把に言うと、各変数の分散や変数間の相関を考慮した上で、データ点がどれくらい「普通でない」かを数値化したものです。単純なユークリッド距離では変数間のスケール差や相関を考慮できませんが、ホテリングの $T^2$ はこれらを適切に扱います。
数学的な定義
$p$ 次元の正常データ $\bm{x}_1, \bm{x}_2, \dots, \bm{x}_n$ が多変量正規分布 $\mathcal{N}(\bm{\mu}, \bm{\Sigma})$ に従うとします。
平均ベクトルと共分散行列を標本から推定すると、
$$ \hat{\bm{\mu}} = \bar{\bm{x}} = \frac{1}{n} \sum_{i=1}^{n} \bm{x}_i $$
$$ \hat{\bm{\Sigma}} = \bm{S} = \frac{1}{n-1} \sum_{i=1}^{n} (\bm{x}_i – \bar{\bm{x}})(\bm{x}_i – \bar{\bm{x}})^T $$
このとき、新しいデータ点 $\bm{x}_{\text{new}}$ に対するホテリングの $T^2$ 統計量は以下で定義されます。
$$ T^2 = (\bm{x}_{\text{new}} – \bar{\bm{x}})^T \bm{S}^{-1} (\bm{x}_{\text{new}} – \bar{\bm{x}}) $$
マハラノビス距離との関係
ホテリングの $T^2$ 統計量は、マハラノビス距離の二乗と本質的に同じものです。マハラノビス距離 $D_M$ は以下で定義されます。
$$ D_M(\bm{x}) = \sqrt{(\bm{x} – \bm{\mu})^T \bm{\Sigma}^{-1} (\bm{x} – \bm{\mu})} $$
したがって、
$$ T^2 = D_M^2 $$
という関係が成り立ちます。マハラノビス距離は共分散行列の逆行列を用いることで、変数間の相関構造を考慮した距離を計算しています。
導出: なぜ共分散行列の逆行列を使うのか
多変量正規分布の確率密度関数を思い出しましょう。
$$ p(\bm{x}) = \frac{1}{(2\pi)^{p/2} |\bm{\Sigma}|^{1/2}} \exp\left( -\frac{1}{2} (\bm{x} – \bm{\mu})^T \bm{\Sigma}^{-1} (\bm{x} – \bm{\mu}) \right) $$
対数を取ると、
$$ \begin{align} \ln p(\bm{x}) &= -\frac{p}{2}\ln(2\pi) – \frac{1}{2}\ln|\bm{\Sigma}| – \frac{1}{2}(\bm{x} – \bm{\mu})^T \bm{\Sigma}^{-1} (\bm{x} – \bm{\mu}) \end{align} $$
$\bm{\mu}$ と $\bm{\Sigma}$ が固定のとき、対数尤度の値は $(\bm{x} – \bm{\mu})^T \bm{\Sigma}^{-1} (\bm{x} – \bm{\mu})$ のみに依存します。この二次形式が大きいほど、データ点 $\bm{x}$ の確率密度は小さく、すなわち「異常」であると判断できます。
閾値の設定方法
異常を判定するためには、$T^2$ の閾値 $T^2_{\alpha}$ を設定する必要があります。正常データが $p$ 次元の多変量正規分布に従うとき、$T^2$ 統計量はF分布に従います。
訓練データに含まれないテストデータに対して、
$$ \frac{n-p}{p(n+1)} T^2 \sim F(p, n-p) $$
が成り立ちます。ここで $n$ は訓練データのサンプル数、$p$ は次元数です。
有意水準 $\alpha$（例えば $\alpha = 0.01$）に対する閾値は、
$$ T^2_{\alpha} = \frac{p(n+1)}{n-p} F_{\alpha}(p, n-p) $$
で計算されます。$T^2 > T^2_{\alpha}$ のデータ点を異常と判定します。
Pythonでの実装
ホテリングの $T^2$ 統計量を用いた異常検知をPythonでスクラッチ実装します。
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
np.random.seed(42)
# 正常データの生成（2次元多変量正規分布）
n_normal = 200
mean_normal = np.array([2.0, 3.0])
cov_normal = np.array([[1.0, 0.5],
[0.5, 1.5]])
X_normal = np.random.multivariate_normal(mean_normal, cov_normal, n_normal)
# 異常データの生成
n_anomaly = 10
X_anomaly = np.random.multivariate_normal([6.0, 7.0], np.eye(2) * 0.3, n_anomaly)
# 全データ
X_all = np.vstack([X_normal, X_anomaly])
labels = np.array([0] * n_normal + [1] * n_anomaly) # 0:正常, 1:異常
# --- ホテリングのT^2統計量の計算 ---
# 正常データから平均と共分散行列を推定
mu_hat = X_normal.mean(axis=0)
S_hat = np.cov(X_normal.T)
S_inv = np.linalg.inv(S_hat)
# 全データに対してT^2を計算
def hotelling_t2(X, mu, S_inv):
"""ホテリングのT^2統計量を計算"""
diff = X - mu
t2 = np.array([d @ S_inv @ d for d in diff])
return t2
t2_scores = hotelling_t2(X_all, mu_hat, S_inv)
# 閾値の設定（F分布に基づく）
p = X_normal.shape[1] # 次元数
n = n_normal
alpha = 0.01
f_critical = stats.f.ppf(1 - alpha, p, n - p)
threshold = p * (n + 1) / (n - p) * f_critical
print(f"次元数 p = {p}, 訓練データ数 n = {n}")
print(f"有意水準 alpha = {alpha}")
print(f"F臨界値 = {f_critical:.4f}")
print(f"T^2 閾値 = {threshold:.4f}")
# 異常判定
predictions = (t2_scores > threshold).astype(int)
print(f"\n異常と判定されたデータ数: {predictions.sum()}")
print(f"実際の異常データ数: {labels.sum()}")
# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
# 散布図
ax1 = axes[0]
normal_mask = predictions == 0
anomaly_mask = predictions == 1
ax1.scatter(X_all[normal_mask, 0], X_all[normal_mask, 1],
c='blue', alpha=0.5, label='Normal', s=20)
ax1.scatter(X_all[anomaly_mask, 0], X_all[anomaly_mask, 1],
c='red', alpha=0.8, label='Anomaly', s=50, marker='x')
ax1.set_xlabel('$x_1$')
ax1.set_ylabel('$x_2$')
ax1.set_title("Hotelling's $T^2$ Anomaly Detection")
ax1.legend()
ax1.grid(True, alpha=0.3)
# T^2スコアのプロット
ax2 = axes[1]
ax2.bar(range(len(t2_scores)), t2_scores, color=['blue' if l == 0 else 'red' for l in labels], alpha=0.6)
ax2.axhline(y=threshold, color='green', linestyle='--', linewidth=2, label=f'Threshold = {threshold:.2f}')
ax2.set_xlabel('Sample Index')
ax2.set_ylabel('$T^2$ Score')
ax2.set_title('$T^2$ Scores with Threshold')
ax2.legend()
ax2.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
上記のコードでは、2次元の正常データから平均ベクトルと共分散行列を推定し、テストデータに対して $T^2$ 統計量を計算しています。F分布に基づく閾値を超えたデータ点を異常として判定しています。
まとめ
本記事では、ホテリングの $T^2$ 統計量を用いた異常検知について解説しました。
- ホテリングの $T^2$ は、多変量正規分布を仮定した上での統計的な距離指標であり、マハラノビス距離の二乗と等価
- 共分散行列の逆行列を用いることで、変数間の相関を考慮した異常度を計算できる
- 閾値はF分布を利用して統計的に設定できる
- シンプルな手法だが、正常データが多変量正規分布に近い場合に高い性能を発揮する
次のステップとして、以下の記事も参考にしてください。