SVMのカーネル設計:Mercerの定理・グラム行列の正定値性とRandom Fourier Featuresの理論とPython実装

SVMのカーネルトリックをMercerの定理から掘り下げ、numpy.linalg.eigvalshでグラム行列の正定値性を検証し多項式カーネルの明示的特徴写像との一致を実測。RBFカーネルが無限次元特徴空間に対応することをRandom Fourier Features(Rahimi & Recht)でnumpy実装し近似誤差の収束を確認、sklearn.svm.SVCでハイパーパラメータgammaの過学習・未学習遷移を検証します。

はじめに

https://yuhi-sa.github.io/posts/20260226_svm/1/ ではSVMの基礎理論とカーネルトリックを紹介し、線形・RBF・多項式カーネルを一覧表で示すにとどめていました。しかし「なぜ好き勝手な関数をカーネルとして使えないのか」「RBFカーネルはなぜ無限次元の特徴写像に対応すると言われるのか」という核心には踏み込んでいませんでした。本記事では、Mercerの定理によるカーネルの妥当性判定条件をグラム行列の正定値性から導出し、多項式カーネルの明示的特徴写像を構成して数値的に検証します。さらに、RBFカーネルの無限次元性をRandom Fourier Features(Rahimi & Recht, 2007)という実用的な近似手法で体感し、ハイパーパラメータ \(\gamma\) が汎化性能に与える影響を実験で確認します。

カーネルトリックの復習

https://yuhi-sa.github.io/posts/20260226_svm/1/ で導出したSVMの双対問題では、データは内積 \(\mathbf{x}_i \cdot \mathbf{x}_j\) を通じてのみ現れます。この内積をカーネル関数 \(K(\mathbf{x}_i, \mathbf{x}_j) = \phi(\mathbf{x}_i) \cdot \phi(\mathbf{x}_j)\) に置き換えることで、高次元の特徴写像 \(\phi\) を明示的に計算せずに非線形分類ができます。しかし、任意の2変数関数 \(K(\mathbf{x}, \mathbf{y})\) をこの用途に使えるわけではありません。

Mercerの定理:カーネルの妥当性条件

関数 \(K(\mathbf{x}, \mathbf{y})\) が何らかの特徴写像 \(\phi\) に対して \(K(\mathbf{x}, \mathbf{y}) = \phi(\mathbf{x}) \cdot \phi(\mathbf{y})\) と書ける(=妥当なカーネルである)ための必要十分条件は、任意の有限個のデータ点 \(\{\mathbf{x}_1, \ldots, \mathbf{x}_n\}\) に対して、グラム行列

\[ K_{ij} = K(\mathbf{x}_i, \mathbf{x}_j) \tag{1} \]

が半正定値(すべての固有値が0以上)であることです(Mercerの定理)。直感的には、グラム行列が半正定値であれば、そのグラム行列自体を「特徴ベクトルの内積行列」として解釈できる座標系(固有値分解 \(K = V \Lambda V^\top\) から \(\phi(\mathbf{x}_i) = \sqrt{\Lambda} V_i\) )が必ず存在するためです。

数値実験1:グラム行列の正定値性を直接検証

3次元の乱数データ50点に対し、RBFカーネル・多項式カーネル・そして妥当でないカーネル(\(K(\mathbf{x},\mathbf{y}) = -\|\mathbf{x}-\mathbf{y}\|\) )のグラム行列を計算し、最小固有値を確認しました。

import numpy as np

def rbf_kernel(X, gamma=0.5):
    sq = np.sum(X**2, axis=1)
    D2 = sq[:, None] + sq[None, :] - 2 * X @ X.T
    return np.exp(-gamma * D2)

def poly_kernel(X, degree=2, gamma=1.0, coef0=1.0):
    return (gamma * (X @ X.T) + coef0) ** degree

def invalid_kernel(X):
    sq = np.sum(X**2, axis=1)
    D2 = sq[:, None] + sq[None, :] - 2 * X @ X.T
    return -np.sqrt(np.maximum(D2, 0))  # 負の距離:一般には妥当なカーネルではない


rng = np.random.default_rng(42)
X = rng.normal(size=(50, 3))
K_rbf, K_poly, K_invalid = rbf_kernel(X), poly_kernel(X), invalid_kernel(X)

for name, K in [("RBF", K_rbf), ("Poly(degree2)", K_poly), ("Invalid", K_invalid)]:
    eig = np.linalg.eigvalsh(K)
    print(f"{name}: min_eig={eig.min():.6f} max_eig={eig.max():.6f} n_negative={(eig < -1e-8).sum()}")

実行結果:

RBF: min_eig=0.000181 max_eig=15.644640 n_negative=0
Poly(degree2): min_eig=-0.000000 max_eig=254.878312 n_negative=0
Invalid: min_eig=-98.802080 max_eig=22.587837 n_negative=1
カーネルグラム行列の最小固有値
RBF0.000181 (≥0、妥当)
多項式(次数2)-0.000000 (≥0、妥当)
\(-\|\mathbf{x}-\mathbf{y}\|\) (不正)-98.802080(負、不正)

RBFと多項式カーネルの最小固有値はほぼゼロ以上(浮動小数点誤差の範囲)で半正定値性が確認できました。一方、負の距離をそのままカーネルとして使うと最小固有値が大きく負になり、Mercerの条件を満たさないことが数値的に確認できます(実際、負のユークリッド距離は一般には正定値核ではないことが知られています)。この3つのカーネルの固有値スペクトルを可視化すると、不正なカーネルだけが負の領域に突出したバーを持つことが一目で分かります。

Gram行列の固有値スペクトル:妥当なカーネル(RBF)と不正なカーネルの比較

RBFカーネルの50個の固有値はすべて0以上(最小0.000181)である一方、不正なカーネルは49個が正であるにもかかわらず、たった1個の固有値が**-98.80**という大きな負の値を取っています。この「たった1つの負の固有値」が、SVMの最適化を実際に壊しうることを次の実験で確認します。

数値実験1b:非正定値カーネルは凸最適化ソルバーを実際に破綻させる

SVMの双対問題は、\(H_{ij} = y_i y_j K(\mathbf{x}_i,\mathbf{x}_j)\) を用いて

\[ \max_{\boldsymbol{\alpha}} \sum_i \alpha_i - \frac12 \boldsymbol{\alpha}^\top H \boldsymbol{\alpha} \quad \text{s.t. } 0 \le \alpha_i \le C,\ \sum_i \alpha_i y_i = 0 \]

という二次計画問題(QP)です。この目的関数が凹関数(=最大化問題として解ける)であるためには \(H\) が半正定値である必要があり、\(H\) の半正定値性は \(K\) の半正定値性と同値です(\(H = \text{diag}(y)\,K\,\text{diag}(y)\) で \(\text{diag}(y)\) は直交行列なので固有値は \(K\) と一致します)。つまり、Mercerの条件を満たさないカーネルを使うと、SVMの双対問題はもはや凸最適化問題ではなくなります

これを実際に確かめるため、標準的な内点法QPソルバー cvxopt.solvers.qp(半正定値性を前提にコレスキー分解でKKT系を解く)と、sklearn.svm.SVC(libsvmのSMOという発見的解法を使い、半正定値性を事前にチェックしない)の両方に同じ3つのカーネルを与えて比較しました。

import numpy as np
from sklearn.svm import SVC
from cvxopt import matrix, solvers

solvers.options["show_progress"] = False


def solve_svm_dual_cvxopt(K, y, C=1.0):
    n = len(y)
    P = matrix(np.outer(y, y) * K)
    q = matrix(-np.ones(n))
    G = matrix(np.vstack([-np.eye(n), np.eye(n)]))
    h = matrix(np.hstack([np.zeros(n), np.full(n, C)]))
    A = matrix(y.reshape(1, -1))
    b = matrix(np.zeros(1))
    return solvers.qp(P, q, G, h, A, b)


rng = np.random.default_rng(42)
X = rng.normal(size=(50, 3))
y = (X[:, 0] > 0).astype(float) * 2 - 1  # +-1ラベル
K_rbf, K_poly, K_invalid = rbf_kernel(X), poly_kernel(X), invalid_kernel(X)  # 数値実験1で定義した関数を再利用

for name, K in [("rbf", K_rbf), ("poly", K_poly), ("invalid", K_invalid)]:
    clf = SVC(kernel="precomputed", C=1.0).fit(K, y)
    print(f"{name}: SVC train_acc={clf.score(K, y):.3f}  (警告なしで完走)")
    try:
        sol = solve_svm_dual_cvxopt(K, y)
        print(f"      cvxopt: status={sol['status']}")
    except Exception as e:
        print(f"      cvxopt: EXCEPTION {type(e).__name__}: {e}")
カーネルsklearn.svm.SVC(libsvm/SMO)cvxopt.solvers.qp(内点法)
RBFtrain_acc=0.980, SV=29/50, 警告なしstatus=optimal, dual obj=-14.358952, gap=2.11e-06
多項式(次数2)train_acc=0.980, SV=13/50, 警告なしstatus=optimal, dual obj=-6.463952, gap=2.82e-07
\(-\|\mathbf{x}-\mathbf{y}\|\) (不正)train_acc=0.980, SV=19/50, 警告なしで完走(危険)ValueError: Rank(A) < p or Rank([P; A; G]) < n

結果は対照的でした。RBF・多項式カーネルでは両ソルバーとも問題なく解け、cvxopt は双対ギャップ \(10^{-6}\) 程度で最適解に到達します。しかし不正なカーネルでは、cvxopt は \(H\) が半正定値でないことを検出してKKT系の分解に失敗し、明示的な例外を送出して停止しました。一方 sklearn.svm.SVC は半正定値性を事前チェックしないSMO発見的解法を使っているため、警告もエラーも出さずに「訓練精度98.0%」というもっともらしい数値を返します。しかしこの数値は、目的関数が凹関数であることを前提にした停止条件(KKT条件を満たしたら終了)を、非凹関数に対して形式的に適用しているだけであり、大域最適性はおろか、解の一意性すら理論的に保証されません。つまり「エラーが出ないこと」は「カーネルが妥当であること」を意味しない、という実務上重要な教訓が得られます。カスタムカーネルを設計する際は、必ず事前にグラム行列の半正定値性を検証すべきです。

数値実験2:多項式カーネルの明示的特徴写像

次数2の多項式カーネル \(K(\mathbf{x},\mathbf{y}) = (\gamma\,\mathbf{x}\cdot\mathbf{y} + c)^2\) (\(\mathbf{x},\mathbf{y}\in\mathbb{R}^2\) )を多項定理で展開すると、明示的な6次元特徴写像

\[ \phi(\mathbf{x}) = \left(\gamma x_1^2,\ \gamma x_2^2,\ \sqrt{2}\gamma x_1 x_2,\ \sqrt{2c\gamma}x_1,\ \sqrt{2c\gamma}x_2,\ c\right) \tag{2} \]

が得られます。この写像を使って内積 \(\phi(\mathbf{x})\cdot\phi(\mathbf{y})\) を直接計算し、カーネル関数の値と比較しました。

def explicit_feature_map_deg2(x, gamma, coef0):
    x1, x2 = x
    return np.array([
        gamma * x1**2,
        gamma * x2**2,
        np.sqrt(2) * gamma * x1 * x2,
        np.sqrt(2 * coef0 * gamma) * x1,
        np.sqrt(2 * coef0 * gamma) * x2,
        coef0,
    ])
K(x,y) via kernel formula : 0.0018919198
K(x,y) via <phi(x),phi(y)>: 0.0018919198
abs diff: 3.70e-16

カーネル関数の直接計算と、明示的特徴写像の内積が機械精度で完全一致しました。これは「カーネルトリックが本当に高次元内積のショートカットになっている」ことの直接的な証拠です。

RBFカーネルは無限次元特徴写像に対応する

RBFカーネル \(K(\mathbf{x},\mathbf{y}) = \exp(-\gamma\|\mathbf{x}-\mathbf{y}\|^2)\) は、多項式カーネルのような有限次元の明示的特徴写像を持ちません。これは指数関数をテイラー展開すると

\[ \exp(2\gamma\,\mathbf{x}\cdot\mathbf{y}) = \sum_{k=0}^{\infty} \frac{(2\gamma\,\mathbf{x}\cdot\mathbf{y})^k}{k!} \tag{3} \]

という無限次数の多項式カーネルの正の重み付き和になっているためです(各項 \((\mathbf{x}\cdot\mathbf{y})^k\) はそれ自身が妥当なカーネルであり、正の係数の和もまた妥当なカーネルです)。RBFカーネル自体は \(K(\mathbf{x},\mathbf{y}) = \exp(-\gamma\|\mathbf{x}\|^2)\exp(2\gamma\,\mathbf{x}\cdot\mathbf{y})\exp(-\gamma\|\mathbf{y}\|^2)\) と分解でき、正の対角スケーリングは半正定値性を保つため、RBFカーネル全体も妥当です。

数値実験3:Random Fourier Featuresで無限次元性を体感する

無限次元の特徴写像を厳密に計算することはできませんが、Random Fourier Features(RFF、Rahimi & Recht, 2007)という手法で有限次元にランダム近似できます。Bochnerの定理により、RBFカーネルはランダムな周波数 \(\mathbf{w}\sim\mathcal{N}(0, 2\gamma I)\) とランダムな位相 \(b\sim\text{Uniform}(0,2\pi)\) を使って

\[ K(\mathbf{x},\mathbf{y}) \approx \frac{1}{D}\sum_{i=1}^{D} z_{\mathbf{w}_i,b_i}(\mathbf{x})\, z_{\mathbf{w}_i,b_i}(\mathbf{y}), \qquad z_{\mathbf{w},b}(\mathbf{x}) = \sqrt{2}\cos(\mathbf{w}^\top\mathbf{x}+b) \tag{4} \]

と近似できます(\(D\) 個のランダム特徴を使った \(D\) 次元の明示的特徴写像)。

w = rng.normal(0, np.sqrt(2 * gamma), size=(D, d))
b = rng.uniform(0, 2 * np.pi, size=D)
def rff(X):
    proj = X @ w.T + b
    return np.sqrt(2.0 / D) * np.cos(proj)

5次元データ200組のペアで真のRBFカーネル値と近似値の RMSE を \(D\) を変えて測定しました。

ランダム特徴数 \(D\)RMSE
100.30419
1000.09632
1,0000.03050
10,0000.00966

\(D\) を10倍にするごとにRMSEが約 \(1/\sqrt{10} \approx 0.32\) 倍になっており(モンテカルロ近似の標準的な収束率 \(O(1/\sqrt{D})\) 通り)、有限個のランダム特徴で無限次元カーネルに収束していく様子が定量的に確認できます。RFFはこの原理を使って、カーネルSVMを線形SVM(sklearn.svm.LinearSVCsklearn.linear_model.SGDClassifier)で近似することで、大規模データに対する計算量を \(O(n^2)\) から \(O(n)\) に削減する実用手法(sklearn.kernel_approximation.RBFSampler)としても使われています。

収束率 \(O(1/\sqrt{D})\) を対数対数プロットで定量検証する

上の表は \(D\) を4点しか振っていないため、収束率をより厳密に確認するため \(D\) を10から20,000まで対数的に11段階振り、各 \(D\) でシードを変えた10試行のRMSEの平均・標準偏差を測定しました。

import numpy as np


def true_rbf(x, y, gamma):
    return np.exp(-gamma * np.sum((x - y) ** 2))


def rff_features(X, w, b, D):
    proj = X @ w.T + b
    return np.sqrt(2.0 / D) * np.cos(proj)


gamma, d, n_pairs, n_trials = 0.5, 5, 200, 10
master_rng = np.random.default_rng(0)
Xp = master_rng.normal(size=(n_pairs, d))
Yp = master_rng.normal(size=(n_pairs, d))
K_true = np.array([true_rbf(Xp[i], Yp[i], gamma) for i in range(n_pairs)])

Ds = [10, 20, 50, 100, 200, 500, 1000, 2000, 5000, 10000, 20000]
rmse_means = []
for D in Ds:
    trial_rmses = []
    for t in range(n_trials):
        rng = np.random.default_rng(1000 * D + t)  # (D, 試行) ごとに一意な再現可能シード
        w = rng.normal(0, np.sqrt(2 * gamma), size=(D, d))
        b = rng.uniform(0, 2 * np.pi, size=D)
        ZX, ZY = rff_features(Xp, w, b, D), rff_features(Yp, w, b, D)
        K_approx = np.sum(ZX * ZY, axis=1)
        trial_rmses.append(np.sqrt(np.mean((K_approx - K_true) ** 2)))
    rmse_means.append(np.mean(trial_rmses))

slope, intercept = np.polyfit(np.log10(Ds), np.log10(rmse_means), 1)
print(f"fitted slope: {slope:.4f}  (theory: -0.5)")
fitted slope: -0.4998  (theory: -0.5)

Random Fourier Features近似誤差の収束(対数対数プロット)

対数対数空間で最小二乗フィットした傾きは -0.4998 となり、理論値 \(-1/2\) (\(O(D^{-1/2})\) )と小数点3桁まで一致しました。11点のデータが対数対数プロット上でほぼ完全な直線に乗っており、\(D\) が2桁変化する広い範囲でモンテカルロ収束率の理論がそのまま成り立つことが視覚的にも定量的にも確認できます。

カーネルの構成則:カスタムカーネルの作り方

妥当なカーネルは以下の演算で閉じています(新しいカーネルを作っても妥当性が保たれます)。

  • 正のスカラー倍:\(K\) が妥当なら \(cK\) (\(c>0\) )も妥当
  • :\(K_1, K_2\) が妥当なら \(K_1 + K_2\) も妥当
  • :\(K_1, K_2\) が妥当なら \(K_1 \cdot K_2\) (要素積)も妥当

これを使い、RBFカーネルと多項式カーネルの単純な重み付き和 \(0.5K_{\text{RBF}} + 0.5K_{\text{poly}}\) のグラム行列を検証すると、最小固有値は \(0.000045\) (ほぼゼロ以上)となり、妥当なカーネルであることが確認できました。この構成則は、時系列データに「周期性を捉えるカーネル + トレンドを捉えるカーネル」を組み合わせるなど、ドメイン知識に応じたカスタムカーネル設計の理論的根拠になります(https://yuhi-sa.github.io/posts/20260502_gaussian_process/1/ のガウス過程回帰でも同じ構成則が使われます)。

ハイパーパラメータ γ の効果:未学習から過学習まで

sklearn.datasets.make_moons(ノイズ0.25)で生成した非線形データに対し、\(C=1.0\) を固定して \(\gamma\) を変化させ、sklearn.svm.SVC で訓練・テスト精度とサポートベクター数を測定しました。

\(\gamma\)訓練精度テスト精度サポートベクター数
0.010.8100.833147 / 210
0.10.8380.87895 / 210
1.00.9520.97861 / 210
10.00.9570.98994 / 210
100.00.9810.933189 / 210
1000.01.0000.556210 / 210 (全点)

\(\gamma\) が小さい(0.01)と決定境界が滑らかすぎて訓練データにすら十分適合できない未学習、\(\gamma=10\) 付近でテスト精度が最大(98.9%)となる最適点、\(\gamma=1000\) では訓練精度100%にもかかわらずテスト精度が55.6%まで崩壊する過学習(全210点がサポートベクターになり、各データ点の周りに孤立した鋭いガウス山を作って丸暗記している状態)が明確に観察できます。RBFカーネルの \(\gamma\) は「特徴空間での類似度の減衰の速さ」を制御しており、大きすぎると各点が独立した特徴になり汎化しなくなることが、この実験で定量的に裏付けられます。

関連記事

おすすめ書籍

パターン認識と機械学習 上(C.M.ビショップ、丸善出版)

カーネル法・サポートベクターマシンを含む機械学習の理論的背景を数式で丁寧に解説した定番書です。

※ 上記は Amazon アソシエイトのリンクです。

参考文献

  • Vapnik, V. N. (1998). Statistical Learning Theory. Wiley.
  • Schölkopf, B., & Smola, A. J. (2002). Learning with Kernels. MIT Press.
  • Rahimi, A., & Recht, B. (2007). Random features for large-scale kernel machines. Advances in Neural Information Processing Systems (NeurIPS).
  • Bishop, C. M. (2006). Pattern Recognition and Machine Learning. Springer. Chapter 6 (Kernel Methods).
  • He, M., He, F., Liu, F., & Huang, X. (2024). Random Fourier features for asymmetric kernels. Machine Learning, 113, 8459–8485. 非対称カーネル(正定値でも半正定値でもないカーネルを含む)へRFFを一般化し、PDカーネルと不定値カーネルの両方に適用可能な一様収束保証を与えています。
  • Ma, Z., Yang, J., & Yang, Y. (2025). On the generalization properties of learning the random feature models with learnable activation functions. arXiv:2510.15327. 学習可能な活性化関数を持つランダム特徴モデルの汎化性能を解析し、リーク重み付きサンプリングにより必要な特徴数 \(D\) の理論的下界を大幅に改善できることを示しています。