ドキュメンテーション

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

制約なし非線形最適化アルゴリズム

制約なしの最適化の定義

制約なしの最小化は、スカラー関数 f(x) の局所的最小値ベクトル x を見つける問題です。

minxf(x)

"制約なし" は x の範囲に制限がないことを意味します。

fminunc trust-region アルゴリズム

非線形最小化に対する信頼領域法

Optimization Toolbox™ のソルバーで使用される多くのメソッドは、"信頼領域法" を基にしています。信頼領域法はシンプルなものですが最適化では重要な概念です。

最適化の信頼領域法のアプローチを理解するために制約なし最小化問題を考え、f(x) を最小化します。ここで関数はベクトル引数を取り、スカラーを出力します。n 空間の点 x を想定し、より小さい関数値へ移動して最小化を行う場合を考えてみましょう。基本的な概念は、シンプルな関数 q で f を近似することです。この関数は、点 x の近傍 N で関数 f の挙動をよく表すものです。この近傍が信頼領域です。テスト ステップ s が N における最小化 (または近似最小化) によって計算されます。信頼領域の部分問題を次に示します。

mins{q(s), sN}.(6-1)

f(x + s) < f(x) の場合、現在の点が x + s に更新されます。そうでない場合は現在の点は変更されず、信頼領域 N は縮小され、テスト ステップの計算が繰り返されます。

f(x) を最小化するための信頼領域を決める上で重要な問題は、近似 q (現在の点 x で定義) の選択および計算方法、信頼領域 N の選択および変更方法、信頼領域の部分問題を解く精度です。この節では制約なしの問題に焦点をあてます。変数上に制約があるために複雑度が増す問題については、後節で取り扱います。

標準的な信頼領域法 ([48]) では二次近似 q は、x における F のテイラー近似のはじめの 2 項によって決められます。近傍 N は球形または楕円体です。数学的に、信頼領域の部分問題は、次のように表現できます。

min{12sTHs+sTg  such that  DsΔ},(6-2)

ここで g は現在の点 x における f の勾配です。H はヘッセ行列 (2 次導関数の対称行列) です。D は対角スケーリング行列であり、Δ は正のスカラーです。∥ .∥ は 2 ノルムです。式 6-2 を解くために適したアルゴリズムがあります ([48] を参照してください)。このようなアルゴリズムは、完全な固有システムと以下の永年方程式に適用されるニュートン過程の計算を一般的に含みます。

1Δ1s=0.

このようなアルゴリズムにより 式 6-2 の正確な解が与えられます。しかし、H の因子分解に比例して時間が必要になります。そのために、大規模問題に対しては種々のアプローチが必要となります。式 6-2 に基づく近似とヒューリスティックな方法は文献 ([42][50]) で提案されています。Optimization Toolbox のソルバーがサポートする近似アプローチとして、信頼領域法の部分問題を 2 次元の部分空間 S ([39][42]) に制限します。部分空間 S が計算されると、(部分空間では問題は 2 次元になるため) 完全な次元の固有値 / 固有ベクトルの情報が必要な場合でも 式 6-2 を解く作業は簡単になります。ここでの主な仕事は、部分空間を決定することになります。

2 次元の部分空間 S は、以下に示す前処理を使用した共役勾配法過程を用いて決められます。ソルバーは s1 と s2 で決まる線形空間として S を定義します。ここで、s1 は勾配 g の方向にあります。s2 は近似ニュートン方向、すなわち以下の解

Hs2=g,(6-3)

または負の曲率の方向です。

s2THs2<0.(6-4)

このように S を選択する背景の考え方は、最急降下方向または負の曲率方向にグローバルな収束を進めることと、ニュートン ステップが存在する場合は、これを介して迅速にローカルな収束を達成することです。

信頼領域法の概念を使用した制約なしの最小化の概要は以下になります。

  1. 2 次元の信頼領域の部分問題の定式化

  2. テスト ステップ s を決めるため、式 6-2 を解きます。

  3. f(x + s) < f(x) の場合、x = x + s とします。

  4. Δ を調節します。

これら 4 つのステップは、収束するまで繰り返されます。信頼領域の大きさ Δ は標準的な規則に基づいて調整されます。特に、使用するステップが適用されない場合、すなわち f(x + s) ≥ f(x) の場合、内点集合は小さくなります。詳細は [46][49] を参照してください。

Optimization Toolbox のソルバーは特定の関数 f の重要かつ特別なケースをいくつか扱います。非線形最小二乗法、二次関数、線形最小二乗法を考えてみましょう。しかし、根底に存在するアルゴリズムは、一般的な場合と同じです。これらの特別な場合は後続の節で説明します。

前処理付き共役勾配法

線形方程式系 Hp = –g の大きな対称正定システムを解く一般的な方法は、前処理付き共役勾配法 (PCG) です。この反復法は、形式 H·v の行列ベクトル積を計算する機能を必要とします。ここで v は任意のベクトルです。対称正定行列 M は H の "前提条件子" です。すなわち、M = C2 です。ここで C–1HC–1 は、条件数の良い行列または分類分けされた固有値をもつ行列です。

最小化の過程で、ヘッセ行列 H が対称と仮定します。しかし、H は、強い最小化子の近傍の中でのみ正定であることが保証されます。PCG アルゴリズムは、負または 0 の曲率の勾配をもつ場合に終了します (dTHd ≤ 0)。PCG 出力方向 p は、負の曲率方向またはニュートン システム Hp = –g への近似解のどちらかです (tol は近似方法を制御)。いずれにせよ、信頼領域法のアプローチで使用される 2 次元部分空間を定義するために p が使用されます (非線形最小化に対する信頼領域法を参照)。

fminunc quasi-newton アルゴリズム

制約なしの最適化の基礎

制約なしの最適化に対して、種々の方法が存在していますが、微係数情報を使うか、使わないかにより大きく分類されます。関数評価のみを使う探索法 (たとえば、Nelder and Mead [30] のシンプレックス探索) は、非線形な問題や不連続な問題に非常に有効です。勾配法は、一般に最小化される関数の 1 次導関数が連続型であれば非常に有効です。ニュートン法のような高次の項を使う方法は、2 次導関数の情報が既に用意されていたり、簡単に計算できるときにのみ有効です。これは、2 次導関数の情報を得るのに数値微分を使い、非常に多くの計算を必要とするためです。

勾配法は、最小値が存在すると考えられる探索の方向を示すために勾配の情報を使います。この中で最も簡単な方法は、最急降下法で ∇f(x) を目的関数の勾配とすると –∇f(x) の方向で探索を行うものです。この方法は最小化される関数が長く狭い谷をもつとき、あまり効果的なものではありません。たとえば、Rosenbrock 関数の場合です。

f(x)=100(x2x12)2+(1x1)2.(6-5)

この関数の最小値は f(x) = 0 となる x = [1,1] です。この関数の等高線図を下図に示します。点 [-1.9,2] の位置を初期推定値として最急降下法を使った場合の最小値を求めるまでの解の方向も示します。最適化は、1000 回の反復の後終了しますが、最小値の位置とかなり離れています。図で黒色で表される部分は、谷の片側へ行ったり来たりのジグザグ運動を連続的に行っていることを示しています。プロット図の中央に向かって、ある点が谷の中心に正確に到達するとき、より多くのステップ数が生じることに注意してください。

図 6-1. Rosenbrock 関数における最急降下法

このような関数は原点回りで曲率が歪んでいるために、制約なしの例の中で悪名高いものでバナナ関数として知られています。Rosenbrock 関数はさまざまな最適化手法の使い方を示すために、この節全体で使われています。等高線は、U 字谷の急斜面のために指数的な刻み幅でプロットされています。

この図をアニメーションで表示するには、MATLAB® コマンド ラインで「bandem」と入力します。

準ニュートン法

勾配情報を使う方法の中で、最も良く使われるものは準ニュートン法です。この方法は、各反復で曲率情報を構築し、次の二次モデルの問題を定式化します。

minx12xTHx+cTx+b,(6-6)

ここでヘッセ行列 H は正定対称行列で、c は定数ベクトル、b は定数です。この問題の最適解は x の偏導関数をゼロに近づけたときに求められます。すなわち以下のようになります。

f(x*)=Hx*+c=0.(6-7)

最適解 x* は次のように記述することができます。

x*=H1c.(6-8)

ニュートン タイプの方法は (準ニュートン法と違って) 直接 H を計算し、何回かの反復後、最小値に到達するために減少する方向に進みます。数値的に H を計算するには、非常に多くの計算を必要とします。準ニュートン法は観測された f(x) と ∇f(x) の動作から適切な更新法を使って H を近似するための曲率情報を構築することにより、計算時間が長くならないようにします。

多くのヘッシアン更新法が開発されてきました。しかしながら、Broyden [3]、Fletcher [12]、Goldfarb [20]、Shanno [37] の式 (BFGS 法) は、一般的な問題に使用するのに最も効果的と考えられています。

BFGS 法は次のように与えられます。

Hk+1=Hk+qkqkTqkTskHkskskTHkTskTHksk,(6-9)

ここで、

sk=xk+1xk,qk=f(xk+1)f(xk).

最適化を始める出発点 H0 は、対称な正定行列、たとえば、単位行列 I を設定してください。ヘッシアン H の逆行列を計算することを避けるため、更新法を導入します。この方法は各更新でヘッシアンの逆行列 H–1 の近似を求める式を使うことにより、H の逆行列を直接計算することを避けるものです。良く知られた方法として、Davidon [7]、Fletcher and Powell [14] の DFP 式があります。これは BFGS メソッド (式 6-9) と同じ式を使用しますが、qk が sk の代わりに使用される点が異なります。

勾配の情報は、解析的に計算された勾配から得られるか、または有限差分を使った数値微分法により偏微係数として求めるかのどちらかです。これは各設計変数 x に個々の変動を含ませ、目的関数中で変化率を計算するものです。

主要な各反復の k 番目でライン探索は次の方向になります。

d=Hk1f(xk).(6-10)

準ニュートン法は、Figure 6-2, Rosenbrock 関数の BFGS 法の適用にあるRosenbrock 関数の解を求める経路に示されています。方法は、谷の型に従って、有限差分勾配のみを使って 140 回の関数評価の後、最小値に収束します。

図 6-2. Rosenbrock 関数の BFGS 法の適用

この図をアニメーションで表示するには、MATLAB コマンド ラインで「bandem」と入力します。

ライン探索法

"ライン探索" は、大規模最適化アルゴリズムの一部として使用される探索法です。メイン アルゴリズムの各ステップで、ライン探索法はメイン アルゴリズムが決めるベクトルである、"探索方向" と平行に、現在の点 xk を含むラインに沿って探索します。すなわちこの方法は、次式の次の反復 xk+1 を見つけます。

xk+1=xk+α*dk,(6-11)

ここで xk は現在の反復を示し、dk は探索方向、α* はスカラー ステップ長のパラメーターです。

ライン探索法は、目的関数の多項式の内挿を反復最小化することにより、ライン xk + α*dk に沿って目的関数を減少させようと試みます。ライン探索法には、2 つのメイン ステップがあります。

  • bracketing 段階は探索されるライン xk+1=xk+α*dk 上の点の範囲を決めます。bracket は、α の値の範囲を指定する区間に相当します。

  • sectioning ステップは bracket を部分区間に分け、部分区間上で目的関数の最小が多項式内挿により近似されます。

結果のステップ長 α は Wolfe 条件を満たします。

f(xk+αdk)f(xk)+c1αfkTdk,(6-12)
f(xk+αdk)Tdkc2fkTdk,(6-13)

ここで、c1 と c2 は 0 < c1 < c2 < 1 を満たす定数です。

第 1 の条件 (式 6-12) は、αk が目的関数を十分減少させることを要求します。第 2 の条件 (式 6-13) は、ステップ長が小さすぎないことを保証します。両方の条件 (式 6-12式 6-13) を満たす点は、"許容点" と呼ばれます。

ライン探索法は [13] の 2-6 節で述べるアルゴリズムの実装です。ライン探索の詳細は [31] を参照してください。

ヘッシアンの更新

最適化関数の多くは、BFGS 法 (式 6-9) を使用して、各反復でヘッセ行列を更新することにより探索方向を決めます。関数 fminunc準ニュートン法 に与えられた DFP 法を使用するオプションも提供します (DFP 法の選択は optionsHessUpdate'dfp' に設定します)。ヘッシアン H は探索方向 d が常に減少する方向になるように正定に保存されます。これは、方向 d で任意のステップ α に対して目的関数の大きさが必ず小さくなることを意味しています。H が正定値であることは H を正定値で初期化し、その後 qkTsk (式 6-14 から導出) が常に正であることにより保証されます。qkTsk の項は、ライン探索ステップ長を表すパラメーター αk と、探索方向 d と過去および現在の勾配評価を組み合わせたものの積として表せます。

qkTsk=αk(f(xk+1)Tdf(xk)Td).(6-14)

qkTsk が正である条件は、十分な精度でライン探索を行うことにより必ず可能になります。これは探索方向 d は減少する方向であるために、αk と負の勾配 –∇f(xk)Td は常に正です。したがって –∇f(xk+1)Td の項が負の値を取る場合は、ライン探索の精度増加と同程度の小ささにしなければなりません。

fminsearch アルゴリズム

fminsearch は、Lagarias らの [57] で説明されている Nelder-Mead のシンプレックス アルゴリズムを使用します。このアルゴリズムは n 次元のベクトル x に対して、n + 1 点のシンプレックスを使用します。このアルゴリズムはまず始めに、x0(i) の要素の 5% を x0 に加え、これらの n ベクトルをx0 に付け加えてシンプレックスの要素として使用し、初期推定 x0 のまわりにシンプレックスを作成します(x0(i) = 0 の場合、i の要素として 0.00025 を使用します)。次にアルゴリズムは以下の手順に従い、シンプレックスを繰り返し変更します。

    メモ:   fminsearch 反復表示のためのキーワードは、手順の説明後に [ ] 付きで表示されます。

  1. x(i) が現在のシンプレックスの点のリスト i = 1,...,n+1 を示すとします。

  2. 最小の関数値 f(x(1)) から最大の関数値 f(x(n+1)) へシンプレックスの点を並べます。反復の各ステップにおいて、アルゴリズムは現在の最悪の点 x(n+1) を破棄し、他の点をシンプレックスにします(あるいは、次の手順 7 の場合に、f(x(1)) より大きい値をもつすべての n 点が変更されます)。

  3. 反射点の生成

    r = 2m – x(n+1),

    ここで、

    m = Σx(i)/n、i = 1...n

    そして、f(r) を計算します。

  4. f(x(1)) ≤ f(r) < f(x(n)) の場合、r を採用し、この反復を終了します。[反射]

  5. f(r) < f(x(1)) の場合、拡張点 s を計算し、

    s = m + 2(m – x(n+1)),

    そして、f(s) を計算します。

    1. f(s) < f(r) の場合、s を採用してこの反復を終了します。[展開]

    2. そうでない場合は r を採用し、この反復を終了します。[反射]

  6. f(r) ≥ f(x(n)) の場合、x(n+1) と r の良い方の値と m の間で縮小を実行します。

    1. f(r) < f(x(n+1)) (すなわち r が x(n+1)) より良い場合は、

      c = m + (r – m)/2

      と、f(c) を計算します。f(c) < f(r) の場合、c を採用してこの反復を終了します。[外側に縮小] そうでない場合は、手順 7 (収縮) を続行します。

    2. f(r) ≥ f(x(n+1)) の場合、

      cc = m + (x(n+1) – m)/2

      そして、f(cc) を計算します。f(cc) < f(x(n+1)) の場合、cc を採用してこの反復を終了します。[内側に縮小] そうでない場合は、手順 7 (収縮) を続行します。

  7. n 点を計算します

    v(i) = x(1) + (x(i) – x(1))/2

    そして、f(v(i)), i = 2,...,n+1 を計算します。次の反復でシンプレックスは x(1), v(2),...,v(n+1) です。[収縮]

以下の図は、fminsearch が新しい各シンプレックスと共にステップで計算するであろう点を示します。元のシンプレックスは太い外側の線で示されます。反復は停止条件を満たすまで行われます。

この情報は役に立ちましたか?