Tunny Icon
TunnyDocs

The next-gen Grasshopper optimization tool.

ガウス過程

ガウス過程回帰（Gaussian Process regression、別名 Kriging）は、ガウス過程（GP）を使った非パラメトリック回帰手法です。 訓練点の間を確率論的に補間し、予測の不確実性も同時に定量化できます。 予測が 1 つの値ではなく分布として得られるため、「この領域はデータが少ないので予測が当てにならない」ということまで曲面から読み取れます。 Tunny Dashboard は egobox-gp クレート（v0.36、Apache-2.0）を使用しており、SMT サロゲートモデリングツールキットの Rust 移植版です。

2 つのスパース GP バリアントが選択できます: GP-FITC と GP-VFE。 どちらもアーキテクチャ、ハイパーパラメータ探索、計算量はまったく同一で、学習時に使う周辺尤度の下界の種類のみが異なります（FITC と VFE の違いを参照）。

概要

GP-FITC と GP-VFE はどちらも egobox の Sparse GP（FITC または VFE 近似）を内部で使用します。

オプション 近似手法 誘導点数 M
GP-FITC FITC min(N, 100) — 下記参照
GP-VFE VFE min(N, 100) — 下記参照

誘導点の選択: N ≤ 100 のとき誘導点 = 訓練点（Z = X）となり、FITC/VFE はノイズ推定付きの厳密 GP と数学的に等価になります。 N > 100 のとき M = 100 個の誘導点を k-means セントロイド（決定論的シード）で選択します。 セントロイドを使うことで、訓練点の分布に偏りがあっても誘導点が空間全体に行き渡ります。 モデルは全 N 点で学習し、サブサンプリングは行いません。 多目的サロゲート最適化では、誘導点をパレートフロント側に偏らせ、非劣 trial に集中させることで、フロントを改善する領域でサロゲートがもっとも正確になるようにできます。

実測学習時間（N = 10,000）: GP-FITC ≈ 2.4 秒、GP-VFE ≈ 2.0 秒（release ビルド）。 これらは開発環境での計測による参考値であり、環境やデータによって変動します。

ガウス過程の基礎

定義

ガウス過程は、任意の有限点集合 {x1,…,xN}{x_1, \ldots, x_N} 上の関数値が多変量正規分布に従うことを仮定する確率過程です:

f(x)∼GP(m(x),k(x,x′))f(x) \sim \mathcal{GP}(m(x), k(x,x'))

  • m(x)m(x): 平均関数（Tunny Dashboard では m(x)=0m(x)=0 を仮定）
  • k(x,x′)k(x, x'): カーネル関数（データ点間の類似度を表す共分散関数）

予測（事後分布）

訓練データ (X,y)(X, y) が与えられたとき、新しい点 x∗x^* に対する事後分布は:

μ(x∗)=k(x∗,Z) w,w=KZZ−1KXZ⊤(Q+Λ)−1y\mu(x^) = k(x^,Z),w, \qquad w = K_{ZZ}^{-1}K_{XZ}^\top(Q+\Lambda)^{-1}y

σ2(x∗)=k(x∗,x∗)−k(x∗,Z)KZZ−1k(Z,x∗)\sigma^2(x^) = k(x^,x^) - k(x^,Z)K_{ZZ}^{-1}k(Z,x^*)

予測（平均、分散）はバッチで計算されます。 95% 信頼帯 = 平均 ± 1.96·σ です。 カーネルは距離が離れるほど値が小さくなるため、訓練点から遠い領域ほど σ2(x∗)\sigma^2(x^*) が大きく、信頼帯も広く出ます。

ARD Matérn 5/2 カーネル

カーネル関数

k(x1,x2)=σf2(1+5 r+5r23)exp⁡(−5 r)k(x_1,x_2) = \sigma_f^2\left(1+\sqrt{5},r+\frac{5r^2}{3}\right)\exp(-\sqrt{5},r)

r2=∑d(x1,d−x2,dld)2r^2 = \sum_d \left(\frac{x_{1,d}-x_{2,d}}{l_d}\right)^2

パラメータ 意味
σfσ_f シグナル標準偏差（関数の振幅スケール）
ldl_d 次元 d の長さスケール（大きいほど滑らか）
σnσ_n 観測ノイズ標準偏差

ARD（Automatic Relevance Determination）: 次元ごとに独立な長さスケール ldl_d を持つことで、重要な次元（小さい ldl_d）と無関係な次元（大きい ldl_d）を自動的に識別する仕組みです。 たとえば目的関数にほとんど影響しないパラメータは ldl_d が大きく推定されるため、$r$ への寄与が小さくなり、その方向には曲面がほぼ平坦になります。

RBF ではなく Matérn 5/2 を選ぶ理由

カーネル 滑らかさ 特徴
RBF（Gaussian） C∞（無限微分可能） データから遠い領域で不確実性を過小評価しやすい
Matérn 5/2 C²（2 回微分可能） 工学・ハイパーパラメータ最適化の現実的な滑らかさに合致

工学の目的関数は完全に滑らかとは限らないので、無限回微分可能な RBF を仮定すると、実際より確信を持った予測になりがちです。 滑らかさを C² にとどめることで、データから離れた領域での過信を避けています。

スパース GP: 誘導点

N 個の訓練点をすべてカーネル求解に使う代わりに、M 個（M ≪ N）の代表点 Z={z1,…,zM}Z = {z_1, \ldots, z_M} を媒介変数として導入します:

u=f(Z)∼GP(0,KZZ)u = f(Z) \sim \mathcal{GP}(0, K_{ZZ})

訓練点をそのまま使うと N×N のカーネル行列を逆行列にかけることになり、N が大きいと計算が現実的でなくなるためです。

FITC は誘導変数 uu を条件としたとき各訓練点が独立であると仮定します:

p(f(X)∣u)≈∏ip(f(xi)∣u)p(f(X)\mid u) \approx \prod_i p(f(x_i)\mid u)

主要な行列:

行列 サイズ 内容
K_ZZ M × M 誘導点間のカーネル行列
K_XZ N × M 訓練点と誘導点間のカーネル行列

Q 行列（低ランク近似）:

QXX≈KXZKZZ−1KXZ⊤Q_{XX} \approx K_{XZ}K_{ZZ}^{-1}K_{XZ}^\top

FITC ダイアゴナル行列 Λ:

Λ=diag⁡(σf2−Qdiag)+σn2I\Lambda = \operatorname{diag}(\sigma_f^2 - Q_{\mathrm{diag}}) + \sigma_n^2 I

Woodbury 恒等式により、高コストな N×N 逆行列を M×M 演算に削減できます:

(Q+Λ)−1=Λ−1−Λ−1KXZΣ−1KXZ⊤Λ−1(Q+\Lambda)^{-1} = \Lambda^{-1} - \Lambda^{-1}K_{XZ}\Sigma^{-1}K_{XZ}^\top\Lambda^{-1}

Σ=KZZ+KXZ⊤Λ−1KXZ\Sigma = K_{ZZ} + K_{XZ}^\top\Lambda^{-1}K_{XZ}

主要コスト: O(N · M²)。 M を 100 で頭打ちにしているので、訓練点が増えてもコストは N に対して線形にしか伸びません。

FITC と VFE の違い

FITC と VFE はどちらも同じ M×M 計算形式に帰着しますが、ハイパーパラメータ学習時に最大化する目的関数が異なります:

基準 目的関数 ノイズ推定の傾向
FITC FITC 周辺尤度近似 データに密着した積極的なノイズ推定
VFE Variational Free Energy（ELBO） 真の GP 周辺尤度の真の下界；やや保守的なノイズ推定 → 滑らかなフィット

VFE の目的関数:

LVFE=log⁡N(y; 0, QXX+σn2I)−12σn2tr⁡(KXX−QXX)\mathcal{L}\text{VFE} = \log \mathcal{N}(y;,0,,Q{XX}+\sigma_n^2 I) - \frac{1}{2\sigma_n^2}\operatorname{tr}(K_{XX}-Q_{XX})

トレース項は誘導点近似によって失われた情報量へのペナルティであり、FITC にはこの項がありません。 このペナルティを課すぶん VFE のフィットは控えめになるため、次のベンチマークでも R² は GP-FITC より低く出ています。

実用的なベンチマーク（ノイズなし関数、N = 100、M = 100；開発環境での計測による参考値であり、環境やデータによって変動します）: GP-FITC R² ≈ 0.88、GP-VFE R² ≈ 0.76 です。 VFE はわずかにフィットを犠牲にして、より原理的な下界を実現します。 ノイズのあるデータでは両者はほぼ一致します。 そのため R² の数字だけを見て GP-VFE を切り捨てるのは早計で、曲面がスパイク状になっていないかとあわせて判断してください。

使い分けの目安:

  • GP-FITC：デフォルト。クリーンまたは軽度ノイズデータで最良のフィット。まずこちらを選びます。
  • GP-VFE：GP-FITC の曲面が過学習気味、スパイク状に見える場合。VFE の保守的なノイズ推定がより滑らかな曲面を生み出します。

観測ノイズ分散

等分散ノイズ分散 σn2\sigma_n^2 はカーネルハイパーパラメータと同時に推定されます。 下限は正規化 y 単位で 1e-6 です（数値失敗時は 1e-3 で再試行します）。 ノイズ分散に下限を設けることで共分散行列の正定値性が保たれ、Cholesky 分解が数値的に失敗しにくくなります。

ハイパーパラメータ最適化

egobox は選択した目的関数（FITC 尤度または VFE ELBO）を、勾配不要の COBYLA オプティマイザと 10 点マルチスタート（決定論的: 固定グリッド、固定シード）を用いて最大化します。

最適化変数は θ=[log⁡l1,…,log⁡lD,log⁡σf,log⁡σn]\theta = [\log l_1,\ldots,\log l_D,\log\sigma_f,\log\sigma_n] です。 対数空間で最適化することで、すべてのパラメータは明示的な制約なしに正値を保ちます。 マルチスタートの初期点も固定グリッドと固定シードで決めているので、同じデータからはいつも同じモデルが得られます。

データ正規化

GP はハイパーパラメータの初期値 log⁡ls=0\log l_s = 0（長さスケール =1=1）を前提としているため、$x$ と yy が [0,1][0,1] 程度のスケールでないと最適化が適切に収束しません。

  • X: 各次元を [0,1] にスケーリング（min/max 正規化）
  • Y: Z スコア正規化（平均 0、標準偏差 1）
  • グリッド予測を逆変換し、元のスケールに戻す

グリッド予測は逆変換して元のスケールに戻すので、正規化が行われていることは表示側からは見えません。

計算量

処理 計算量 N = 10,000, M = 100 の概算
FITC/VFE カーネル行列 O(N·M²) 1×10⁸ ops
Cholesky 分解（M×M） O(M³) 1×10⁶ ops
グリッド予測（50×50） O(2500·M) 2.5×10⁵ ops

目標: 10,000 ms 以内（release ビルド）。 実測（N = 10,000）: GP-FITC ≈ 2.4 秒、GP-VFE ≈ 2.0 秒（開発環境での計測による参考値であり、環境やデータによって変動します）。

R² の見方

R² 対応
≥ 0.8 良好なフィットです。曲面は信頼できます。
< 0.5 フィット不良。Random Forest または LightGBM に切り替えるか、N を増やします。

モデルを切り替えても R² が上がらない場合は、試行数そのものが足りていない可能性があります。

特性・限界

強み:

  • データが少なくても高品質な補間（N = 20 程度でも機能します）
  • 滑らかな曲面（Matérn 5/2 の C² 連続性）
  • ARD により各次元の重要度を自動推定
  • 全 N 点で学習し、データのサブサンプリングは行いません（egobox-gp バックエンド）
  • 95% 信頼帯による不確実性の可視化

弱み:

  • N > 100 では M = 100 の誘導点上限によりコストを抑えますが近似誤差が生じます
  • GP-FITC はノイズが多いデータで過学習しやすい（その場合は GP-VFE を使用）
  • COBYLA マルチスタートが局所最適解に収束する場合があります

ですが、GP の価値は点予測の精度そのものよりも、信頼帯を通じて「どこがまだ分かっていないか」を示せる点にあります。 信頼帯が広い領域は試行が足りていない領域なので、次にどこを探索すべきかの手がかりとして読めます。

使用場面の目安

目的関数の形が...

  線形に近い ───────────────────────────────────→ Ridge（最速）
  非線形・不連続 ──────────────────────────────→ LightGBM または Random Forest
  滑らかな非線形（デフォルト） ─────────────────→ GP-FITC（最高品質・デフォルト）
  滑らかな非線形・曲面が過学習気味 ────────────→ GP-VFE（より滑らか・保守的）
  不連続・多領域の応答 ────────────────────────→ GP-MOE（gaussian-process-moe.md 参照）
  • N < 50 の Study では GP-FITC がもっとも信頼性の高い補間を提供します
  • R2>0.8R^2 > 0.8 なら曲面はデータをよく説明できています

参考文献

← 前へサロゲートモデル: Random Forest次へ →サロゲートモデル: GP-MOE