ブロイデン・フレッチャー・ゴールドファーブ・シャノアルゴリズム

数値 最適化において、ブロイデン・フレッチャー・ゴールドファーブ・シャノ（BFGS）アルゴリズムは、制約のない非線形最適化問題を解く反復法です。^[1]関連するデイビッドン・フレッチャー・パウエル法と同様に、BFGSは曲率情報で勾配を前処理することで降下方向を決定します。これは、一般化セカント法による勾配評価（または近似勾配評価）のみから得られる損失関数のヘッセ行列の近似値を徐々に改善することで行われます。^[2]

BFGS曲率行列の更新には逆行列が必要ないため、ニュートン法と比較して計算量はだけです。また、一般的に使用されているのはL-BFGSです。これはBFGSのメモリ制限バージョンで、非常に多数の変数（たとえば、> 1000）の問題に特に適しています。BFGS-Bバリアントは、単純なボックス制約を処理します。^[3]。BFGS行列はコンパクトな表現も許容するため、大規模な制約付き問題に適しています。 ${\displaystyle {\mathcal {O}}(n^{2})}$${\displaystyle {\mathcal {O}}(n^{3})}$

このアルゴリズムは、チャールズ・ジョージ・ブロイデン、ロジャー・フレッチャー、ドナルド・ゴールドファーブ、デビッド・シャノにちなんで名付けられました。^{[4] [5] [6] [7]}

最適化問題は、 を最小化することです。ここで、は のベクトルであり、 は微分可能なスカラー関数です。が取ることができる値には制約はありません。 ${\displaystyle f(\mathbf {x} )}$${\displaystyle \mathbf {x} }$${\displaystyle \mathbb {R} ^{n}}$${\displaystyle f}$${\displaystyle \mathbf {x} }$

アルゴリズムは最適値の初期推定から始まり、各段階でより良い推定値を得るために反復的に進行します。 ${\displaystyle \mathbf {x} _{0}}$

ステージkでの検索方向 p _{k は}、ニュートン方程式の類似体の解によって与えられます。

${\displaystyle B_{k}\mathbf {p} _{k}=-\nabla f(\mathbf {x} _{k}),}$

ここで、は各段階で繰り返し更新される におけるヘッセ行列の近似値であり、はx kで評価さ_れる関数の勾配である_。次に、p _k方向の直線探索を使用して、_スカラー${\displaystyle B_{k}}$${\displaystyle \mathbf {x} _{k}}$${\displaystyle \nabla f(\mathbf {x} _{k})}$${\displaystyle f(\mathbf {x} _{k}+\gamma \mathbf {p} _{k})}$${\displaystyle \gamma >0.}$

の更新に課される準ニュートン条件は ${\displaystyle B_{k}}$

${\displaystyle B_{k+1}(\mathbf {x} _{k+1}-\mathbf {x} _{k})=\nabla f(\mathbf {x} _{k+1})- \nbla f(\mathbf {x} _{k}).}$

とすると、は次を満たす 。${\displaystyle \mathbf {y} _{k}=\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}$${\displaystyle \mathbf {s} _{k}=\mathbf {x} _{k+1}-\mathbf {x} _{k}}$${\displaystyle B_{k+1}}$

${\displaystyle B_{k+1}\mathbf {s} _{k}=\mathbf {y} _{k}}$、

これは正割方程式です。

が正定値となるためには曲率条件が満たされている必要があり、これは割線方程式を であらかじめ乗算することで確認できます。関数が強く凸でない場合は、直線探索を使用して曲率条件を伴うWolfe条件を満たす点x _{k +1}を見つけるなどして、条件を明示的に強制する必要があります。 ${\displaystyle \mathbf {s} _{k}^{\top }\mathbf {y} _{k}>0}$${\displaystyle B_{k+1}}$${\displaystyle \mathbf {s} _{k}^{T}}$

点における完全なヘッセ行列を として計算する必要はなく、ステージkにおける近似ヘッセ行列は2 つの行列を追加することで更新されます。 ${\displaystyle \mathbf {x} _{k+1}}$${\displaystyle B_{k+1}}$

${\displaystyle B_{k+1}=B_{k}+U_{k}+V_{k}.}$

と はどちらも対称ランク 1 行列ですが、それらの合計はランク 2 更新行列になります。BFGS とDFP更新行列はどちらも、ランク 2 行列によってその前身と異なります。もう 1 つのより単純なランク 1 法は対称ランク 1法として知られており、これは正定値性を保証しません。 の対称性と正定値性を維持するために、更新形式を として選択できます。セカント条件 を課します。 と を選択すると、次が得られます。^[8]${\displaystyle U_{k}}$${\displaystyle V_{k}}$${\displaystyle B_{k+1}}$${\displaystyle B_{k+1}=B_{k}+\alpha \mathbf {u} \mathbf {u} ^{\top }+\beta \mathbf {v} \mathbf {v} ^{\top }}$${\displaystyle B_{k+1}\mathbf {s} _{k}=\mathbf {y} _{k}}$${\displaystyle \mathbf {u} =\mathbf {y} _{k}}$${\displaystyle \mathbf {v} =B_{k}\mathbf {s} _{k}}$

${\displaystyle \alpha ={\frac {1}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}},}$

${\displaystyle \beta =-{\frac {1}{\mathbf {s} _{k}^{T}B_{k}\mathbf {s} _{k}}}.}$

最後に、とをに代入して、の更新方程式を取得します。 ${\displaystyle \alpha }$${\displaystyle \beta }$${\displaystyle B_{k+1}=B_{k}+\alpha \mathbf {u} \mathbf {u} ^{\top }+\beta \mathbf {v} \mathbf {v} ^{\top }}$${\displaystyle B_{k+1}}$

${\displaystyle B_{k+1}=B_{k}+{\frac {\mathbf {y} _{k}\mathbf {y} _{k}^{\mathrm {T} }}{\mathbf {y} _{k}^{\mathrm {T} }\mathbf {s} _{k}}}-{\frac {B_{k}\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} }B_{k}^{\mathrm {T} }}{\mathbf {s} _{k}^{\mathrm {T} }B_{k}\mathbf {s} _{k}}}.}$

非線形目的関数 である次の制約のない最適化問題を考えます 。$${\displaystyle {\begin{aligned}{\underset {\mathbf {x} \in \mathbb {R} ^{n}}{\text{minimize}}}\quad &f(\mathbf {x} ),\end{aligned}}}$$${\displaystyle f:\mathbb {R} ^{n}\to \mathbb {R} }$

初期推定値とヘッセ行列の初期推定値から、解に収束する まで次の手順が繰り返されます。${\displaystyle \mathbf {x} _{0}\in \mathbb {R} ^{n}}$${\displaystyle B_{0}\in \mathbb {R} ^{n\times n}}$${\displaystyle \mathbf {x} _{k}}$

1. を解いて方向を取得します。${\displaystyle \mathbf {p} _{k}}$${\displaystyle B_{k}\mathbf {p} _{k}=-\nabla f(\mathbf {x} _{k})}$
2. 1 次元の最適化 (直線探索) を実行して、最初のステップで見つかった方向の許容可能なステップ サイズを見つけます。正確な直線探索を実行すると、 になります。実際には、通常は不正確な直線探索で十分であり、許容可能な はWolfe 条件を満たします。${\displaystyle \alpha _{k}}$${\displaystyle \alpha _{k}=\arg \min f(\mathbf {x} _{k}+\alpha \mathbf {p} _{k})}$${\displaystyle \alpha _{k}}$
3. 設定して更新します。${\displaystyle \mathbf {s} _{k}=\alpha _{k}\mathbf {p} _{k}}$${\displaystyle \mathbf {x} _{k+1}=\mathbf {x} _{k}+\mathbf {s} _{k}}$
4. ${\displaystyle \mathbf {y} _{k}={\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}}$。
5. ${\displaystyle B_{k+1}=B_{k}+{\frac {\mathbf {y} _{k}\mathbf {y} _{k}^{\mathrm {T} }}{\mathbf {y} _{k}^{\mathrm {T} }\mathbf {s} _{k}}}-{\frac {B_{k}\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} }B_{k}^{\mathrm {T} }}{\mathbf {s} _{k}^{\mathrm {T} }B_{k}\mathbf {s} _{k}}}}$。

収束は勾配のノルムを観察することで判定できます。ある が与えられた場合、 のときにアルゴリズムを停止できます。がで初期化された場合、最初のステップは勾配降下法と同等になりますが、それ以降のステップはヘッセ行列の近似値 である によってますます洗練されていきます。${\displaystyle \epsilon >0}$${\displaystyle ||\nabla f(\mathbf {x} _{k})||\leq \epsilon .}$${\displaystyle B_{0}}$${\displaystyle B_{0}=I}$${\displaystyle B_{k}}$

アルゴリズムの最初のステップは行列の逆行列を使って実行される。これは、アルゴリズムのステップ5に シャーマン・モリソンの公式を適用することで効率的に得られ、${\displaystyle B_{k}}$

${\displaystyle B_{k+1}^{-1}=\left(I-{\frac {\mathbf {s} _{k}\mathbf {y} _{k}^{T}}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}}\right)B_{k}^{-1}\left(I-{\frac {\mathbf {y} _{k}\mathbf {s} _{k}^{T}}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}}\right)+{\frac {\mathbf {s} _{k}\mathbf {s} _{k}^{T}}{\mathbf {y} _{k}^{T}\mathbf {s} _{k}}}.}$

これは、が対称であり、およびがスカラーであることを認識し、次のような展開を使用して、 一時的な行列を使わずに効率的に計算できます。${\displaystyle B_{k}^{-1}}$${\displaystyle \mathbf {y} _{k}^{\mathrm {T} }B_{k}^{-1}\mathbf {y} _{k}}$${\displaystyle \mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}}$

${\displaystyle B_{k+1}^{-1}=B_{k}^{-1}+{\frac {(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}+\mathbf {y} _{k}^{\mathrm {T} }B_{k}^{-1}\mathbf {y} _{k})(\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} })}{(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k})^{2}}}-{\frac {B_{k}^{-1}\mathbf {y} _{k}\mathbf {s} _{k}^{\mathrm {T} }+\mathbf {s} _{k}\mathbf {y} _{k}^{\mathrm {T} }B_{k}^{-1}}{\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}}}.}$

したがって、行列の反転を避けるために、ヘッセ行列そのものではなく、ヘッセ行列の逆行列を近似することができる。 ^[9]${\displaystyle H_{k}{\overset {\operatorname {def} }{=}}B_{k}^{-1}.}$

初期推定値と近似逆ヘッセ行列から、解に収束する まで次の手順を繰り返します。${\displaystyle \mathbf {x} _{0}}$${\displaystyle H_{0}}$${\displaystyle \mathbf {x} _{k}}$

1. を解いて方向を取得します。${\displaystyle \mathbf {p} _{k}}$${\displaystyle \mathbf {p} _{k}=-H_{k}\nabla f(\mathbf {x} _{k})}$
2. 1 次元の最適化 (直線探索) を実行して、最初のステップで見つかった方向の許容可能なステップ サイズを見つけます。正確な直線探索を実行すると、 になります。実際には、通常は不正確な直線探索で十分であり、許容可能な はWolfe 条件を満たします。${\displaystyle \alpha _{k}}$${\displaystyle \alpha _{k}=\arg \min f(\mathbf {x} _{k}+\alpha \mathbf {p} _{k})}$${\displaystyle \alpha _{k}}$
3. 設定して更新します。${\displaystyle \mathbf {s} _{k}=\alpha _{k}\mathbf {p} _{k}}$${\displaystyle \mathbf {x} _{k+1}=\mathbf {x} _{k}+\mathbf {s} _{k}}$
4. ${\displaystyle \mathbf {y} _{k}={\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}}$。
5. ${\displaystyle H_{k+1}=H_{k}+{\frac {(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}+\mathbf {y} _{k}^{\mathrm {T} }H_{k}\mathbf {y} _{k})(\mathbf {s} _{k}\mathbf {s} _{k}^{\mathrm {T} })}{(\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k})^{2}}}-{\frac {H_{k}\mathbf {y} _{k}\mathbf {s} _{k}^{\mathrm {T} }+\mathbf {s} _{k}\mathbf {y} _{k}^{\mathrm {T} }H_{k}}{\mathbf {s} _{k}^{\mathrm {T} }\mathbf {y} _{k}}}}$。

統計的推定問題（最大尤度法やベイズ推定など）では、解の信頼区間または信頼区間は、最終的なヘッセ行列の逆行列から推定することができます^[要出典]。しかし、これらの量は技術的には真のヘッセ行列によって定義され、BFGS近似は真のヘッセ行列に収束しない可能性があります。^[10]

BFGS 更新式は、曲率が厳密に正で、ゼロから制限されていることに大きく依存しています。この条件は、凸ターゲットに対して Wolfe 条件で直線探索を実行すると満たされます。ただし、実際のアプリケーション (逐次二次計画法など) では、負またはほぼゼロの曲率が定期的に生成されます。これは、非凸ターゲットを最適化するときや、直線探索の代わりに信頼領域アプローチを採用するときに発生することがあります。ターゲットのノイズにより、誤った値が生成される可能性もあります。 ${\displaystyle \mathbf {s} _{k}^{\top }\mathbf {y} _{k}}$

このような場合には、より堅牢な更新を得るために 、いわゆる減衰BFGS更新の1つを使用することができる（^[11]を参照） 。${\displaystyle \mathbf {s} _{k}}$${\displaystyle \mathbf {y} _{k}}$

注目すべきオープンソース実装は次のとおりです。

• ALGLIBはBFGSとそのメモリ制限バージョンをC++とC#で実装します。
• GNU Octave はfsolve、信頼領域拡張を備えたBFGS の形式を関数内で使用します。
• GSLはBFGSをgsl_multimin_fdfminimizer_vector_bfgs2として実装している。^[12]
• Rでは、BFGSアルゴリズム（およびボックス制約を許可するL-BFGS-Bバージョン）は、基本関数optim()のオプションとして実装されています。^[13]
• SciPyでは、scipy.optimize.fmin_bfgs関数がBFGSを実装しています。^[14]また、パラメータLを非常に大きな数に設定することで、L-BFGSアルゴリズムのいずれかを使用してBFGSを実行することも可能です。
• Juliaでは、Optim.jlパッケージがBFGSとL-BFGSをoptimize()関数のソルバーオプションとして実装しています（他のオプションも含まれています）。^[15]

注目すべき独自の実装には次のものがあります。

• 大規模非線形最適化ソフトウェアArtelys Knitro は、BFGS アルゴリズムと L-BFGS アルゴリズムの両方を実装しています。
• MATLAB Optimization Toolboxでは、問題のサイズが「中規模」に設定されている場合に、 fminunc 関数は 3 次直線探索による BFGS を使用します。
• Mathematica にはBFGS が含まれています。
• LS-DYNA は暗黙的な問題を解決するために BFGS も使用します。

• アヴリエル、モーデカイ（2003）、非線形計画法：分析と方法、ドーバー出版、ISBN 978-0-486-43227-4
• Bonnans, J. Frédéric; Gilbert, J. Charles; Lemaréchal, Claude ; Sagastizábal, Claudia A. (2006)、「ニュートン法」、数値最適化: 理論的および実践的側面(第 2 版)、ベルリン: Springer、pp. 51–66、ISBN 3-540-35445-X
• フレッチャー、ロジャー（1987）、最適化の実践的方法（第2版）、ニューヨーク：ジョン・ワイリー＆サンズ、ISBN 978-0-471-91547-8
• Luenberger, David G. ; Ye, Yinyu (2008)、線形および非線形計画法、International Series in Operations Research & Management Science、第 116 巻 (第 3 版)、ニューヨーク: Springer、pp. xiv+546、ISBN 978-0-387-74502-2、MR  2423726
• ケリー、CT（1999）、反復法による最適化、フィラデルフィア：産業応用数学協会、pp. 71-86、ISBN 0-89871-433-8