electricpockyの日記

電気の人です。生まれ変わったらピカチュウになりたい。嘘です。

中心差分近似によるヘッセ行列の計算方法

1.はじめに

最適化演算をニュートン法を適用して解く場合、目的関数の2階微分で構成されたヘッセ行列を使って解の更新ベクトルを算出する必要がある。本記事では、このヘッセ行列の各要素を差分近似の一つである中心差分近似を使って求める方法を記載してみる。

 

2.中心差分による実装のメリット

ニュートン法等のヘッセ行列を使った最適化演算のプログラムを実装する際、ヘッセ行列の計算を差分近似で実装すれば、目的関数に関係なく各行列要素を求めることが可能。

波動方程式ポアソン方程式といった偏微分方程式変分法で解いている資料は数多くあるが、変分法の対象として扱われている簡易なパターンだけを考えれば、これらの方程式は∂x^2や∂t^2といった形であり、ヘッセ行列で必要とされる∂x1∂x2のような複数の変数による微分への考慮は不要である。

このような背景から、中心差分近似を使ったヘッセ行列の計算の実装方法を形に残しておこうと思った次第です。

 

3.中心差分による実装方法

まず、前進差分、後退差分、中心差分の基礎的な解説は、先人が数多くの媒体で解説しているので割愛します。

あくまでプログラム実装に使える記事とするべく、結論(算出式)を先出しし、算出式の導出は次点に設ける構成にしてみました。

 

f:id:electricpocky:20210131014514j:plain

f:id:electricpocky:20210131014526j:plain

f:id:electricpocky:20210131014539j:plain

f:id:electricpocky:20210131014551j:plain

f:id:electricpocky:20210131014609j:plain

f:id:electricpocky:20210131014619j:plain

4.注意事項

 この内容はプログラムの実装や直感的な理解を目的に記載しているので、一部数学的な厳密さを欠いた記載となっております。そのため、数学的な議論もふまえた厳密な理解を目指す人や、ヘッセ行列を用いた新たなアルゴリズムを検討している人は、ちゃんとした文献を読んで勉強してください。また、本記事の内容に起因して発生した損害については一切責任を負いかねます。

 

信頼領域法についてわかったこと【Dogleg法ほか】

信頼領域法について調べたので、忘れ防止でメモ。

もう少し深堀したいが、まずはアウトプットし、

適宜修正する方向でやっていきたい。

 

1.信頼領域法とは

  • 暫定解を中心とした半径を定義(以後、信頼半径とする)し、この信頼半径内を二次モデル関数で近似することで、探索方向を決定する最適化手法である。
  • 決定した探索方向に対し、二次モデル関数の減少量と、目的関数の減少量の比較により、信頼半径を修正しながら、最適解を求める。

  • 最急降下法ニュートン法と同様、反復計算により解を求める

 

2.信頼領域法を実装するための手法(2/9 一部修正)

  • 二次モデル関数の最小化問題は正確に解くことが難しいため、近似的に解を求める「dogleg法」等の手法が主に使われる。
  • 「dogleg法」の他には、「Conjugated Gradient Steihaug法」と制約なし非線形最小二乗問題の解法である「Levenberg-Marquardt法」を用いて解く方法もある。
  • 特に「Levenberg-Marquardt法」を用いて解く方法は「Nearly exact trust-region法」と呼ばれている。
  • Conjugated Gradient Steihaug法とNearly exact trust-region法は、信頼半径の修正による反復に加えて、二次モデル関数の最適化計算自体の反復計算も必要な手法であり、目的関数や初期値に応じて計算時間が増大するリスクあり。そのため、計算時間で考えると「dogleg法」で実装する方が無難。
  • 信頼領域法はpythonのscipy.optimize.minimizeにより実装可能

 

3.信頼領域法のアルゴリズム

f:id:electricpocky:20210209094423p:plainf:id:electricpocky:20210209094435p:plainf:id:electricpocky:20210209094450p:plainf:id:electricpocky:20210209094502p:plainf:id:electricpocky:20210209094515p:plainf:id:electricpocky:20210209094527p:plainf:id:electricpocky:20210209094651p:plain
f:id:electricpocky:20210209094609p:plainf:id:electricpocky:20210209094619p:plainf:id:electricpocky:20210209094631p:plainf:id:electricpocky:20210209094725p:plainf:id:electricpocky:20210209094835p:plainf:id:electricpocky:20210209094747p:plainf:id:electricpocky:20210209094759p:plain















 

 

 

 

 

 

 

 

 



4.注意事項

 この内容はプログラムの実装や直感的な理解を目的に記載しているので、一部数学的な厳密さを欠いています。そのため、数学的な議論もふまえた厳密な理解を目指す人や、信頼領域法を応用したアルゴリズムを検討している人は、ちゃんとした文献を読んで勉強してください。また、本記事の内容に起因して発生した損害については一切責任を負いかねます。

 

5.参考文献

[1] Jorge Nocedal and Stephen J. Wright "Numerical Optimization" Springer-Verlag New York, Inc. 1999.

[2] 矢部 博「工学基礎 最適化とその応用」 数理工学社、2006年。

[3] http://www.ams.jhu.edu/~abasu9/AMS_553-761/lecture06_handout4.pdf

 (ジョンズホプキンス大学の講義資料 2021年2月3日閲覧)

 

【2/5 追記】

・lm法によるλの更新式の導出を追記

・誤記を修正

・数値例の計算誤りを修正

 

【2/9 追記】

・Conjugated Gradient Steihaug法の例題を追記

・誤記、および一部表現を修正

 

表現の修正は適宜実行するとして、以後の内容の拡充は別記事にて扱います。

気が向いたらマージするかも。

 

はじめまして【Numpyを使ったベクトルと行列の演算】

初投稿。

どうもelectricpockyです。電気を生業としてる人です。

数学とかデータサイエンスに興味があり、自身の勉強の成果を備忘録として残すため、ブログを作りました。

見やすい、すぐに内容を把握できる、すぐに使える備忘録を作るべく奮闘しようと思います。

 

自己紹介はおいおい実施するとして、pythonの拡張モジュールであるNumpyを使ったベクトルと行列の演算について、自分なりに整理してみました。

 

学生時代にC言語でプログラムを書いてたとき、行列演算をfor文使ってせこせこ実装してたけど、pythonだとNumpy使えばすごく楽に書けるんですよね。

np.dotを使えば大体の演算ができる感じ。

いちいち調べるの面倒なので、一つの画像に整理。

Input、Outputとしてpythonのコードを記載することで、数値例との結びつきを意識してみました。

何か困ったときの字引的な感じで使えると良いなと思う。

 

f:id:electricpocky:20210115235415j:plain