カーネルリッジ回帰の理論と実装をわかりやすく解説

カーネルリッジ回帰(Kernel Ridge Regression)は、リッジ回帰にカーネルトリックを適用することで、非線形な回帰問題を解く手法です。明示的に高次元の特徴量空間を構成することなく、カーネル関数を通じて非線形な関係を捉えることができます。

本記事では、リッジ回帰の双対表現の導出からカーネルトリックの適用、そしてPythonでの実装までを丁寧に解説します。

本記事の内容

  • リッジ回帰の双対表現の導出
  • カーネルトリックの適用
  • 代表的なカーネル関数
  • Pythonでのスクラッチ実装と可視化

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

リッジ回帰の復習

リッジ回帰の解析解は以下の通りです。

$$ \hat{\bm{w}} = (\bm{X}^T\bm{X} + \lambda \bm{I}_d)^{-1}\bm{X}^T\bm{y} $$

ここで $\bm{X} \in \mathbb{R}^{n \times d}$ は入力データ行列、$\bm{y} \in \mathbb{R}^n$ は目的変数ベクトル、$d$ は特徴量の次元、$n$ はサンプル数です。

双対表現の導出

リッジ回帰の解を別の形で表現します。行列の公式(Woodburyの恒等式)を利用すると、

$$ (\bm{X}^T\bm{X} + \lambda \bm{I}_d)^{-1}\bm{X}^T = \bm{X}^T(\bm{X}\bm{X}^T + \lambda \bm{I}_n)^{-1} $$

が成り立ちます。これを用いると、

$$ \hat{\bm{w}} = \bm{X}^T(\bm{X}\bm{X}^T + \lambda \bm{I}_n)^{-1}\bm{y} $$

新しいデータ $\bm{x}_*$ に対する予測は、

$$ \begin{align} \hat{y}_* &= \bm{x}_*^T \hat{\bm{w}} \\ &= \bm{x}_*^T \bm{X}^T(\bm{X}\bm{X}^T + \lambda \bm{I}_n)^{-1}\bm{y} \\ &= \bm{k}_*^T (\bm{K} + \lambda \bm{I}_n)^{-1}\bm{y} \end{align} $$

ここで $\bm{K} = \bm{X}\bm{X}^T$ はグラム行列で、$K_{ij} = \bm{x}_i^T \bm{x}_j$(内積)です。$\bm{k}_* = \bm{X}\bm{x}_*$ はテスト点と各訓練データの内積ベクトルです。

重要な点は、予測値が内積 $\bm{x}_i^T \bm{x}_j$ のみで表現されていることです。

カーネルトリックの適用

内積を一般的なカーネル関数 $k(\bm{x}_i, \bm{x}_j)$ に置き換えます。

$$ K_{ij} = k(\bm{x}_i, \bm{x}_j) $$

$$ \hat{y}_* = \bm{k}_*^T (\bm{K} + \lambda \bm{I}_n)^{-1}\bm{y} $$

ここで $k_{*i} = k(\bm{x}_*, \bm{x}_i)$ です。カーネル関数を適切に選ぶことで、暗黙的に高次元の特徴空間での線形回帰を行うことができます。

代表的なカーネル関数

線形カーネル:

$$ k(\bm{x}, \bm{x}’) = \bm{x}^T\bm{x}’ $$

多項式カーネル:

$$ k(\bm{x}, \bm{x}’) = (\bm{x}^T\bm{x}’ + c)^d $$

RBFカーネル(ガウシアンカーネル):

$$ k(\bm{x}, \bm{x}’) = \exp\left(-\frac{\|\bm{x} – \bm{x}’\|^2}{2\sigma^2}\right) $$

RBFカーネルは無限次元の特徴空間に対応し、十分に滑らかな任意の関数を近似できます。

Pythonでの実装

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)

# --- カーネル関数 ---
def rbf_kernel(X1, X2, sigma=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 np.exp(-sq_dist / (2 * sigma**2))

def polynomial_kernel(X1, X2, degree=3, c=1.0):
    """多項式カーネル行列"""
    return (X1 @ X2.T + c) ** degree

# --- カーネルリッジ回帰 ---
class KernelRidgeRegression:
    def __init__(self, kernel_func, lam=1.0, **kernel_params):
        self.kernel_func = kernel_func
        self.lam = lam
        self.kernel_params = kernel_params

    def fit(self, X, y):
        self.X_train = X
        K = self.kernel_func(X, X, **self.kernel_params)
        n = K.shape[0]
        self.alpha = np.linalg.solve(K + self.lam * np.eye(n), y)
        return self

    def predict(self, X):
        K_star = self.kernel_func(X, self.X_train, **self.kernel_params)
        return K_star @ self.alpha

# --- データ生成 ---
n_train = 30
X_train = np.sort(np.random.uniform(-3, 3, n_train)).reshape(-1, 1)
y_true_func = lambda x: np.sin(x) + 0.5 * np.cos(2 * x)
y_train = y_true_func(X_train.ravel()) + np.random.randn(n_train) * 0.2

X_test = np.linspace(-3.5, 3.5, 200).reshape(-1, 1)
y_test_true = y_true_func(X_test.ravel())

# --- 異なるカーネルとパラメータでの比較 ---
fig, axes = plt.subplots(2, 2, figsize=(12, 10))

# RBFカーネル: sigma の影響
for i, sigma in enumerate([0.3, 1.0]):
    ax = axes[0, i]
    model = KernelRidgeRegression(rbf_kernel, lam=0.01, sigma=sigma)
    model.fit(X_train, y_train)
    y_pred = model.predict(X_test)

    ax.scatter(X_train, y_train, c='red', s=30, zorder=5, label='Data')
    ax.plot(X_test, y_test_true, 'k--', alpha=0.5, label='True')
    ax.plot(X_test, y_pred, 'b-', label='Prediction')
    ax.set_title(f'RBF Kernel ($\\sigma$={sigma}, $\\lambda$=0.01)')
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.legend(fontsize=8)
    ax.grid(True, alpha=0.3)

# RBFカーネル: lambda の影響
for i, lam in enumerate([0.001, 1.0]):
    ax = axes[1, i]
    model = KernelRidgeRegression(rbf_kernel, lam=lam, sigma=1.0)
    model.fit(X_train, y_train)
    y_pred = model.predict(X_test)

    ax.scatter(X_train, y_train, c='red', s=30, zorder=5, label='Data')
    ax.plot(X_test, y_test_true, 'k--', alpha=0.5, label='True')
    ax.plot(X_test, y_pred, 'b-', label='Prediction')
    ax.set_title(f'RBF Kernel ($\\sigma$=1.0, $\\lambda$={lam})')
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.legend(fontsize=8)
    ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# --- 多項式カーネルとの比較 ---
fig, axes = plt.subplots(1, 2, figsize=(12, 5))

# 多項式カーネル
model_poly = KernelRidgeRegression(polynomial_kernel, lam=0.01, degree=5, c=1.0)
model_poly.fit(X_train, y_train)
y_pred_poly = model_poly.predict(X_test)

axes[0].scatter(X_train, y_train, c='red', s=30, zorder=5, label='Data')
axes[0].plot(X_test, y_test_true, 'k--', alpha=0.5, label='True')
axes[0].plot(X_test, y_pred_poly, 'b-', label='Poly(d=5)')
axes[0].set_title('Polynomial Kernel')
axes[0].legend()
axes[0].grid(True, alpha=0.3)

# RBFカーネル
model_rbf = KernelRidgeRegression(rbf_kernel, lam=0.01, sigma=1.0)
model_rbf.fit(X_train, y_train)
y_pred_rbf = model_rbf.predict(X_test)

axes[1].scatter(X_train, y_train, c='red', s=30, zorder=5, label='Data')
axes[1].plot(X_test, y_test_true, 'k--', alpha=0.5, label='True')
axes[1].plot(X_test, y_pred_rbf, 'b-', label='RBF')
axes[1].set_title('RBF Kernel')
axes[1].legend()
axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

RBFカーネルの幅 $\sigma$ が小さいと局所的にフィットし、大きいと滑らかな曲線になります。正則化パラメータ $\lambda$ は過学習の防止に寄与します。

まとめ

本記事では、カーネルリッジ回帰の理論と実装について解説しました。

  • リッジ回帰の双対表現により、予測値は内積(カーネル関数)のみで計算可能
  • カーネルトリックにより、暗黙的に高次元特徴空間での線形回帰を実現
  • RBFカーネルのパラメータ $\sigma$ と正則化パラメータ $\lambda$ が予測の滑らかさを制御
  • ガウス過程回帰と密接な関係があり、GPRの予測平均はカーネルリッジ回帰の解と一致する

次のステップとして、以下の記事も参考にしてください。