Solving Fredholm Integral Equations of the Second Kind via Wasserstein Gradient Flows

1. 問題設定 (Problem)

フレドホルム積分方程式の第 2 種は、以下の形式で定義されます。
$$\pi(x) = \phi(x) + \lambda \int_{\mathbb{R}^d} k(x, y)\pi(y)dy$$
ここで、$\pi$ は未知の関数（解）、$\phi$ は既知のforcing項、$k$ は積分核、$\lambda$ は定数です。

• 課題: 従来の数値解法（基底関数による離散化や数値積分）は低次元では有効ですが、次元 $d$ が高い場合や、解の定義域が有界でない（unbounded）場合、計算コストが爆発的に増大するか、精度が低下するという問題があります。
• 目的: 解 $\pi$ が確率分布（確率密度関数）である場合を想定し、非有界領域における第 2 種のフレドホルム方程式を効率的に解く手法の開発。特に、核 $k$ による積分演算子が非可逆または条件数が悪い（ill-posed）場合の安定した解法が求められています。

2. 手法 (Methodology)

著者らは、解 $\pi$ を求める問題を、特定の汎関数（Functional）の最小化問題として定式化し、その**ワッセルシュタイン勾配流（Wasserstein Gradient Flow）**を粒子系で近似するアプローチを採用しました。

2.1 最小化問題の定式化

解 $\pi$ を求めるために、以下の正則化付き汎関数 $F_\alpha$ を最小化する問題を考えます。
$$F_\alpha(\pi) = \text{KL}\left(\pi \bigg\| \phi + \lambda \int_{\mathbb{R}^d} k(\cdot, y)\pi(y)dy \right) + \alpha \text{KL}(\pi | \pi_0)$$

• KL 発散: 第一項は、解 $\pi$ と方程式の右辺で定義される分布との間のカルバック・ライブラー（KL）発散です。
• 正則化項: 第二項は、参照測度 $\pi_0$ に対する KL 発散で、正則化パラメータ $\alpha > 0$ によって制御されます。これにより、解の一意性と安定性が保証されます。
• 正則化の必要性: $\alpha > 0$ かつ $\pi_0$ が適切に選ばれることで、この汎関数は狭義凸（strictly convex）かつ強制的（coercive）となり、一意な最小解が存在することが証明されます。

2.2 ワッセルシュタイン勾配流と McKean-Vlasov SDE

汎関数の勾配流を導出すると、確率測度の時間発展を表す偏微分方程式（Fokker-Planck 方程式）が得られます。これに対応する確率微分方程式（SDE）は、McKean-Vlasov 型 SDE となります。
$$dX_t = \left( \int_{\mathbb{R}^d} b_\eta(X_t, z, \pi_t) d\pi_t(z) - \alpha \nabla U(X_t) \right) dt + \sqrt{2(1+\alpha)} dB_t$$

• 特徴: この SDE のドリフト項は、粒子の分布 $\pi_t$ 自体に依存しており、さらにその分布に対する積分（期待値）がネスト（入れ子）構造になっています。
  • 第 1 種の方程式に対する既存手法 [16] と異なり、核 $k$ を $\pi$ に対して 2 回積分する必要があり、ドリフトの計算が複雑化しています。
  • $\eta > 0$ は分母の安定化のために導入されたハイパーパラメータです。

2.3 相互作用粒子系と数値アルゴリズム

解析的な SDE の解は直接計算できないため、**相互作用粒子系（Interacting Particle System）**を用いて近似します。

• 粒子系: $N$ 個の粒子 $\{X^i_t\}_{i=1}^N$ を用意し、各粒子の運動方程式を離散化（Euler-Maruyama 法）してシミュレーションします。
• 計算コストの工夫: ドリフト項のネスト構造により、単純な実装では計算コストが $O(N^3)$ になりますが、中間結果（$\pi^N_n[k(X^j, \cdot)]$）を事前計算・保存することで、$O(N^2)$ まで削減可能です。さらに並列化により $O(N)$ まで可能です。
• アルゴリズム:
  1. 初期分布から粒子をサンプリング。
  2. 各ステップで粒子間の相互作用（核関数による重み付け）を計算し、ドリフトを評価。
  3. 正則化項（$\nabla U$）と拡散項を加えて粒子を更新。
  4. 収束後に、粒子の分布から滑らかな解 $\hat{\pi}$ を再構成（核密度推定など）。

3. 主要な貢献 (Key Contributions)

1. 第 2 種方程式への拡張: 第 1 種の方程式に対する Wasserstein 勾配流アプローチを、第 2 種（解が両辺に現れる）に拡張し、ネストされた期待値を含む非線形 SDE を導出しました。
2. 理論的保証:
   • 正則化された汎関数の一意な最小解の存在と、正則化パラメータ $\alpha, \eta \to 0$ における収束性を証明しました。
   • 粒子系の**カオス伝播（Propagation of Chaos）**と、離散化されたアルゴリズムの収束性（強収束次数 1/2）を証明し、誤差の上限を導出しました。
   • 粒子数 $N$ と時間ステップ $\gamma$ の最適バランス（$N \propto B^{1/3}$）を示唆し、計算コスト $O(N^3)$ に対して誤差が $O(N^{-1/2})$ で減衰することを示しました。
3. 数値的検証: 複数の例題（ガウス分布の単純な例、Karhunen-Loève 展開、Gaussian Process State Space Models）を通じて、既存手法（可逆ジャンプ MCMC、Nystrom 法）と比較し、特に高次元や不安定な問題において優れた性能を示しました。

4. 結果 (Results)

• 正則化の効果: 参照測度 $\pi_0$ を適切に設定することで、解の安定性が大幅に向上し、特に $\lambda$ が 1 に近い（積分演算子が特異に近い）場合でも安定した解が得られました。
• 既存手法との比較:
  • Reversible Jump MCMC: $\lambda \to 1$ で von Neumann 級数が不安定になるため失敗しますが、提案手法は $\lambda=1$（固有値問題）でも機能しました。
  • Nystrom 法: 離散化格子に依存するため高次元では困難ですが、提案手法は適応的な粒子配置により、同等の計算コストでより高精度な近似（最大 5 倍の精度向上）を実現しました。
• 収束性: 粒子数 $N$ の増加に伴い、解の誤差が理論的に予測された $O(N^{-1/2})$ のレートで減少することが確認されました。

5. 意義 (Significance)

• 高次元・非有界領域への適用: 従来の格子ベースの手法が直面する「次元の呪い」や境界条件の扱いの難しさを克服し、確率粒子法を用いることで柔軟な空間探索を可能にしました。
• 逆問題としての安定性: 積分方程式は本質的に ill-posed な逆問題ですが、KL 正則化と Wasserstein 勾配流の組み合わせにより、安定した数値解を提供する枠組みを確立しました。
• 応用範囲の広がり: 経済計量学、光輸送モデリング、強化学習、ガウス過程状態空間モデルなど、確率測度を解とする多様な分野での応用可能性が開かれました。
• 理論的深さ: ネストされた期待値を持つ McKean-Vlasov SDE に対する収束解析は、確率数値解析の分野において新たな知見を提供しています。

総じて、この論文は、積分方程式の数値解法において、確率的粒子法と最適輸送理論（Wasserstein 幾何）を融合させることで、従来の決定論的離散化手法の限界を突破する強力なアプローチを提示したものです。