統計学輪講 第24回

説明変数が何らかの微分方程式で拘束されている目的関数の微分を計算することは、ハミルトニアンモンテカルロ法・ニューラルネットワーク・４次元変分法データ同化などの統計関連分野に限らず、航空・工学分野のトポロジー最適化など、様々な科学分野で広く必要とされる。特にそのような目的関数の２階微分行列であるヘッセ行列は、非線形共役勾配法やニュートン法などの高速な最適化法に利用されるほか、逆行列要素は最適解の不確実性の評価にも直結するため重要であるが、計算コストが大きいため、大自由度系に対してもヘッセ行列を高速かつ高精度に計算する枠組みが希求されている。これまでに我々は、ヘッセ行列およびその逆行列要素の高速・高精度計算を目的として、４次元変分法データ同化で用いられる 2nd-order adjoint (SOA) 法に基づくアルゴリズム開発およびその応用研究を推進してきた [1,2]。SOA 法で利用される SOA モデルは一般に常微分方程式の形式で与えられ、ヘッセ行列評価には SOA モデルの数値積分が必要になるが、その数値積分法の選び方によっては必要なメモリが増大するだけでなく、ヘッセ行列の計算精度が著しく低下しそれに基づいて計算される不確実性などの結果が信頼できないものになる可能性があるなどの問題があり、省メモリ化およびヘッセ行列計算精度の担保のためにどのような数値積分法を選択すべきかの指針を数理学的な観点から与えることは解決すべき積年の課題であった。
このような背景から我々は、近年提案された Sanz-Serna (2016) [3] の理論に基づいて、必要なメモリを最低限に抑えさらにヘッセ行列の数値誤差を計算機誤差まで抑えることを可能にする SOA モデルの最適な数値積分法の選択法を提案した。本手法では SOA 法に登場する微分方程式群に内在する保存量を離散化後も保存するような数値積分法を構築することで高精度なヘッセ行列計算を可能にしている。反応拡散系や波動方程式系の初期値推定問題やパラメータ推定問題などを通じて本手法を検証したところ、本手法から提案される数値積分法は従来用いられてきた数値積分法に比べて、ヘッセ行列に含まれる数値誤差を劇的に抑えることができていることが確かめられた。講演では、これら適用事例結果を本手法の導出と合わせて紹介する。

※ 本研究は、理化学研究所脳神経科学研究センターの松田孟留先生、大阪大学サイバーメディアセンターの宮武勇登先生との共同研究です [4]。

参考文献：
[1] SI, H. Nagao, A. Yamanaka, Y. Tsukada, T. Koyama, M. Kano and J. Inoue, Data assimilation for massive autonomous systems based on second-order adjoint method, Phys. Rev. E 94, 043307 (2016).
[2] SI, H. Nagao, T. Kasuya, and J. Inoue, Grain growth prediction based on data assimilation by implementing 4DVar on multi-phase-field model, Sci. Tech. Adv. Mater. 18:1, 857-869 (2017).
[3] J. M. Sanz-Serna, Symplectic Runge--Kutta schemes for adjoint equations, automatic differentiation, optimal control, and more, SIAM Rev. 58(1), 3–33 (2016).
[4] SI, T. Matsuda, and Y. Miyatake, Adjoint-based exact Hessian computation, BIT Numerical Mathematics (2020), in press. (arXiv:1910.06524)