論文「Simultaneous model discovery and state estimation under high data corruption」の技術的サマリー

この論文は、不完全でノイズの多いデータから常微分方程式（ODE）を発見し、同時に状態変数を推定するための新しいスパース回帰戦略を提案しています。従来の手法が抱えるノイズへの感度や不完全データへの対応不足を克服し、統計的な尤度関数に基づいたハイブリッド回帰アプローチを採用しています。

以下に、問題定義、手法、主要な貢献、結果、および意義について詳細にまとめます。

---

1. 問題定義と背景

データ駆動型のモデル発見（Model Discovery）は、物理法則と観測データを組み合わせる重要な課題ですが、以下の現実的な制約に直面しています。

• 高ノイズデータ: 従来の残差回帰（例：SINDy アルゴリズム）は、データ微分に依存するためノイズに極めて敏感です。
• 不完全なデータ: 観測データがドメイン全体に存在しない、または一部の状態変数が測定されていない（潜在変数）場合、従来の手法は機能しません。
• 状態推定の難しさ: カオス系などでは、発見されたモデルを時間積分しても初期条件のノイズにより真の軌道から急速に乖離します。
• パラメータ非線形性: 多くの既存手法はパラメータに対して線形な項のライブラリを前提としていますが、より一般的な非線形パラメータを含むモデルには対応できません。

本研究は、これらの課題に対処し、モデルのパラメータと状態変数を同時に推定する枠組みを構築することを目的としています。

2. 提案手法：ハイブリッドスパース回帰

2.1 統計的定式化（尤度関数）

提案手法は、データとモデルの両方の不確実性を考慮した**混合回帰（Hybrid Regression）**アプローチを採用しています。
尤度関数は、データ誤差とモデル誤差（残差）のバランスを取る形で定義されます。

$$L(\theta, u) = \min_{u} \left( \|N(u; \theta)\|^2 + \lambda \|\hat{u} - u|_D\|^2 \right)$$

ここで、

• $N(u; \theta) = 0$: 微分方程式（離散化された形式）。
• $\hat{u}$: 観測データ（ノイズあり、欠損あり）。
• $u$: 推定される真の状態（滑らかな解）。
• $\lambda$: モデルの不確実性とデータの不確実性の比率（分散比）を表す重み。
• $D$: データが観測されている部分集合。

この定式化により、$\lambda$ が小さい場合は古典的なパラメータ推定（最小二乗法）に、$\lambda$ が大きい場合は残差最小化（SINDy 的アプローチ）に近づきます。本研究では中間的な値を用いることで、ノイズに対して頑健な状態推定とパラメータ推定を同時に行います。

2.2 スパース性の確保とモデル選択

モデルのスパース性（不要な項を排除）を確保するために、以下の戦略を採用しています。

• スパース化ペナルティ: $l_0$ ノルム（非ゼロ成分の数）の滑らかな近似（Mohimani et al. による指数関数型ペナルティ）を使用し、勾配ベースの最適化を可能にします。
• モデル選択アルゴリズム: 硬い閾値ではなく、**ベイズ情報量基準（BIC）**を用いてモデルを評価・選択します。
  • 適応的剪定（Adaptive Pruning）: 初期モデルから不要な項を BIC が改善する限り順次削除します。
  • 前方・後方リファインメント: 削除しすぎた場合や、最適なモデルが見つからない場合に、項を戻す（後方）またはさらに削除する（前方）動的なプロセスを実行します。これにより、局所解に陥りにくく、最適なモデル構造を探索します。

2.3 最適化手法

高次元かつ非凸な最適化問題に対して、以下の技術を用いて効率的に解いています。

• 第 2 次最適化（Levenberg-Marquardt 法）: 勾配法や準ニュートン法ではなく、真のヘッセ行列（Hessian）を用いた第 2 次変分法を採用。
• 自動微分（Automatic Differentiation）: 勾配とヘッセ行列の計算を自動化。
• グラフ彩色アルゴリズム: ヘッセ行列のスパース性を利用し、必要なヘッセ行列 - 積（HVP）の回数を最小化することで、大規模な計算を可能にしています。

3. 主要な貢献

1. ノイズと欠損データへの頑健性: 直接的なフィルタリングや微分を行わず、状態変数 $u$ を推定することで、高ノイズ（50% 以上）や連続したデータ欠損（30% 以上）がある場合でも、正確なモデルと状態を復元可能であることを示しました。
2. パラメータ非線形性の対応: パラメータに対して線形なライブラリに限定されない一般化された定式化を提供。指数関数項など、パラメータ非線形な項を含むモデル（Colpitts 発振器など）の発見に成功しました。
3. 効率的な第 2 次最適化の実装: 自動微分とグラフ彩色アルゴリズムを組み合わせることで、システム同定において第 2 次最適化手法を実用的な計算コストで適用可能にしました。
4. BIC によるモデル選択: 従来の SINDy などの硬い閾値削除ではなく、統計的基準（BIC）に基づく適応的剪定アルゴリズムを導入し、過学習を防ぎつつ最適なモデルを特定します。

4. 実験結果

提案手法は、Van der Pol 振動子、Lorenz 系、Lorenz-96 系、Colpitts 発振器など、多様なシステムで検証されました。比較対象として、現在の最先端であるWSINDy（弱形式の SINDy）が用いられました。

• ノイズ耐性: 50% のノイズレベルでも、WSINDy よりもはるかに高いパラメータ推定精度とモデル構造の復元率（True Positive Rate）を達成しました。
• 不完全データ: 連続したデータ欠損やランダムな欠損がある場合でも、状態変数の軌道を高精度に再構成できました。
• 非線形項の発見: 指数関数項を含む Colpitts 発振器において、ノイズ下でも正しい非線形項とそのパラメータを正確に発見しました。
• 計算効率: 提案された LM 法は、他の最適化アルゴリズム（L-BFGS, 共役勾配法など）と比較して、10 倍以上高速に収束することが確認されました。

5. 意義と結論

この研究は、データ駆動型科学（Data-Driven Science）における重要な課題である「不完全でノイズの多いデータからのモデル発見」に対する強力な解決策を提供します。

• 実用性: 物理的な制約や専門家の知識を前提とせず、一般的な多項式ライブラリからでも、統計的に正当化された手法で高精度なモデルを自動発見できます。
• 拡張性: 常微分方程式だけでなく、偏微分方程式や微分代数方程式への拡張も可能であり、自動微分の利用により実装の柔軟性が高いです。
• 将来展望: 現在の手法は複雑な最適化ランドスケープに課題を残していますが、ハイパーパラメータの自動推定や、弱形式との組み合わせなど、さらなる発展の余地があります。

総じて、この論文は、観測データの質が低い現実的な環境下でも、物理法則（微分方程式）を正確に発見し、システムの挙動を再現できる画期的なフレームワークを提示しています。