› Kriging サロゲートモデル シミュレーター 戻る
サロゲート / 機械学習

Kriging サロゲートモデル シミュレーター — ガウス過程回帰

CFD や FEM の高コスト評価を少数の訓練点で代替する Kriging（ガウス過程回帰）。SE / Matérn カーネルと相関長 ℓ・ナゲット σ_n² を変えながら、予測平均と95%信用区間・対数周辺尤度・RMSE の振る舞いを直感的に学べます。

パラメータ設定
カーネル関数
真関数の滑らかさの事前仮定
トレーニング数 N
高コスト評価点（CFD/FEM/実験）の数
相関長 ℓ
カーネル幅。大きいほど滑らかな関数を仮定
信号分散 σ_f²
関数の振幅スケール
ナゲット σ_n²
観測ノイズの分散。決定論的 CAE では 1e-6〜1e-3
真関数振幅 A
f(x) = A·sin(2πx) の振幅
計算結果
—
トレーニング数 N
—
平均サンプル間隔
—
予測標準偏差 σ_pred
—
95% 信用区間幅
—
RMSE 予測誤差
—
対数周辺尤度
—
Kriging 補間 — 真関数・訓練点・予測±95%CI

白：真関数 A·sin(2πx)、赤丸：訓練点（観測ノイズ付）、青線：Kriging 予測平均、青帯：95% 信用区間。観測ノイズが動的に揺らぐ様子も表示します。

Kriging 予測 + 信用区間（Chart.js）
対数周辺尤度 logML vs 相関長 ℓ
理論・主要公式

$$\hat\mu(x_*) = k_*^T K^{-1} y,\qquad \hat\sigma^2(x_*) = k(x_*,x_*) - k_*^T K^{-1} k_*$$

K = カーネル行列（N×N、K_{ij}=k(x_i,x_j)+σ_n²δ_{ij}）、k_* = 予測点と訓練点の共分散ベクトル。事後平均 μ̂ と分散 σ̂² を同時に提供する。

$$k_{SE}(r) = \sigma_f^2 \exp\!\left(-\frac{r^2}{2\ell^2}\right)$$

平方指数（RBF）カーネル。r=|x−x'|、ℓ=相関長、σ_f²=信号分散。無限回微分可能な「最も滑らかな」関数族を仮定する。

$$k_{M5/2}(r) = \sigma_f^2\left(1+\frac{\sqrt{5}\,r}{\ell}+\frac{5r^2}{3\ell^2}\right)\exp\!\left(-\frac{\sqrt{5}\,r}{\ell}\right)$$

Matérn 5/2 カーネル（Bayesian Optimization のデフォルト）。2回微分可能で、SE より急な変化も表現できる。

$$\log p(y\mid X,\theta) = -\tfrac12 y^T K^{-1} y - \tfrac12 \log|K| - \tfrac{N}{2}\log(2\pi)$$

対数周辺尤度。第1項=データ適合、第2項=モデル複雑さペナルティ。これを θ=(ℓ,σ_f²,σ_n²) で最大化してハイパーパラメータを推定する。

Kriging サロゲートモデル（ガウス過程回帰）

🙋
「サロゲートモデル」って、本物の CFD の代わりに使うってことですよね？でも、近似なら多項式回帰でも良くないですか？なんで Kriging なんですか？
🎓
いい質問だね。Kriging が他の回帰と決定的に違うのは、「予測値だけじゃなく、予測の不確かさ（分散）も同時に出してくれる」点なんだ。多項式や RBF 補間は「中央値らしき値」だけを返すけど、Kriging は各点 x_* で「平均 μ̂(x_*)」と「標準偏差 σ̂(x_*)」を一緒に返す。だから「ここはサンプル点が近いから自信ある」「ここはスカスカだから当てにならない」が定量的に分かる。これが Bayesian Optimization で次の評価点を選ぶ Expected Improvement の基礎になるんだ。
🙋
なるほど。じゃあ左の「相関長 ℓ」と「ナゲット σ_n²」って、どう決めればいいんですか？適当に動かしてみると RMSE が大きく変わって、何が正解か分からなくなります…
🎓
そこが Kriging の腕の見せどころで、答えは「対数周辺尤度（logML）を最大化する ℓ を選ぶ」だよ。右下のチャートが logML を ℓ の関数として描いてくれているから、ピークの位置を読めばいい。実務では BFGS や L-BFGS で数値最適化するけど、よく局所最適にハマるから、ℓ=0.05, 0.5, 5 みたいに複数の初期値から開始する「マルチスタート」が標準。GPyOpt、scikit-learn の GaussianProcessRegressor、SMT、Dakota などのライブラリは全部この最適化を内部でやってくれる。
🙋
カーネルを SE と Matérn で切り替えると予測曲線の感じが結構変わりますね。SE が一番きれいに見えるけど、Matérn を使う場面ってあるんですか？
🎓
SE は「無限回微分可能」を仮定するから、本当に滑らかな現象には強い。たとえば翼型の揚力係数を AoA に対して回帰するなら SE で十分。でも、座屈・接触・相転移みたいに「ある点で挙動が急変する」物理を扱うと、SE は無理に滑らかに通そうとして逆に振動するんだ。そういうケースは Matérn 3/2（1回微分可能）や Matérn 5/2（2回微分可能）の方が誠実。Bayesian Optimization の世界では Matérn 5/2 が事実上のデフォルトになっていて、迷ったらまず 5/2 から始めるのが定石だよ。
🙋
N=8 のデフォルトだと RMSE が 1.07 もあって判定が NG ですが、N を増やせば改善しますよね？どのくらい必要なんですか？
🎓
そう、N を 20 くらいまで増やすと sin(2πx) なら RMSE が 0.2 を切る感じになる。経験則として「1次元なら 10〜20点、2次元なら 30〜50点、5次元なら 100〜200点」あれば実用的な精度が出る。次元 d に対して必要点数は概ね 10^d でスケールするので、これが「次元の呪い」と呼ばれる現象。実務で次元が 10 を超えると、素の Kriging から PCE や Deep Kernel Learning に乗り換える判断ポイントになるよ。
🙋
最後にもうひとつ。CFD 結果みたいに「ノイズがない決定論的な出力」を補間するときも、ナゲット σ_n² を入れた方がいいんでしょうか？
🎓
実は数値的にも入れた方がいい。決定論的なら理論上は σ_n²=0 で良いんだけど、それだと K がほぼ特異行列になって K⁻¹ が爆発する。だから「jitter」と呼ばれる小さな値（1e-6〜1e-3）を必ず加える。これは「ナゲット」と同じ働きで、数値的安定性のための正則化項。一方、実験データなら測定誤差の分散をそのまま入れる。CAE 結果と実験データを混ぜる Multi-Fidelity Kriging では、各データセットに別々のナゲットを与える設計が標準だよ。

よくある質問

本質的には同じ手法です。Kriging は1950年代に南アフリカの鉱山技師 Krige が地下鉱床の品位を補間するために提案した地質統計学の用語で、機械学習コミュニティでは同じ式が Gaussian Process Regression（GPR）と呼ばれます。両者とも「観測点の重み付き和で予測値を作り、予測分散も同時に出す」線形ベイズ推定で、カーネル関数 k(x,x') によって滑らかさを規定します。実装上の違いは、Kriging はバリオグラム γ(h) で書くのに対し、GPR は共分散関数 k(r) で書く流儀ぐらいです。CAE・最適化文脈では Kriging、機械学習・UQ 文脈では GPR と呼ぶことが多いと考えてください。
SE（平方指数・RBF）カーネルは無限回微分可能な「最も滑らかな」関数を仮定します。空力解析や熱問題のように真の応答が滑らかな場合に第一選択です。Matérn 3/2 は1回しか微分できない関数（角がある／不連続な微分を含む）を表現でき、構造物の座屈や接触問題のように応答が急変する場合に向きます。Matérn 5/2 はその中間で、Bayesian Optimization の事実上のデフォルトです。迷ったら Matérn 5/2 で開始し、残差プロットや交差検証で SE と比較するのがおすすめです。
通常は対数周辺尤度（log marginal likelihood）を最大化してハイパーパラメータを推定します。本ツールの2番目のグラフがこの logML を ℓ の関数として表示しており、ピークがある ℓ が最尤推定値です。実務では BFGS や L-BFGS で数値最適化しますが、観測点が少ない場合は局所最適に陥りやすいため、複数の初期値（ℓ=0.05, 0.5, 5 など）から開始するマルチスタートが標準です。ナゲット σ_n² は観測ノイズの分散で、CAE 結果のように決定論的な出力には小さな値（1e-6〜1e-3）を、実験データには測定誤差の分散をそのまま入れます。
標準的な GPR は計算量 O(N³)・メモリ O(N²) なので、訓練点数 N が 2,000 程度を超えると一般の PC では厳しくなります。次元数（入力変数の数）は ARD（Automatic Relevance Determination）で各次元別の相関長を学習することで20次元程度まで扱えますが、それ以上では「次元の呪い」で性能が落ちます。高次元・大規模データでは Sparse GP（Inducing Points）・SVGP・Deep Kernel Learning・PCE（Polynomial Chaos）・ニューラルネットなどへ移行します。CAE 設計最適化の典型ケース（5〜10次元・100〜500点）であれば、素の Kriging が最良の選択です。

実世界での応用

Bayesian Optimization（ベイズ最適化）：1回の評価に数時間かかる CFD や FEM の設計最適化で、Kriging サロゲートと Expected Improvement（EI）獲得関数を組み合わせ、「次にどこを評価すれば最も改善が期待できるか」を提案します。航空機翼形状最適化、ロケットノズル設計、化学プラント運転条件最適化など、数十〜数百回の高コスト評価で大域最適に到達する代表的な手法。Google Vizier、Facebook Ax、scikit-optimize などのオープンソース実装が広く使われます。

CAE 結果の代理モデルとロバスト設計：自動車衝突解析・電磁界解析・流体熱連成解析など、1ケース数十分〜数日かかるシミュレーションを Kriging で代替し、寸法・材料定数・荷重の不確かさを Monte Carlo で1万ケース流す UQ（不確かさ定量化）に使います。直接 CAE で1万回回すと現実的でないが、サロゲートで0.1秒/ケースになれば一晩で完了。Dakota、SMT、UQLab などのツールが標準。

地質統計学・資源探査：本来の Kriging が生まれた領域。坑井ボーリングで得られた数十点の鉱石品位や地下水位データから、3次元的に補間して鉱床のグレード分布図を作ります。石油・天然ガス探査、レアアース埋蔵量推定、地下水汚染拡散予測など、地球科学の標準ツール。バリオグラム解析と組み合わせ、空間相関構造の同定もセットで行うのが業界慣行。

機械学習の Uncertainty Quantification：自動運転の障害物検知・医療画像診断など「予測が間違っているかもしれない領域で警告を出したい」場面で、GPR を補助モデルとして用います。Deep Kernel Learning（DKL）はニューラルネットの最終層を GP に置き換える手法で、画像のような高次元入力に対しても予測分散を提供できます。NASA の宇宙船航法、Toyota の運転支援、Genentech の創薬スクリーニングなど、安全性が問われる領域で導入が進んでいます。

よくある誤解と注意点

最大の落とし穴が、「Kriging の信用区間を絶対視する」こと。GP の事後分散 σ̂² は「仮定したカーネルと観測データが正しい場合の」分散であり、カーネルが真の関数族と一致していなければ、信用区間は系統的にずれます。たとえば、本当は不連続のあるステップ関数を SE カーネル（無限回微分可能）で当てはめると、不連続点周辺で大きく外れているのに信用区間は「自信あり」を返してしまう。実務では交差検証で「実測値が95%信用区間に入る頻度」が公称値（95%）に近いか必ず確認してください（calibration test）。

次に、「ハイパーパラメータを固定したまま使い続ける」こと。設計探索の初期段階で推定した ℓ が、設計空間の他の領域でも適切とは限りません。新しい評価点が加わるたびにハイパーパラメータを再最適化するのが正しい運用です。逆にやりすぎて毎反復で大きく動くようなら、観測ノイズが過小か、初期サンプリングが不適切（LHS や Sobol 系列を使わず均等格子で取ってしまった、など）の兆候。Latin Hypercube Sampling と D-optimal design の組み合わせで初期点を取るのが定石です。

最後に、「Kriging は外挿に強い」と誤解すること。実は逆で、訓練データのバウンディングボックス外（外挿領域）では予測平均が事前平均（通常はゼロ）に急速に戻り、分散は σ_f² に発散します。つまり「外では何も知らない」と正直に答えるのが GP の特性で、これ自体は健全ですが、最適化アルゴリズムが境界の外を提案してしまうと無意味な評価を繰り返すことになります。Bayesian Optimization では設計領域の境界をハードに enforce し、必要なら境界に近づくほど EI を減衰させる工夫を入れます。また、PCE のように「グローバルな多項式基底」を選ぶ手法と組み合わせるハイブリッド（PC-Kriging）も有力な解です。

使い方ガイド

  1. トレーニングデータ点数（numSamples）を1～50の範囲で設定し、高コスト解析（有限要素法・流体解析など）の評価点数を指定
  2. カーネル関数（SquaredExponential/Matérn3-2/Matérn5-2）を選択し、長さスケール（lengthScale）を0.1～5.0の範囲で調整して関数の滑らかさを制御
  3. 信号分散（signalVariance）とナゲット分散（nuggetVariance）を入力してモデルの不確実性を定量化し、予測値・95%信用区間幅・RMSE・対数周辺尤度を確認

具体的な計算例

航空機翼の曲げ剛性最適化で、12個の有限要素解析評価点（スパン方向剛性材厚さ1.2～4.5mm）からKrigingモデルを構築。長さスケール=0.8、信号分散=850、ナゲット分散=2.5の設定で、厚さ2.8mmにおける予測曲げ剛性=245 N·m²、95%信用区間幅=18.5 N·m²、RMSE=3.2 N·m²、対数周辺尤度=-8.6を得た。これにより未評価点での応答予測と不確実性を同時に把握でき、次の評価点の最適配置（ラテン超方形計画など）を決定可能

実務での注意点

  1. ナゲット分散は計算誤差や製造ばらつきを反映；0.5～5.0の範囲で小値は過学習を招くため慎重に設定
  2. Matérn5-2カーネルは自動車ボディ応力解析など滑らかな応答曲面で高精度；Matérn3-2は複合材板のバケリング荷重など中程度の滑らかさに適用
  3. トレーニング数が少ない場合（N<8）は信用区間が広がり、設計空間全域では信頼度が低いため追加サンプリングが必須
  4. 対数周辺尤度が-15以下の場合はカーネルパラメータが過度に調整されており、正規化された入力データ（平均0、標準偏差1）の再確認が必要