制約条件のない最適化:局所的最小化の方法
はじめに
ほとんどのメソッドの本質は,次のステップを決定するために使われる下の式
で採用される局所的な二次モデルにある.Wolfram言語のFindMinimum関数では本質的に異なる5種類のモデルが選べ,これらはメソッドオプションで制御される.FindMaximumとFindFitでも同様にこれらのメソッドが使われている.
"Newton" | 記号的な導関数が計算できないときには,厳密なヘッセ行列あるいは有限差分近似を用いる |
"QuasiNewton" | 過去のステップによる更新で構築されたヘッセ行列に準ニュートンBFGS近似を用いる |
"LevenbergMarquardt" | 最小2乗問題についてのガウス・ニュートン法である.ヘッセ行列は で近似される.ここで は剰余関数のヤコビ行列である |
"ConjugateGradient" | 線形系を解くための共役勾配法の非線形版である.モデルのヘッセ行列は明示的には形成されない |
"PrincipalAxis" | これは勾配にさえ導関数を使わず,過去のステップからの値を保つことで作用する.各変数にそれぞれ2つの初期条件が必要である |
FindMinimumで選択可能な基本的メソッド
Method->Automaticとすると,問題が構造的に平方和になっていない限り,Wolfram言語は「準ニュートン」法を使う.平方和になっている場合は「ガウス・ニュートン」法のLevenberg–Marquardt法が使われる.各変数に2つの初期条件が与えられると「主軸」法が用いられる.
ニュートン法
Wolfram言語の利点のひとつは,導関数を記号的に計算できる点である.Method->"Newton"と指定し,関数が明示的に微分可能であれば,自動的に記号的な導関数の計算が行われるのである.一方で,関数が明示的に近似できる形式になっていなければ,Wolfram言語は有限差分近似を用いてヘッセ行列を計算する.この際,必要な評価回数を最小にする構造的な情報を用いるか,あるいは変数の数値的値を持つヘッセ行列を返すWolfram言語の式を指定することができる.
この関数は変数の数値的な値で評価するためだけに定義されているので,この関数の導関数を記号的に求めることはできない.
勾配とヘッセ行列の両方が有限差分を使って計算できた場合,ヘッセ行列の誤差が極めて大きい可能性があるので,別のメソッドを使った方がよい.この場合FindMinimumは極めて正確に最小値を見付けるが,導関数の情報が十分ではないので,それも絶対に確かとはいえない.また,関数と勾配の評価数は記号的導関数を自動的に計算した例題におけるよりもはるかに大きくなる.これは,勾配とヘッセ行列のそれぞれを近似するための余分な評価が必要だからである.
勾配(あるいは自動的に計算できる関数)を与えることが可能であるなら,この方法は一般にはるかにうまくいく.勾配はGradientオプションを使って与えることができる.このオプションには「導関数を指定する」ことのできる方法がいくかある.
ヘッセ行列を与えるプログラムがあれば,それも与えるとよい.ヘッセ行列はニュートン法でしか使われないので,"Newton"のオプションとして与える.
原則として,ニュートン法は記号微分を評価するか有限差分を使うかして計算されたヘッセ行列を使う.このようにして計算されたヘッセ行列は常に定正値行列であり,メソッドの収束は凸関数に依存する.この条件が破られた位置から探索が始まることもよくあるので,アルゴリズムはこのことを考慮に入れる必要がある.
ヘッセ行列が定正値行列ではない場合だけにニュートンのステップ式を適用する場合は,関数値の減少を導かないステップ方向にすることができる.
ヘッセ行列が定正値行列ではない場合だけにニュートンのステップ式を適用する場合は,関数値の減少を導かないステップ方向にすることができる.
直線探索法の収束では,正定値の二次形式を使ってその方向を計算することが重要である.そこから導かれる探索過程と収束結果が,十分に降下した方向に依存しているためである.「直線探索法」を参照されたい.Wolfram言語は が定正値行列になれるほど大きい要素を持つ対角行列 によってヘッセ行列を修正する.このようなメソッドは修正ニュートン法と呼ばれることがある.密なヘッセ行列でも疎なヘッセ行列でも, に対する修正は [GMW81]で説明してあるようにコレスキー分解の計算過程のどこかで行われる.この修正はが定正値行列ではないときにのみ行われる.この分解法を独自に使いたければ,LinearSolveからアクセスすることができる.
(修正)ニュートン法は,その強力さの他,収束率も主要な側面といえる.探索が極小値に十分近付くと,収束は二次()収束といわれる.これは が極小値の点であるならば,の定数について以下の式
一般にステップのほとんどは極小値に近付くために使われるので,機械精度では,二次収束により常に重大な相違が生まれるわけではない.しかし,根を非常に高精度で求めたいのであれば,収束が速いという理由から,ニュートン法がたいていの場合最善の選択である.
オプションWorkingPrecision->prec が使われると,AccuracyGoalとPrecisionGoalのデフォルトは prec/2になる.このため,この例題では少なくとも5万桁までの最小値を見付けることになる.
ステップ数よりも関数の評価回数の方が多く必要とされた理由は,Wolfram言語が初期値の精度から,要求された最高のWorkingPrecisionまで適応的に精度を上げるためである.使用される一連の精度は,点が最小値に収束しつつあるという仮定のもとで,最も高価な最終精度での評価回数をできるだけ少なくするように選ばれる.Wolfram言語が精度を変更したときに,より高い精度で関数を再評価しなければならないこともある.
一般的にいって,精度は誤差スケールの約2倍である.ニュートン法では,ステップが計算されるとき,誤差スケールが二次収束に従って実質的に2倍になるので,これが当てはまる.
FindMinimumは常に指定された初期値の精度から始まる.そのため,FindMinimumで適応的な精度制御を使いたくなければ,厳密値,あるいは少なくとも最大のWorkingPrecisionを持つ値から始めることができる.
この方法では関数の評価回数は減るが,すべて最高の精度で行われるので,適応精度を使った方が時間がかなり節約できる.例えば,先程のコマンドで適応精度を使わなかったら,機械精度で始めた場合に比べて50倍以上もの時間がかかる.
ニュートン法には,「直線探索法」と「信頼領域法」の両方の刻み幅制御が実装されている.デフォルトは,上記の例にも用いられている直線探索法であるが,これらの多くは信頼領域法のアプローチで行われる.このアプローチでは一般に1ステップつきより多くの数値線形代数計算が必要になる.しかし,ステップがうまく制御されるので,少ない反復回数で収束する可能性がある.
上記2つのプロットの比較から,信頼領域法では探索が谷に沿って行われ,結果として直線探索法よりも少ない関数評価回数で収束しているため,信頼領域法の方がステップをよく制御していることが分かる.
オプション名
|
デフォルト値
| |
"Hessian" | Automatic | ヘッセ行列式の計算に使う式 |
"StepControl" | "LineSearch" | ステップの制御方法."LineSearch","TrustRegion",Noneのいずれか |
Method->"Newton"のメソッドオプション
準ニュートン法
準ニュートン法にはいろいろな変形がある.そのすべてに共通する概念は,下の二次モデル
の行列 が,関数から構築されたヘッセ行列式の近似,および,以前に取られたステップの一部あるいはすべての勾配値に基づくというものである.
上記の経路でまず目に付くのはこれが間違った方向に向かって始まっている点である.最初のステップではすべての方法がその勾配に従わなければならないために,最急降下のこの方向が選ばれたのである.しかし,後続のステップでは,ヘッセ行列の近似モデルの構築のために取られたステップにおける関数値および勾配値からの情報が組込まれている.
Wolfram言語は,BFGS (Broyden–Fletcher–Goldfarb–Shanno)更新法を使い,大規模システムでは,メモリ節約型のL-BFGS法を使う.L-BFGS法ではモデル は明示的には保存されず,過去のステップから保存された勾配とステップ方向によって が計算される.
BFGSメソッドは,各ステップでヘッセ行列の近似 を形成する代りに,回の操作で 変数の問題の系 [DS96]が解けるように, を計算するコレスキー因子 を実装している.
大規模で疎な問題の場合,BFGS法には少々問題がある.一般にコレスキー因子(あるいはヘッセ行列の近似 またはその逆行列)は密なので,のメモリと操作の必要条件が,疎性を利用するアルゴリズムに比べて抑制的になるからである.L-BFGSアルゴリズム[NW99]は,最後の ステップに基づいてヘッセ行列の逆行列の近似行列を形成し,保存する.ヘッセ行列の近似はそれほど完全ではないかもしれないが, 変数の問題の操作メモリと順序は に限られている.Wolfram言語では,250以上の変数を持つ問題になると自動的にアルゴリズムがL-BFGSに切り換えられる.これは,メソッドオプションの"StepMemory"->m を使って制御することができる.のときは,完全なBFGSメソッドが常に使われる. の適切な値を選択することで,収束速度と各ステップでの作業量との間にトレードオフの関係が生じる.の場合は,「共役勾配」アルゴリズムを使った方が結果がよくなることが多い.
準ニュートン法は,一般に極めて速く,記号計算と数値評価の両方でかなり高くつくヘッセ行列の計算を必要としないので,Wolfram言語はデフォルトでこのメソッドを使う.準ニュートン法で適切な「直線探索」を行うと,ヘッセ行列が正定値行列となる極小値に超直線的に収束する[NW99] .つまり,以下の式
が成り立つ.言い換えれば,ステップがどんどん小さくなっている.しかし,非常に高精度の場合は,「ニュートン」法の -二次収束率とは比較にならない.
ニュートン法は二次収束率が高いので,同じ桁数をより少ないステップで10倍速く見付けることができる.しかし,準ニュートン法は誤差の割合が明かにゼロに近付いているために,収束がなおも超線形的なのが分かる.
オプション名
|
デフォルト値
| |
"StepMemory" | Automatic | ヘッセ行列の近似で「記憶」される有効なステップ数.正の整数あるいはAutomatic |
"StepControl" | "LineSearch" | 刻み幅制御法."LineSearch"かNone |
Method->"QuasiNewton"のメソッドオプション
ガウス・ニュートン法
では,問題の特殊構造を使うとよい場合がよくある.剰余関数 とその導関数であるヤコビ行列 を計算することで,時間と努力を節約することができる.ガウス・ニュートン法はこれを行う洗練された方法である.ガウス・ニュートン法は,二次形式として完全な二階のヘッセ行列を使う代りに,ステップ が以下の公式
から計算されるように,式(1)で を使う.ここで 等となる.これは完全なヘッセ行列 の近似であることに注意のこと.の点が最小値となったり, が最小点近くの線形方程式のように変化したりする場合,剰余がゼロというが,この場合ヘッセ行列の近似はきわめてよく,「ニュートン法」の二次収束がよく観察される.
目的関数は平方和であるが,極めて一般的である.実際に,FindFitがNormFunctionオプションのデフォルト値で用いられるときは,目的関数の形式をしている.ガウス・ニュートン法は,最小二乗の問題の点から見ることもできる.ガウス・ニュートンステップを解くことは線形最小二乗問題を解くことと同じなので,ガウス・ニュートン法を適用することは,実質的に非線形関数に線形の最小二乗フィットを適用することである.このように見てみると,このメソッドがFindFitが行う非線形フィットのようなものに特に適しているのも理解できる.
Wolfram言語が,Rosenbrockの例のような明示的に平方和である問題,あるいはベクトルとそのべクトル自身の内積であるような関数に出会うと,自動的にガウス・ニュートン法が使われる.
同じ例を「ニュートン法」で行ったものと比べると,このガウス・ニュートン法の方がステップ数も評価回数も少ないのが分かる.これは,ガウス ・ニュートン法が問題の特殊構造を利用しているからである.最小値で剰余がゼロとなるので,最小値付近の収束率はニュートン法と同じくらいよい.
Levenberg–Marquardt法は「信頼領域」刻み幅制御を使ったガウス・ニュートン法である(この方法が初めて提唱されたのは,一般的な信頼領域概念が開発される前である).このメソッドはFindMinimumオプションのMethod->"LevenbergMarquardt"あるいはMethod->"GaussNewton"によって使うことができる.
時に関数を明示的に平方和,あるいはベクトルとそのべクトル自身の内積になるように表現するのが厄介なことがある.そのような場合には"Residual"メソッドオプションを使って剰余を直接指定することもできる.同様に,剰余の導関数を"Jacobian"メソッドオプションで指定することができる."Residual"メソッドオプションで剰余を指定した場合は,FindMinimumの第1引数との整合性のチェックは行われない.返された値はオプションで与えられた値に依存する.
オプション名
|
デフォルト値
| |
"Residual" | Automatic | となるように剰余 を直接指定することを許可する |
"EvaluationMonitor" | Automatic | 剰余が評価されるたびに評価される式 |
"Jacobian" | Automatic | 剰余の(行列)導関数の指定を許可する |
"StepControl" | "TrustRegion" | "TrustRegion"でなければならないが,メソッドオプションを使ってパラメータの制御を変更することを許可する |
Method->"LevenbergMarquardt"のメソッドオプション
Wolfram言語で自然に平方和の問題を設定するための別の方法に,FindFitを使うものがある.
この例では,FindFitが内部的に剰余関数とヤコビ行列を形成している.これらはガウス・ニュートン法で平方和,あるいは非線形最小二乗フィットを求めるのに使われる.無論,FindFitは他のメソッドで使うこともできるが,速く評価できる剰余関数が構築できるので,通常は他のメソッドを使うよりも速い.
非線形共役勾配法
非線形共役勾配法の基礎は,実質的に線形共役勾配法を適用し,剰余を勾配で置き換えるものである.モデル二次関数は決して明示的には形成されないので,常に「直線探索法」と一緒に使われる.
最初の非線形共役勾配法はFletcherとReevesによって提唱された次のようなものである.ステップの方向 が与えられると,直線探索を使って が成り立つ を求めて,次に以下
を計算する. を求めるための直線探索は,強いWolfe条件を満たさなければならない.このことは,方向 が確実に降下方向 [NW99]となるために必要である.
これに代る方法で,実際には大抵の場合で(常にではないが)これよりうまくいく方法は,PolaktとRibiereによるものである.ここでは式(2)が以下で置き換えられる.
式(3)では が負になることも可能である.そのような場合には,Wolfram言語が を使って変更されたアルゴリズムを使用する.Wolfram言語ではデフォルトの共役勾配法はPolak–Ribiereであるが,メソッドオプション
Method->{"ConjugateGradient",Method->"FletcherReeves"}.
を使ってFletcher–Reevesメソッドを選ぶこともできる.共役勾配法の利点は,大規模問題でも比較的少量のメモリしか使わず,数値線形代数は必要としないので,各ステップが非常に速いという点である.不利な点は一般に「ニュートン法」や「準ニュートン法」に比べて収束に時間がかるという点である.また,長さに対するステップのスケールが貧弱なので,「直線探索」アルゴリズムで許容し得るステップを求めようとするたびに必要な反復回数が多くなることがある.
非線形共役勾配法における問題点のひとつは,いつ再スタートするかということである.探索が動くにつれて,関数の局所二次近似の性質が根本的に変化する可能性がある.このメソッドの局所的収束は,二次関数が一定となる線形共役勾配法の収束に依存している.n 変数の一定の二次関数と厳密な直線探索で,線形アルゴリズムが n 回以下の反復回数で収束する.再スタートする(で最急降下ステップを取る)ことで,現行の探索点における局所二次モデルには無関係な前のステップからの情報を除去できることがしばしばある.例題を注意深く見ると,メソッドが再スタートされ,最急降下ステップが取られた点が分かる.ひとつのオプションとして k 回(ここで k<=n)の反復毎に再スタートするというものがある.これはメソッドオプションの"RestartIterations"->k で指定することができる.あるいは,下記のテスト
の結果,連続する勾配が十分に直角ではない場合に再スタートすることも可能である.閾値 ν は0から1までの間である.これはメソッドオプションの"RestartThreshold"->ν を使って指定することができる.
次の表は,共役勾配法で使うことができるオプションのまとめである.
オプション名
|
デフォルト値
| |
"Method" | "PolakRibiere" | 非線形共役勾配法は"PolakRibiere"あるいは"FletcherReeves"である |
"RestartThreshold" | 1/10 | それ以下では再スタートする勾配直角の閾値 ν |
"RestartIterations" | ∞ | 次の再スタートまでの反復回数 |
"StepControl" | "LineSearch" | "LineSearch"でなければならないが、これを使って直線探索法を指定することができる |
Method->"ConjugateGradient"のメソッドオプション
バージョン4におけるFindMinimumのデフォルトオプションは,ほぼ厳密な直線探索を使う共役勾配法であった点に注意のこと.これはレガシーとして残っており,FindMinimumオプションのMethod->"Gradient"でアクセスすることができる.一般にこちらの方が新しいMethod->"ConjugateGradient"よりも関数や勾配の評価回数が多い.この新しい方でも,Wolfram言語が現在デフォルトで使うメソッドよりもはるかに多くのメソッドを使っている.
主軸法
「ガウス・ニュートン法」と「共役勾配法」には導関数が用いられる.Wolfram言語が記号的な導関数を計算できないときは,有限差分が使われる.有限差分を用いて導関数を計算するのは,場合によっては非常に大きなコストがかかり,導関数の妥当性にも影響し,究極的には最小値へどれほどまで近似できるかに影響する.記号的な導関数が使えない関数については,近似モデルが関数の評価による値のみを使って構築される,導関数を使わないアルゴリズムを使うこともできる.
Wolfram言語では,導関数を使わないアルゴリズムとしてブレント[Br02]の主軸法を使う. 変数の問題では,一連の探索方向 と点 を取る. が から 方向に沿って を最小化する点だとする(つまり,から「直線探索」を行う).次に を で置換する.最後に, を で置換する.理想的には,新たな反復ができるように,新たな は線形に独立でなければならないが,実際にはそうはならない.ブレントのアルゴリズムでは,局所的な二次形式の主方向に を再調整するのに,行列 の特異値分解が使われる(固有分解を行ってもよいが,ブラントは特異値分解の方が有効であることを示している).新たに獲得した の集合を使って,もう一度反復を行うことができる.
ベクトル の大きさを決定するのに使うため,このメソッドでは,それぞれの変数に2つの異なる初期条件が必要である.事実,各変数について2つの初期条件を指定する場合はいつも,FindMinimum,FindMaximum,FindFitはデフォルトで主軸法を使う.
導関数を使わない直線探索アルゴリズムは極めて多数の関数評価を必要とするため,この探索アルゴリズムの基本はこのプロットからよく分かる.まず, 方向に直線探索が行われる.次に,その点から 方向への直線探索が行われ,ステップの方向が決定される.一旦ステップが取られると,ベクトル は局所的な二次形式近似の主方向に適切に再調整され,次のステップが同じように計算される.
このアルゴリズムは収束率という点で効果的である.これはステップについて二次収束するからである.しかし,関数評価の点から見ると,導関数を使わない直線探索が必要となるため極めて高くなっている.直線探索に与えられる方向(特に最初で)は,必ずしも降下方向である必要はないので,直線探索は両方向の探索が行えなければならないという点に注意のこと.多変数の問題では,全方向での別個の直線探索は非常に高くつく.そのため,この方法は変数が多過ぎない問題に適している.