BFGSメソッドまたは最も効果的な最適化メソッドの1つ。 Pythonの実装例

数値最適化の反復法であるBFGSメソッドは、その研究者にちなんで名付けられました： Bロイデン、 Fレッチャー、 Gオールドファーブ、 Sハノ。 いわゆる準ニュートン法のクラスに属します。 ニュートン法とは対照的に、準ニュートン法はヘッセ関数を直接計算しません。 二次偏導関数を見つける必要はありません。 代わりに、ヘッセ行列はこれまでに実行された手順に基づいて近似的に計算されます。



メソッドにはいくつかの変更があります。

L-BFGS （メモリ使用量の制限）-多数の未知の場合に使用されます。

L-BFGS-Bは、多次元キューブでのメモリ使用量が制限された変更です。



この方法は効果的で安定しているため、最適化機能でよく使用されます。 たとえば、Python言語の一般的なライブラリであるSciPyでは、最適化機能はデフォルトでBFGS、L-BFGS-Bを使用します。

アルゴリズム

何らかの機能を与えましょう $$$f（x、y）$ そして最適化問題を解きます： $$$min f（x、y）$ 。

一般的な場合の場所 $$$f（x、y）$ は、連続的な二次導関数を持つ非凸関数です。



ステップ1

開始点を初期化する $$$x_0$ ;

検索の精度を設定します $$$ε$ > 0;

初期近似を決定します $$$H_0 = B_0 ^ {-1}$ どこで $$$B_0 ^ {-1}$ -逆ヘッセ関数;



初期近似の選択方法 $$$H_0$ ？

残念ながら、すべての場合にうまく機能する一般的な公式はありません。 初期近似として、初期点で計算された関数のヘッセ行列を取ることができます $$$x_0$ 。 それ以外の場合、適切に調整された非縮退行列を使用できますが、実際には、多くの場合、単一の行列が使用されます。



ステップ番号2

検索する方向のポイントを見つけると、次のように決定されます。

$$

$$p_k = -H_k * \ nabla f_k$$

ステップ番号3

計算する $$$x_ {k + 1}$ 再帰関係を通じて：

$$

$$x_ {k + 1} = x_k +α_k* p_k$$

係数 $$$α_k$ 線形検索（行検索）を使用して見つける $$$α_k$ ウルフの条件を満たす ：

$$

$$f（x_k +α_k* p_k）\ leq f（x_k）+ c_1 *α_k* \ nabla f_k ^ T * p_k$$

$$

$$\ nabla f（x_k +α_k* p_k）^ T * p_k \ geq c_2 * \ nabla f_k ^ T * p_k$$

定数 $$$s_1$ そして $$$s_2$ 次のように選択します。 $$$0 \ leq c_1 \ leq c_2 \ leq 1$ 。 ほとんどの実装では： $$$c_1 = 0.0001$ そして $$$s_2 = 0.9$ 。



実際、これを見つけます $$$α_k$ 関数の値 $$$f（x_k +α_k* p_k）$ 最低限。



ステップ4

ベクトルを定義します。

$$

$$s_k = x_ {k + 1}-x_k$$

$$

$$y_k = \ nabla f_ {k + 1}-\ nabla f_k$$

$$$s_k$ -反復のアルゴリズムのステップ、 $$$y_k$ -反復で勾配を変更します。



ステップ番号5

次の式に従ってヘッセ関数を更新します。

$$

$$H_ {k + 1} =（I-ρ_k* s_k * y_k ^ T）H_k（I-ρ_k* y_k * s_k ^ T）+ρ* s_k * s_k ^ T$$

どこで $$$ρ_k$

$$

$$ρ_k= \ frac {1} {y_k ^ T s_k}$$

$$$I$ 単位行列です。

発言

式を見る $$$y_k * s_k ^ T$ 2つのベクトルの外積です 。

2つのベクトルを定義します $$$U$ そして $$$V$ その外積は行列積と同等です $$$UV ^ T$ 。 たとえば、平面内のベクトルの場合：

$$

$$UV ^ T = \ begin {pmatrix} U_1 \\ U_2 \ end {pmatrix} \ begin {pmatrix} V_1＆V_2 \ end {pmatrix} = \ begin {pmatrix} U_1V_1＆U_1V_2 \\ U_2V_1＆U_2V_2 \ end {pmatrix }$$

ステップ6

このアルゴリズムは、不等式が真になるまで実行され続けます。 $$$| \ nabla f_k | >ε$ 。

アルゴリズムの概略

例

次の関数の極値を見つけます。 $$$f（x、y）= x ^ 2-xy + y ^ 2 + 9x-6y + 20$

出発点として、 $$$x_0 =（1、1）$ ;

必要な精度 $$$ε= 0.001$ ;



関数の勾配を見つける $$$f（x、y）$ ：

$$

$$\ nabla f = \ begin {pmatrix} 2x-y + 9 \\ -x + 2y-6 \ end {pmatrix}$$

反復番号0 ：

$$

$$x_0 =（1、1）$$

ポイントで勾配を見つける $$$x_0$ ：

$$

$$\ nabla f（x_0）= \ begin {pmatrix} 10 \\ -5 \ end {pmatrix}$$

検索の終了を確認します。

$$

$$| \ nabla f（x_0）| = \ sqrt {10 ^ 2 +（-5）^ 2} = 11.18$$

$$

$$| \ nabla f（x_0）| = 11.18> 0.001$$

不等式は成り立たないため、反復プロセスを継続します。



検索する方向にポイントを見つけます。

$$

$$p_0 = -H_0 * \ nabla f（x_0）=-\ begin {pmatrix} 1＆0 \\ 0＆1 \ end {pmatrix} \ begin {pmatrix} 10 \\ -5 \ end {pmatrix} = \ begin {pmatrix} -10 \\ 5 \ end {pmatrix}$$

どこで $$$H_0 = I$ 単位行列です。



これを探しています $$$α_0$$$$\ geq$ 0から $$$f（x_0 +α_0* p_0）= min（f（x_0 +α_0* p_0））$



分析的に、問題は1つの変数の関数の極値を見つけることで軽減できます。

$$ $$表示$$ x_0 +α_0* p_0 =（1、1）+α_0（-10、5）=（1-10α_0、1 +5α_0）$$表示$$

それから



$$ $インライン$ f（α_0）=（1-10α_0）^ 2-（1-10α_0）（1 +5α_0）+（1 +5α_0）^ 2 + 9（1-10α_0）-6（1 +5α_0）+ 20 $インライン$



式を単純化すると、導関数が見つかります。

$$

$$\ frac {\ partial f} {\ partial x} =350α_0-125 = 0 \右矢印α_0= 0.357$$

次のポイントを見つける $$$x_1$ ：

$$

$$x_1 = x_0 +α_0* p_0 =（-2.571、2.786）$$

$$

$$s_0 = x_1-x_0 =（-2.571、2.786）-（1、1）=（-3.571、1.786）$$

tの勾配の値を計算します。 $$$x_1$ ：

$$

$$\ nabla f（x_1）= \ begin {pmatrix} 1.071 \\ 2.143 \ end {pmatrix}$$

$$

$$y_0 = \ nabla f（x_1）-\ nabla f（x_0）=（1.071、2.143）-（10、-5）=（-8.929、7.143）$$

ヘッセ近似を求めます。

$$

$$H_1 = \ begin {pmatrix} 0.694および0.367 \\ 0.367および0.709 \ end {pmatrix}$$

反復1：

$$

$$x_1 =（-2.571、2.786）$$

$$

$$\ nabla f（x_1）= \ begin {pmatrix} 1.071 \\ 2.143 \ end {pmatrix}$$

検索の終了を確認します。

$$

$$| \ nabla f（x_1）| = 2.396> 0.001$$

不等式は成り立たないため、反復プロセスを続けます。

$$

$$H_1 = \ begin {pmatrix} 0.694および0.367 \\ 0.367および0.709 \ end {pmatrix}$$

$$

$$p_1 = -H_1 \ nabla f（x_1）=-\ begin {pmatrix} 0.694＆0.367 \\ 0.367＆0.709 \ end {pmatrix} \ begin {pmatrix} 1.071 \\ 2.143 \ end {pmatrix} =-\ begin { pmatrix} 1.531 \\ 1.913 \ end {pmatrix}$$

見つける $$$α_1$ > 0：



$$ $インライン$ f（α_1）=-2.296α_1+（-1.913α_1+ 2.786）^ 2-（-1.913α_1+ 2.786）（-1.53​​1α_1-2.571）+（-1.53​​1α_1-2.571）^ 2-19.86 $ inline $

$$

$$\ frac {\ partial f} {\ partialα_1} =6.15α_1-5.74 = 0 \ Rightarrowα_1= 0.933$$

次に：

$$

$$x_2 =（-3.571、0.786）$$

$$

$$s_1 = x_2-x_1 =（-4、1）-（-2.571、2.786）=（-1.429、-1.786）$$

ポイントで勾配の値を計算する $$$x_2$ ：

$$

$$\ nabla f（x_2）= \ begin {pmatrix} 0 \\ 0 \ end {pmatrix}$$

$$

$$g_1 = \ nabla f（x_2）-\ nabla f（x_1）=（0、0）-（1.071、2.143）=（-1.071、-2.143）$$

$$

$$H_2 = \ begin {pmatrix} 0.667および0.333 \\ 0.333および0.667 \ end {pmatrix}$$

反復2 ：

$$

$$x_2 =（-4、1）$$

検索の終了を確認します。

$$

$$| \ nabla f（x_2）| = 0 <0.001$$

不等式が満たされるため、検索は終了します。 関数の極値を分析的に見つけます;それは $$$x_2$ 。

メソッドの詳細

各反復はコストをかけて行うことができます。 $$$O（n ^ 2）$ （さらに、関数の計算と勾配の推定のコスト）。 ここではない $$$O（n ^ 3）$ 線形システムの解法や複雑な数学演算などの操作。 このアルゴリズムは安定しており、超線形収束を備えており、ほとんどの実用的な問題には十分です。 ニュートンの手法がはるかに高速に（2次的に）収束しても、線形システムを解く必要があるため、各反復のコストは高くなります。 アルゴリズムの否定できない利点は、もちろん、二次導関数を計算する必要がないことです。



この方法が実際にうまく機能することは一般に受け入れられていますが、非凸平滑関数の場合のアルゴリズムの収束については、理論的にはほとんど知られていません。



BFGS式には自己修正プロパティがあります。 行列 $$$H$ 関数の曲率を誤って推定し、この不適切な推定によりアルゴリズムの速度が低下した場合、ヘッセ行列の近似により、数ステップで状況を修正しようとします。 アルゴリズムの自己修正プロパティは、対応する線形検索が実装されている場合にのみ機能します（Wolfe条件が満たされている）。

Pythonの実装例

アルゴリズムは、numpyおよびscipyライブラリを使用して実装されます。 線形検索は、scipy.optimizeモジュールのline_search（）関数を使用して実装されました。 この例では、反復回数に制限が課されています。

''' Pure Python/Numpy implementation of the BFGS algorithm. Reference: https://en.wikipedia.org/wiki/Broyden%E2%80%93Fletcher%E2%80%93Goldfarb%E2%80%93Shanno_algorithm ''' #!/usr/bin/python # -*- coding: utf-8 -*- import numpy as np import numpy.linalg as ln import scipy as sp import scipy.optimize # Objective function def f(x): return x[0]**2 - x[0]*x[1] + x[1]**2 + 9*x[0] - 6*x[1] + 20 # Derivative def f1(x): return np.array([2 * x[0] - x[1] + 9, -x[0] + 2*x[1] - 6]) def bfgs_method(f, fprime, x0, maxiter=None, epsi=10e-3): """ Minimize a function func using the BFGS algorithm. Parameters ---------- func : f(x) Function to minimise. x0 : ndarray Initial guess. fprime : fprime(x) The gradient of `func`. """ if maxiter is None: maxiter = len(x0) * 200 # initial values k = 0 gfk = fprime(x0) N = len(x0) # Set the Identity matrix I. I = np.eye(N, dtype=int) Hk = I xk = x0 while ln.norm(gfk) > epsi and k < maxiter: # pk - direction of search pk = -np.dot(Hk, gfk) # Line search constants for the Wolfe conditions. # Repeating the line search # line_search returns not only alpha # but only this value is interesting for us line_search = sp.optimize.line_search(f, f1, xk, pk) alpha_k = line_search[0] xkp1 = xk + alpha_k * pk sk = xkp1 - xk xk = xkp1 gfkp1 = fprime(xkp1) yk = gfkp1 - gfk gfk = gfkp1 k += 1 ro = 1.0 / (np.dot(yk, sk)) A1 = I - ro * sk[:, np.newaxis] * yk[np.newaxis, :] A2 = I - ro * yk[:, np.newaxis] * sk[np.newaxis, :] Hk = np.dot(A1, np.dot(Hk, A2)) + (ro * sk[:, np.newaxis] * sk[np.newaxis, :]) return (xk, k) result, k = bfgs_method(f, f1, np.array([1, 1])) print('Result of BFGS method:') print('Final Result (best point): %s' % (result)) print('Iteration Count: %s' % (k))
      
      

        
        
        
      

    
        
        
        
      
      

        
        
        
      

    
     

私の記事に関心をお寄せいただきありがとうございます。 彼女があなたにとって有益であり、あなたが多くを学んだことを願っています。



FUNNYDMANはあなたと一緒でした。 みなさんさようなら！