Tunny Icon
TunnyDocs

The next-gen Grasshopper optimization tool.

Sobol

概要

Sobol 感度指数は、目的関数の分散分解に基づいてパラメータの重要度を定量化するグローバル感度分析手法です。 線形か非線形かを問わず、パラメータ間の交互作用も含めた影響度を [0,1][0,1] の値として返します。

Tunny Dashboard では 2 種類の指数を提供します:

指数 記号 意味
一次指数 (First-Order) SiS_i パラメータ xix_i の単独寄与率
全効果指数 (Total-Effect) STiST_i xix_i の単独 + 他パラメータとの交互作用を含む総寄与率

理論背景

ANOVA 分散分解

Sobol 分解は、モデル出力 Y=f(X1,,Xp)Y = f(X_1, \ldots, X_p) の分散を各パラメータの寄与に分解します:

Var(Y)=iVi+i<jVij++V1p\mathrm{Var}(Y) = \sum_i V_i + \sum_{i<j} V_{ij} + \cdots + V_{1\ldots p}

  • Vi=Var(E[YXi])V_i = \mathrm{Var}(\mathbb{E}[Y \mid X_i])xix_i のみが変化したときの YY の期待値の分散
  • Vij=Var(E[YXi,Xj])ViVjV_{ij} = \mathrm{Var}(\mathbb{E}[Y \mid X_i, X_j]) - V_i - V_jxix_ixjx_j の交互作用分散

一次指数

Si=ViVar(Y)S_i = \frac{V_i}{\mathrm{Var}(Y)}

x_i の変化だけで説明できる Y の分散の割合です。 S_i の合計は ≤ 1(交互作用がなければ = 1)です。

全効果指数

STi=E[Var(YXi)]Var(Y)=1Var(E[YXi])Var(Y)ST_i = \frac{\mathbb{E}[\mathrm{Var}(Y \mid X_{\sim i})]}{\mathrm{Var}(Y)} = 1 - \frac{\mathrm{Var}(\mathbb{E}[Y \mid X_{\sim i}])}{\mathrm{Var}(Y)}

X_{~i} は「x_i 以外のすべてのパラメータ」です。 ST_i は、x_i の単独効果と x_i が絡むすべての交互作用効果の合計です。 STiSiST_i \ge S_i は常に成立します。 差 STiSiST_i - S_i が大きいほど、 xix_i が他パラメータとの交互作用に強く関与していることを意味します。

有限サンプル推定の注意: Saltelli 一次指数推定量と Jansen 全効果推定量は独立に計算されるため、有限サンプルでは S^i>ST^i\hat{S}_i > \hat{ST}_i が生じることがあります。 Tunny Dashboard では推定後に ST^imax(ST^i,S^i)\hat{ST}_i \leftarrow \max(\hat{ST}_i, \hat{S}_i) を適用してから [0,1][0,1] にクランプし、理論的性質を出力値で保証しています。


Tunny Dashboard での計算方法

全体フロー

直接 Monte Carlo 積分を行うと n × p 倍のシミュレーション実行が必要になり、最適化問題では現実的ではありません。 代わりに二次 Ridge サロゲートモデルを構築し、そのモデル上で Saltelli サンプリングを行います。

  1. 二次 Ridge サロゲートを構築します。
  2. Saltelli 行列 A, B を ChaCha8 で生成します。
  3. A, B, AB_i でサロゲートを評価します。
  4. Jansen 推定量で S_i, ST_i を計算します。

Step 1: 二次 Ridge サロゲートモデル

特徴量構成

p 個の標準化済み線形入力から p(p+3)/2 次元の二次特徴量ベクトルを構築します:

ϕ(x)=[x1,,xp,  x12,,xp2,  x1x2,,xp1xp]\phi(x) = [x_1, \ldots, x_p,; x_1^2, \ldots, x_p^2,; x_1x_2, \ldots, x_{p-1}x_p]

例: p = 3 → 3 + 3 + 3 = 9 次元 (= 3×6/2)

したがって、サロゲートの特徴量次元は

dϕ=p+p+p(p1)2=p(p+3)2d_\phi = p + p + \frac{p(p-1)}{2} = \frac{p(p+3)}{2}

であり、学習サンプル数 nn は少なくとも dϕd_\phi と同程度、実務上は ndϕn \gg d_\phi が望ましいです。

サロゲートの学習

  1. XX の各列を Z スコア標準化

xj=xjμjσj\tilde x_j = \frac{x_j-\mu_j}{\sigma_j}j=σjxjμj

  1. 二次特徴量 Φ\Phi を構築
  2. Φ\Phi の各列を再度 Z スコア標準化
  3. 各目的関数に対して Ridge( α=1.0\alpha=1.0) を適合

βk=(ΦTΦ+I)1ΦT(ykyˉk)\beta_k = (\Phi^T\Phi + I)^{-1}\Phi^T(y_k-\bar y_k)

サロゲートの評価

f^(x)=βTϕ(x)+intercept\hat f(x) = \beta^T \tilde\phi(x) + \mathrm{intercept}(x)+intercept

ここで φ̃(x) は、二次特徴量を訓練データの統計量で標準化したものです。


Step 2: Saltelli サンプリング

疑似乱数生成: ChaCha8

ChaCha8 アルゴリズムに基づき、シードベースの決定論的な一様乱数([0, 1) の範囲)を生成します。

初期シード: 0xDEAD_BEEF_1234_5678(固定・再現性あり)

行列 A, B の生成

各パラメータの値域 [lo_j, hi_j] は、実際のトライアルの min/max から取得します。

A[s,j]=loj+u(hijloj),B[s,j]=loj+u(hijloj)A[s,j] = lo_j + u,(hi_j - lo_j),\qquad B[s,j] = lo_j + u,(hi_j - lo_j)

AB_i 行列

A の第 i 列のみを B の第 i 列で置換した行列です:

ABi[s]=[A[s,0],,A[s,i1],B[s,i],A[s,i+1],,A[s,p1]]AB_i[s] = [A[s,0], \ldots, A[s,i-1], B[s,i], A[s,i+1], \ldots, A[s,p-1]]


Step 3: サロゲート評価

fA[k,s]=surrogate_eval(A[s],objective=k)f_A[k,s] = \mathrm{surrogate_eval}(A[s],,\mathrm{objective}=k)

fB[k,s]=surrogate_eval(B[s],objective=k)f_B[k,s] = \mathrm{surrogate_eval}(B[s],,\mathrm{objective}=k)

fABi[k,s]=surrogateeval(ABi[s],objective=k)f{AB_i}[k,s] = \mathrm{surrogate_eval}(AB_i[s],,\mathrm{objective}=k)


Step 4: Jansen 推定量

一次指数(Saltelli 2010)

SisfB[s](fABi[s]fA[s])NVarYS_i \approx \frac{\sum_s f_B[s],(f_{AB_i}[s]-f_A[s])}{N,\mathrm{Var}_Y}

全効果指数(Jansen 1999)

STis(fA[s]fABi[s])22NVarYST_i \approx \frac{\sum_s (f_A[s]-f_{AB_i}[s])^2}{2N,\mathrm{Var}_Y}

分散の推定

VarY=mean(fA2)mean(fA)2\mathrm{Var}_Y = \mathrm{mean}(f_A^2) - \mathrm{mean}(f_A)^2 VarY0\mathrm{Var}_Y \approx 0 の場合(サロゲートが定数に近い)は Si=STi=0S_i = ST_i = 0 を返します。

クリッピング

推定誤差による範囲外値を防ぐため、最終値を [0,1][0,1] にクリップします:

Si=clamp(Siraw,0,1),STi=clamp(STiraw,0,1)S_i = \mathrm{clamp}(S_i^{\mathrm{raw}}, 0, 1),\qquad ST_i = \mathrm{clamp}(ST_i^{\mathrm{raw}}, 0, 1)


複数目的関数への対応

Importance チャートでは、各パラメータのスコアとして全目的関数に対する指数の平均を表示します:

display_score(pj)=1mkSi[j][k]\mathrm{display_score}(p_j)=\frac{1}{m}\sum_k S_i[j][k]

display_score(pj)=1mkSTi[j][k]\mathrm{display_score}(p_j)=\frac{1}{m}\sum_k ST_i[j][k]


パラメータ設定

設定項目
Saltelli サンプル数 1024
Ridge 正則化強度 α 1.0
乱数シード 0xDEAD_BEEF_1234_5678

必要データ量

  • 最低 2 トライアル(n ≥ 2)、p ≥ 1、目的関数 ≥ 1 が必要です
  • サロゲートの精度は n とともに向上します。目安として n ≥ 10 × p(p+3)/2 を推奨します

パラメータの扱いに関する前提

数値列はそのまま使用し、カテゴリ列は文字列ラベルを出現順の整数 ID(0.0, 1.0, …)へラベル符号化して使用します。

カテゴリパラメータはラベル符号化済みの数値として扱われるため Sobol 計算の対象となりますが、整数 ID は順序や距離の情報を持たないため、Sobol 指数の解釈はカテゴリ変数に対して近似的なものとなります。


特性・限界

強み:

  • 線形、非線形、交互作用をすべて扱えます
  • 値域が [0,1][0,1] で複数パラメータ間の比較が容易です
  • ST_i - S_i で交互作用の強度を把握できます
  • サロゲートを使うため、トライアル数が限られていても計算可能です

弱み:

  • サロゲートモデル(二次 Ridge)の精度に結果が依存します
    • 真の関数が高次非線形または強い不連続性を持つ場合は精度が落ちます
  • p が大きくなると特徴量数 p(p+3)/2 が増加し、サロゲートの学習精度が低下する可能性があります
  • 乱数シードが固定のため再現性はありますが、N の選択で結果が変わります

参考文献

  • Saltelli, A. et al. (2010). Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index. Computer Physics Communications, 181(2), 259–270. https://doi.org/10.1016/j.cpc.2009.09.018
  • Jansen, M. J. W. (1999). Analysis of variance designs for model output. Computer Physics Communications, 117(1), 35–43. https://doi.org/10.1016/S0010-4655(98)00154-4
  • Sobol, I. M. (1993). Sensitivity estimates for nonlinear mathematical models. Mathematical Modelling and Computational Experiments, 1(4), 407–414.
© 2026 hrntsm
Made with Fresh