BFGS法

BFGS法とは

数理最適化におけるBFGS（Broyden–Fletcher–Goldfarb–Shanno）法は、非制約非線形最適化問題を解くための反復アルゴリズムです。この方法は、目的関数の勾配情報のみを用いてヘッセ行列の近似を更新し、探索方向を決定します。DFP法と密接に関連しており、どちらも勾配のプレコンディショニングに曲率情報を用いる点が特徴です。BFGS法では、ヘッセ行列の逆行列を直接計算するのではなく、その推定値を更新することで計算コストを抑えます。

計算複雑度

BFGS法の計算複雑度はO(n^2)であり、ヘッセ行列を直接計算するニュートン法のO(n^3)よりも高速です。これは、BFGS法がヘッセ行列の逆行列の評価を必要としないためです。さらに、メモリ使用量を抑えたL-BFGS法も広く利用されており、多数の変数（1000以上）を扱う問題に適しています。また、BFGS-B法はシンプルなボックス制約に対応できます。

理論的背景

目的関数f(x)を最小化する最適化問題を考えます。ここで、xはn次元ベクトルであり、f(x)は微分可能なスカラー値関数とします。BFGS法は、初期推定値x0から開始し、反復ごとに推定値を更新します。

各反復kにおける探索方向pkは、以下の式で求められます。

Bkp_k = -∇f(x_k)

ここで、Bkはxkにおけるヘッセ行列の推定値です。このBkは、各反復で目的関数の勾配∇f(x_k)を用いて更新されます。

探索方向pkが得られたら、この方向に沿って直線探索を行い、f(x_k + γp_k)を最小化するスカラーγを決定し、次の点x_k+1を求めます。

準ニュートン条件

Bkの更新には、以下の準ニュートン条件が適用されます。

B_k+1(x_k+1 - x_k) = ∇f(x_k+1) - ∇f(x_k)

ここで、

y_k = ∇f(x_k+1) - ∇f(x_k)
s_k = x_k+1 - x_k

と定義すると、Bk+1は次の正割方程式を満たします。

B_k+1s_k = y_k

曲率条件

行列B_k+1が正定値であるためには、曲率条件s_k^T y_k > 0を満たす必要があります。この条件は、目的関数が強凸関数でない場合には明示的に課す必要があり、たとえば、直線探索でウルフ条件を満たす点を選べばこの条件を満足できます。

ヘッセ行列の更新

BFGS法では、各反復ごとにヘッセ行列を再計算する代わりに、前回の推定値に2つの行列を加算することで更新を行います。

B_k+1 = B_k + U_k + V_k

ここで、U_kとV_kはどちらも階数1の対称行列です。この更新により、階数2の対称行列を用いた更新が行われます。より単純な対称ランク1法では、階数1の行列を用いて更新を行いますが、正定値性は保証されません。

正割条件を課すと、更新式は以下のようになります。

B_k+1 = B_k + (y_k y_k^T) / (y_k^T s_k) - (B_k s_k s_k^T B_k^T) / (s_k^T B_k s_k)

アルゴリズム

非線形関数f: R^n → R の非制約最適化問題を考えます。

1. 初期推定解x0と初期ヘッセ行列の推定値B0を設定します。
2. 各反復kにおいて、次のステップを実行します。
1. 降下方向pkをBkp_k = -∇f(x_k)を解いて求めます。
2. 直線探索により、ステップサイズαkを決定します。この際、厳密な直線探索でなく、ウルフ条件を満たす非厳密な直線探索で十分な場合が多いです。
3. s_k = α_k p_k とし、x_k+1 = x_k + s_k により推定解を更新します。
4. y_k = ∇f(x_k+1) - ∇f(x_k)を計算します。
5. ヘッセ行列の推定値を更新します。

B_k+1 = B_k + (y_k y_k^T) / (y_k^T s_k) - (B_k s_k s_k^T B_k^T) / (s_k^T B_k s_k)

3. 勾配のノルム||∇f(x_k)||が十分に小さくなった場合、収束したとみなし、終了します。

最初のステップでは、B0を単位行列Iとすると、最急降下法と等価になります。しかし、その後のステップでは、Bkがヘッセ行列を推定することにより、収束が改善されます。

ヘッセ行列の逆行列の推定

実際には、ヘッセ行列の逆行列H_kを直接推定する方が効率的な場合があります。逆行列の更新にはシャーマン・モリソンの公式を利用できます。

H_k+1 = (I - (s_k y_k^T) / (y_k^T s_k)) H_k (I - (y_k s_k^T) / (y_k^T s_k)) + (s_k s_k^T) / (y_k^T s_k)

この計算は、B_k^-1が対称行列であること、y_k^T B_k^-1 y_kとs_k^T y_kがスカラーであることを利用して、一時的な行列を必要とせずに実行できます。

逆行列の推定ステップ

1. 初期推定解x0とヘッセ行列の逆行列の推定値H0を設定します。
2. 各反復kにおいて、次のステップを実行します。
1. 降下方向pkをp_k = -H_k∇f(x_k)で求めます。
2. 直線探索により、ステップサイズαkを決定します。
3. s_k = α_k p_k とし、x_k+1 = x_k + s_k により推定解を更新します。
4. y_k = ∇f(x_k+1) - ∇f(x_k)を計算します。
5. ヘッセ行列の逆行列の推定値を更新します。

H_k+1 = H_k + ((s_k^T y_k + y_k^T H_k y_k)(s_k s_k^T)) / ((s_k^T y_k)^2) - (H_k y_k s_k^T + s_k y_k^T H_k) / (s_k^T y_k)

発展

BFGS法の更新公式は、曲率s_k^T y_kが常に正で、ゼロから離れた下界があることに依存しています。しかし、実際の問題では、負やほぼゼロの曲率が現れることがあります。このような場合には、減衰BFGS更新など、より頑健な更新式が用いられます。

実装

BFGS法は、様々なライブラリで実装されています。

ALGLIB: C++およびC#用の実装
GNU Octave: 信頼領域を用いたBFGS法の実装
GSL: C言語での実装
R言語: optim()関数で利用可能
SciPy: scipy.optimize.fmin_bfgs関数で利用可能
Julia: Optim.jlパッケージで利用可能
Artelys Knitro: 大規模非線形最適化ソフトウェア
MATLAB Optimization Toolbox: fminunc関数で利用可能
Mathematica: BFGS法を含む
LS-DYNA: 陰解法に利用

出典

Avriel, Mordecai (2003), Nonlinear Programming: Analysis and Methods, Dover Publishing
Bonnans, J. Frédéric; Gilbert, J. Charles; Lemaréchal, Claude; Sagastizábal, Claudia A. (2006), “Newtonian Methods”, Numerical Optimization: Theoretical and Practical Aspects (Second ed.), Berlin: Springer
Fletcher, Roger (1987), Practical Methods of Optimization (2nd ed.), New York: John Wiley & Sons
Luenberger, David G.; Ye, Yinyu (2008), Linear and nonlinear programming, International Series in Operations Research & Management Science, 116 (Third ed.), New York: Springer
Kelley, C. T. (1999), Iterative Methods for Optimization, Philadelphia: Society for Industrial and Applied Mathematics

外部リンク

* Source code of high-precision BFGS

もう一度検索