多出力ガウス過程(Multi-Output Gaussian Process, MOGP)は、複数の出力を同時に予測するガウス過程の拡張です。通常のガウス過程回帰(GPR)は1つの出力に対する回帰ですが、MOGPでは出力間の相関構造を利用することで、各出力の予測精度を向上させることができます。
例えば、複数のセンサーからの信号を同時に予測する場合や、多変量時系列の回帰問題において、出力間の相関をモデルに取り込むことが有効です。
本記事の内容
- MOGPの定義と動機
- ICM(Intrinsic Coregionalization Model)カーネル
- 予測分布の導出
- Pythonでの実装と可視化
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
単出力ガウス過程の復習
通常のガウス過程回帰では、入力 $\bm{x}$ に対して1つのスカラー出力 $f(\bm{x})$ を予測します。
$$ f(\bm{x}) \sim \mathcal{GP}(0, k(\bm{x}, \bm{x}’)) $$
予測分布は以下で与えられます。
$$ \bar{f}_* = \bm{k}_*^T(\bm{K} + \sigma_n^2\bm{I})^{-1}\bm{y} $$
多出力ガウス過程(MOGP)とは
$D$ 個の出力 $f_1(\bm{x}), f_2(\bm{x}), \dots, f_D(\bm{x})$ を同時にモデル化する場合を考えます。MOGPでは、異なる出力間のカーネル関数 $k_{dd’}(\bm{x}, \bm{x}’)$ を定義します。
$$ \text{Cov}[f_d(\bm{x}), f_{d’}(\bm{x}’)] = k_{dd’}(\bm{x}, \bm{x}’) $$
全ての出力を1つのベクトルとして扱うことで、大きなガウス過程として定式化できます。
ICM(Intrinsic Coregionalization Model)
出力間のカーネルを構成する代表的な方法がICMです。ICMでは、カーネルを「入力空間のカーネル」と「出力空間の相関行列」のクロネッカー積で表現します。
$$ k_{dd’}(\bm{x}, \bm{x}’) = B_{dd’} \cdot k(\bm{x}, \bm{x}’) $$
ここで、
- $k(\bm{x}, \bm{x}’)$: 入力空間のベースカーネル（例: RBFカーネル）
- $B_{dd’}$: コリジョナリゼーション行列 $\bm{B} \in \mathbb{R}^{D \times D}$ の $(d, d’)$ 要素
$\bm{B}$ は正定値対称行列で、出力間の相関構造を表します。$B_{dd’}$ が大きいほど、出力 $d$ と $d’$ は強い正の相関を持ちます。
全データに対するカーネル行列は、
$$ \bm{K}_{\text{full}} = \bm{B} \otimes \bm{K}_x $$
と表されます。ここで $\otimes$ はクロネッカー積、$\bm{K}_x$ は入力データ間のカーネル行列です。
予測分布
訓練データのラベルベクトルを $\bm{y} = (\bm{y}_1^T, \bm{y}_2^T, \dots, \bm{y}_D^T)^T$ とし、テスト入力 $\bm{x}_*$ に対する全出力の予測分布は、
$$ \bar{\bm{f}}_* = \bm{K}_{*}^T (\bm{K}_{\text{full}} + \sigma_n^2 \bm{I})^{-1} \bm{y} $$
$$ \text{Var}(\bm{f}_*) = \bm{K}_{**} – \bm{K}_{*}^T (\bm{K}_{\text{full}} + \sigma_n^2 \bm{I})^{-1} \bm{K}_{*} $$
Pythonでの実装
2出力のMOGPをICMカーネルでスクラッチ実装します。
import numpy as np
import matplotlib.pyplot as plt
np.random.seed(42)
# --- カーネル関数 ---
def rbf_kernel(X1, X2, length_scale=1.0, variance=1.0):
"""RBFカーネル行列"""
sq_dist = np.sum(X1**2, axis=1).reshape(-1, 1) + \
np.sum(X2**2, axis=1).reshape(1, -1) - \
2 * X1 @ X2.T
return variance * np.exp(-sq_dist / (2 * length_scale**2))
def icm_kernel(X1, X2, B, length_scale=1.0, variance=1.0):
"""ICMカーネル: B (otimes) K_x"""
K_x = rbf_kernel(X1, X2, length_scale, variance)
return np.kron(B, K_x)
# --- MOGP ---
class MOGP:
def __init__(self, B, length_scale=1.0, variance=1.0, noise=0.1):
self.B = B
self.length_scale = length_scale
self.variance = variance
self.noise = noise
self.D = B.shape[0]
def fit(self, X, Y):
"""
X: (n, d) 入力
Y: (n, D) 出力（D個の出力）
"""
self.X_train = X
self.n = X.shape[0]
# ラベルを(n*D, )に並べ替え: [y1_1,...,y1_n, y2_1,...,y2_n]
self.y_train = Y.T.ravel()
# フルカーネル行列
K_full = icm_kernel(X, X, self.B, self.length_scale, self.variance)
K_full += self.noise**2 * np.eye(self.n * self.D)
self.L = np.linalg.cholesky(K_full)
self.alpha = np.linalg.solve(self.L.T,
np.linalg.solve(self.L, self.y_train))
return self
def predict(self, X_test):
"""テスト点での予測平均と分散"""
n_test = X_test.shape[0]
# テスト-訓練間のカーネル
K_star = icm_kernel(X_test, self.X_train, self.B,
self.length_scale, self.variance)
# テスト-テスト間のカーネル
K_ss = icm_kernel(X_test, X_test, self.B,
self.length_scale, self.variance)
# 予測平均
mu = K_star @ self.alpha
# 予測分散
v = np.linalg.solve(self.L, K_star.T)
var = np.diag(K_ss - v.T @ v)
# (n_test*D,) -> (D, n_test) に整形
mu = mu.reshape(self.D, n_test)
var = var.reshape(self.D, n_test)
return mu, var
# --- データ生成 ---
n_train = 20
X_train = np.sort(np.random.uniform(0, 10, n_train)).reshape(-1, 1)
# 2つの相関のある出力
f1_true = np.sin(X_train.ravel())
f2_true = 0.8 * np.sin(X_train.ravel()) + 0.5 * np.cos(X_train.ravel())
Y_train = np.column_stack([
f1_true + np.random.randn(n_train) * 0.1,
f2_true + np.random.randn(n_train) * 0.1
])
# コリジョナリゼーション行列（出力間の相関を表現）
B = np.array([[1.0, 0.8],
[0.8, 1.0]])
# --- MOGP の学習と予測 ---
mogp = MOGP(B, length_scale=1.5, variance=1.0, noise=0.1)
mogp.fit(X_train, Y_train)
X_test = np.linspace(-0.5, 10.5, 100).reshape(-1, 1)
mu_pred, var_pred = mogp.predict(X_test)
# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
colors = ['blue', 'orange']
labels = ['Output 1', 'Output 2']
true_funcs = [
np.sin(X_test.ravel()),
0.8 * np.sin(X_test.ravel()) + 0.5 * np.cos(X_test.ravel())
]
for d in range(2):
ax = axes[d]
std = np.sqrt(np.clip(var_pred[d], 0, None))
ax.plot(X_test, true_funcs[d], 'k--', alpha=0.5, label='True')
ax.plot(X_test, mu_pred[d], color=colors[d], label='MOGP mean')
ax.fill_between(X_test.ravel(),
mu_pred[d] - 2 * std,
mu_pred[d] + 2 * std,
alpha=0.2, color=colors[d], label='95% CI')
ax.scatter(X_train, Y_train[:, d], c='red', s=30, zorder=5, label='Data')
ax.set_xlabel('x')
ax.set_ylabel('f(x)')
ax.set_title(f'{labels[d]}')
ax.legend(fontsize=8)
ax.grid(True, alpha=0.3)
plt.suptitle('Multi-Output Gaussian Process (ICM Kernel)', fontsize=14)
plt.tight_layout()
plt.show()
print(f"コリジョナリゼーション行列 B:\n{B}")
print(f"出力間の相関: {B[0, 1]:.2f}")
2つの出力間に正の相関 ($B_{12} = 0.8$) を設定しているため、出力1のデータが出力2の予測にも寄与し、相互に予測精度が向上しています。
まとめ
本記事では、多出力ガウス過程(MOGP)の理論と実装について解説しました。
- MOGPは複数の出力間の相関構造をモデル化し、出力を同時に予測するガウス過程の拡張
- ICMカーネルは入力カーネルと出力間の相関行列のクロネッカー積で構成される
- コリジョナリゼーション行列 $\bm{B}$ が出力間の相関の強さを表現する
- 出力間に正の相関がある場合、MOGPは各出力を独立に予測するよりも高い精度を実現できる
次のステップとして、以下の記事も参考にしてください。