母集団最適化アルゴリズム:ネルダー–ミード法、またはシンプレックス(NM)検索法
内容
1.はじめに
ネルダー-ミード法は、1965年にジョン・ネルダーとロジャー・ミードによって開発されました。彼らは、導関数を持たない関数や、導関数の解析方程式を持たない関数に対応できる最適化手法を探していました。また、実装が簡単で、当時の計算機で効率的に使用できる方法を開発したいと考えていました。研究の結果、彼らはシンプレックス(関数パラメータ空間の多面体)を使うというアイデアにたどり着きました。
この方法が生まれた歴史は、オックスフォードにあるコンピューティング研究所でのネルダーとその同僚たちの研究から始まりました。彼らは、解析的導関数を持たない、あるいは計算が複雑すぎる関数の最適化という問題に直面していました。このような場合、勾配法などの従来の最適化手法は適用できませんでした。その代わりに、ネルダーとミードは、関数パラメータの空間で最適解を繰り返し探索することに基づいた新しい方法を提案しました。
ネルダー-ミード法は「シンプレックス法」と呼ばれ、1965年に『The Computer Journal』誌に「A Simplex Method for Function Minimization」という論文で発表されました。この方法は科学界に受け入れられ、関数の最適化を必要とする様々な分野で広く使用されるようになりました。
シンプレックスとは、多面体を形成する点の集合のことで、各点は最適化される関数のパラメータ値の集合です。そのアイデアは、関数の最適値を見つけるために、パラメータ空間内でシンプレックスを変更したり移動させたりすることです。
ネルダー-ミード法(ネルダー-ミードシンプレックス法)は、無条件最適化アルゴリズムのクラスに属します。ネルダー-ミード法は決定論的アルゴリズムであり、関数の導関数を使用する必要がなく、複数の局所最小値を持つ関数を扱うことができます。
2.アルゴリズム
ネルダー-ミード法は母集団法ではありません。というのも、使用される最適化エージェントは1つだけで、シンプレックスを表し、そのシンプレックスは探索空間内の点の集合で構成され、それ自体が母集団と考えることができるからです。ただし、私たちは複数のエージェントを使用し、それぞれが独自のシンプレックスを持つので、この実装は集団アルゴリズムに分類されるかもしれません。
そこで、理解を容易にするために、2つの変数関数を最適化する必要があると仮定します。言い換えれば、空間の次元がz=2であれば、シンプレックスの頂点の数はz+1=3となります。シンプレックスの変更は、これら3つのうち最悪の点に対する操作です。これらの点をBest、Good、Worstとし、GoodはWorstの次にくる2番目の点とします(次元が2より大きい場合、Goodは常に頂点の並び替えリストからWorstに続いた2番目の点になる)。
次に、この3つの頂点を使用して、残りの頂点との相対的なWorstを計算し、調整する必要があります。それ自体を除く頂点の幾何学的中心に対してWorstを移動します。つまり、重心を計算し、その重心に対してWorstを移動します。
ネルダー-ミード法の各反復では、以下のステップを実行します。
1.目的関数の値に従ってシンプレックスの頂点を並び替えます。
f(P0) ⩾ f(P1) ⩾ ... ⩾ f(P(n-1))
ここで
- P0, P1 ...P(n-1) - シンプレックスの点
- n - シンプレックスの頂点の数
2.セントロイド計算:最悪の頂点を除く、シンプレックスのすべての頂点の対応する座標の平均値を計算します。XoベクトルをX(n-1)以外のすべての頂点の重心とします。
Xo = (X0 + X1 + ... + X(n-2)) / (n-1)
ここで
- X0, X1 ...X(n-2) - 最悪の頂点を除くすべての頂点の座標の対応するベクトル
3.反射:セントロイドに対して最も悪い頂点を反射させます。新しい反射頂点のXrベクトルは次のように計算されます。
Xr = Xo + α (Xo - Xw)
ここで
- Xw - 最悪頂点ベクトル,頂点P(n-1)に対応
- α - 反射率、通常は1に等しい
4.新しい頂点の評価:目的関数の値がこの時点で計算されます。新しい頂点がシンプレックスの最悪の頂点よりも優れた目的関数値を持つならば、その頂点はその頂点と置き換わることができます。反射操作は、シンプレックスの重心から反射することによって、パラメータ空間の新しい領域を探索することを可能にします。
5.膨張:反射操作によって最良の頂点よりも優れた頂点が得られた場合に、反射方向にさらに移動することでその位置がさらに改善されるという仮定のもと、パラメータ空間をさらに探索するために使用されます。膨張操作では、反射操作よりもさらに大きなパラメータ空間の領域を探索することができます。これは重心と反射点間の距離を増加させ、目的関数の新たな局所最小値を見つけるのに役立ちます。Xeの膨張ベクトルは以下のように計算できます。
Xe=Xo+γ(Xr-Xo)
ここで
- γ -膨張率。通常は1より大きく、デフォルトは2
5.1.膨張操作がおこなわれた場合、Xe頂点の適応度関数が計算されます。Xr頂点より優れていれば、最悪のXwシンプレックス頂点と置き換えることができます。膨張を終えたら、ステップ1に進みます。
γ膨張比の選択は慎重におこなってください。この値が大きすぎると、シンプレックスがパラメータ空間の広い領域に広がってしまい、収束が遅くなったり、極値に関する情報が失われたりする可能性があります。小さすぎると、パラメータ空間の新しい領域を十分に探索できない可能性があります。
6.収縮:反射操作によってXbの良い点(悪い点の次の点)以下の頂点が生成された場合、Xcの頂点を計算するために使用されます。式は以下の通りです。
Xc=Xo+β(Xw-Xo)
ここで
- β-収縮率、通常は0~1、デフォルトは0.5
6.1. 収縮がおこなわれた場合:その頂点のXc目的関数の値が計算されます。もし新しい頂点が、シンプレックスの最悪の頂点よりも優れた目的関数値を持つならば、その頂点は最悪の頂点を置き換えることができます。膨張を終えたら、ステップ1に進みます。
収縮により、単純点を最悪点に近づけることができ、この点の近辺を探るのに役立ちます。膨張と同様、β収縮率を慎重に選ぶことが重要です。この値が大きすぎると、シンプレックスが圧縮されすぎて、極小に関する情報が失われてしまいます。小さすぎると、最悪点の近辺を探るには圧縮が十分でない可能性があります。
7.全収縮:ネルダー-ミード法の最後の操作は、前の操作(反射、拡大、収縮)のどれもが目的関数値の改善につながらない場合に実行されます。全収縮はシンプレックスのサイズを縮小し、検索領域を狭めてスイートスポット周辺に焦点を絞ります。
見てわかるように、著者らはシンプレックスを使用して4つの操作をおこなっています。最適化を開始するには、すべてのシンプレックス点について適応度関数を計算する必要があります。シンプレックス点の数はz+1に等しく、ここでzは探索空間の次元です。例えば、探索次元が1000であれば、適合度関数を1001回評価し、シンプレックス演算を開始しなければなりません。50エージェントの母集団を持つ典型的な母集団アルゴリズムは20エポックを実行するかもしれませんが、このアルゴリズムは1エポックしか実行せず、限られた反復回数のほとんどを消費することになります。
全収縮では、すべての点を最適なものにシフトします。その後、適応度関数を1000回計算する必要があります。つまり、全収縮はリソースを多量に消費します。著者らの考えによれば、この操作は、アルゴリズムが行き詰まり、シンプレックスの最悪点を移動しても解が改善されない場合に呼び出されるべきです。しかし、私の実験によると、この操作はコンピューティングリソースの点で非常に高価であり、さらに、シンプレックスエージェントのすべての点が1つの極値に収束するということは、行き詰まって探索空間の探索が停止することを意味するため、最適化アルゴリズムにとっては無意味で自殺行為であることがわかりました。そのため、アルゴリズムの実装にあたっては、全収縮を断念することにしました。
ここでは、シンプレックスを使用して最初の3つの操作をおこないます。より分かりやすくするため、図1に表示しました。
図1:ネルダー-ミード法の3つの主要操作
シンプレックス探索法は決定論的であるため、初期頂点の選択に失敗する可能性があり、すればその最適化結果は満足のいくものではないかもしれません。このアルゴリズムを使った実験では、限られた座標空間でしか機能しないという懸念が確認されただけでした。
収束の問題があるため、シンプレックスの頂点は単純に固定値だけシフトされます。竹馬の上に立ってベンチに座ろうとしているところを想像してみてください。ベンチが低すぎます。変化するシンプレックスでもまったく同じ状況が生じます。アルゴリズムが完璧に機能するためには、シンプレックスの初期頂点が正しい位置に来るよう、運を味方につける必要があります。これは竹馬に似ています。安全に座るためには、高いベンチが必要なのです。このアルゴリズムは、滑らかな単極端関数、例えば放物線上では比較的よく収束します。ただし、現実的な問題はもっと複雑で、特に離散的な性質を持つ取引問題では、局所的な「罠」に満ちているのが普通です。
そこで疑問が生じます。シンプレックスが極値で「永遠に」立ち往生してしまった場合、どうすればいいのでしょうか。ネルダー-ミードアルゴリズムの論理は、この「罠」から抜け出す方法を提供していません。収縮であれ膨張であれ、各操作の後、アルゴリズムは最悪点を克服しようと反射に戻ります。視覚的には「点滅」のように見えます。この問題を解決するために、私は、シンプレックスに、局所的な「罠」を離れて新たな空間で探索を続ける能力を与えようとしました。これを実現するために、以前の記事で紹介したように、場合によっては「復活」に役立ったレヴィフライト を使用しました。ただし、レヴィフライトが常に有用であるとは限らず、その適用は最適化アルゴリズムのロジックに依存することは注目に値します。私たちの場合、これは実験であり、改善結果は保証されていませんでした。
さて、最も興味深い部分に移りましょう。コードを書きます。
各反復の後、シンプレックスの最悪の頂点に対してどの操作がおこなわれたかを知る必要があります。この作業には列挙子が便利です。何らかの方法でアルゴリズムを補完することにした場合に備えて、シンプレックス頂点の3つの主要な操作に加えて、もう1つ「なし」を追加しました。列挙子をE_SimplexOperationと呼ぶことにしましょう。次のフィールドがあります。
- none:操作なし
- reflection
- expansion
- contraction
//—————————————————————————————————————————————————————————————————————————————— enum E_SimplexOperation { none = 0, reflection = 1, expansion = 2, contraction = 3 }; //——————————————————————————————————————————————————————————————————————————————シンプレックス頂点を記述するには、以下のフィールドを含むS_Point構造体が必要です。
- Init(intcoords):点の初期化。coords引数を取る
- c[]:多次元空間における点の座標の配列。配列のサイズはInitメソッドで決定
- f:頂点適応度関数値
//—————————————————————————————————————————————————————————————————————————————— struct S_Point { void Init (int coords) { ArrayResize (c, coords); f = -DBL_MAX; } double c []; double f; }; //——————————————————————————————————————————————————————————————————————————————
シンプレックスの場合、S_Simplex 構造体を使用します。これは各エージェントの単純なシンプレックスであり、次のフィールドが含まれます。
- Init(intcoords):シンプレックスの初期化メソッドで、座標の数を指定するcoords引数を取る
- S_Pointp[]:シンプレックス頂点の配列。配列サイズはInitメソッドで決定
- c[]:単純重心座標の配列。配列サイズはInitメソッドで決定
- S_Point Xr:反射点
- S_Point Xe:膨張点
- S_Point Xc:収縮点
- indG:シンプレックスのGoodポイントインデックス
- indW:シンプレックスのWorstポイントインデックス
- E_SimplexOperation:シンプレックスの操作タイプ
//—————————————————————————————————————————————————————————————————————————————— struct S_Simplex { void Init (int coords) { ArrayResize (p, coords + 1); for (int i = 0; i <= coords; i++) { p [i].Init (coords); } ArrayResize (c, coords); Xr.Init (coords); Xe.Init (coords); Xc.Init (coords); operation = none; } S_Point p []; double c []; //centroid S_Point Xr; //reflection point S_Point Xe; //expansion point S_Point Xc; //expansion point int indG; //index of good point int indW; //index of worst point E_SimplexOperation operation; //type of simplex operation }; //——————————————————————————————————————————————————————————————————————————————
最適化エージェントのS_エージェント構造体を定義します。この構造体には以下のフィールドが含まれます。
- Init (int coords):エージェントの初期化メソッドで、座標の数を指定するcoords引数を取る。ArrayResize関数を使用して配列cのサイズをcordsに変更する。fフィールドは、エージェント評価関数の初期値である-DBL_MAXに設定される。次に、エージェントのシンプレックスを表すsフィールドに対してInitメソッドが呼ばれる。
- c[]:多次元空間におけるエージェント座標の配列。配列サイズはInitメソッドで決定
- f:エージェントの適応度関数値
- s:S_Simplex構造体で表現されるエージェントシンプレックス。sフィールドのInitメソッドを呼び出すことで初期化される
//—————————————————————————————————————————————————————————————————————————————— struct S_Agent { void Init (int coords) { ArrayResize (c, coords); f = -DBL_MAX; s.Init (coords); } double c []; //coordinates double f; //fitness S_Simplex s; //agent simplex }; //——————————————————————————————————————————————————————————————————————————————
シンプレックス検索メソッドのアルゴリズムクラスC_AO_NMmを定義します。これには次のフィールドが含まれます。
- cB[]:最適座標の配列
- fB:最適座標の適応度関数値
- a[]:S_Agent構造体で表されるエージェントの配列
- rangeMax[]:各座標の検索範囲の最大値の配列
- rangeMin[]:各座標の検索範囲の最小値の配列
- rangeStep[]:各座標の検索ステップの配列
- Init (const int coordsP, const int popSizeP, const double reflectionCoeffP, const double expansionCoeffP, const double contractionCoeffP):アルゴリズムの初期化メソッド。座標の数,母集団のサイズ,反射,膨張,収縮操作の比率を定義する引数を取り、すべてのクラスフィールドの初期化をおこなう
- Moving():アルゴリズム内でエージェントの移動を実行するメソッド
- Revision():アルゴリズム内でエージェントの移動を実行するメソッド
- Sorting、CalcCentroid、Reflection、Expansion、Contraction、Flying:シンプレックスを使った操作を実行する。
- SeInDiSpメソッド、RNDfromCIメソッド、Scaleメソッドは、値を生成し、与えられたステップで有効な範囲に変換するために実行される
//—————————————————————————————————————————————————————————————————————————————— class C_AO_NMm { //---------------------------------------------------------------------------- public: double cB []; //best coordinates public: double fB; //FF of the best coordinates public: S_Agent a []; //agent public: double rangeMax []; //maximum search range public: double rangeMin []; //manimum search range public: double rangeStep []; //step search public: void Init (const int coordsP, //coordinates number const int popSizeP, //population size const double reflectionCoeffP, //Reflection coefficient const double expansionCoeffP, //Expansion coefficient const double contractionCoeffP); //Contraction coefficient public: void Moving (); public: void Revision (); //---------------------------------------------------------------------------- private: int coords; //coordinates number private: int popSize; //population size private: int simplexPoints; //simplex points private: double reflectionCoeff; //Reflection coefficient private: double expansionCoeff; //Expansion coefficient private: double contractionCoeff; //Contraction coefficient private: bool revision; private: S_Point pTemp []; private: int ind []; private: double val []; private: void Sorting (S_Point &p []); private: void CalcCentroid (S_Simplex &s, int indW); private: void Reflection (S_Agent &agent, int indW); private: void Expansion (S_Agent &agent); private: void Contraction (S_Agent &agent, int indW); private: void Flying (S_Agent &agent, int indW); private: double SeInDiSp (double In, double InMin, double InMax, double Step); private: double RNDfromCI (double min, double max); private: double Scale (double In, double InMIN, double InMAX, double OutMIN, double OutMAX, bool revers); }; //——————————————————————————————————————————————————————————————————————————————
アルゴリズムを初期化するには、C_AO_NMmクラスのInitメソッドを使用します。このメソッドは、クラスのすべてのフィールドを初期化し、最適化アルゴリズムが動作するために必要な配列を作成します。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Init (const int coordsP, //coordinates number const int popSizeP, //population size const double reflectionCoeffP, //reflection coefficient const double expansionCoeffP, //expansion coefficient const double contractionCoeffP) //contraction coefficient { MathSrand ((int)GetMicrosecondCount ()); // reset of the generator fB = -DBL_MAX; revision = false; coords = coordsP; popSize = popSizeP; simplexPoints = coords + 1; reflectionCoeff = reflectionCoeffP; expansionCoeff = expansionCoeffP; contractionCoeff = contractionCoeffP; ArrayResize (pTemp, simplexPoints); ArrayResize (ind, simplexPoints); ArrayResize (val, simplexPoints); ArrayResize (a, popSize); for (int i = 0; i < popSize; i++) { a [i].Init (coords); } ArrayResize (rangeMax, coords); ArrayResize (rangeMin, coords); ArrayResize (rangeStep, coords); ArrayResize (cB, coords); } //——————————————————————————————————————————————————————————————————————————————
通常、探索空間内のエージェントの移動にはMovingメソッドを用いますが、NMアルゴリズムでは、探索空間内のエージェントシンプレックス頂点の初期ランダム配置の機能のみを残します。
シンプレックスのランダムな頂点を作成するために、RNDfromCIを使用し、SeInDiSpメソッドを使用して必要なステップで必要な範囲に持っていきます。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Moving () { if (!revision) { int cnt = 0; for (int i = 0; i < popSize; i++) { for (int p = 0; p < simplexPoints; p++) { cnt++; for (int c = 0; c < coords; c++) { a [i].s.p [p].c [c] = RNDfromCI (rangeMin [c], rangeMax [c]); a [i].s.p [p].c [c] = SeInDiSp (a [i].s.p [p].c [c], rangeMin [c], rangeMax [c], rangeStep [c]); } } } } } //——————————————————————————————————————————————————————————————————————————————
エージェントのシンプレックス頂点を移動させる主要なロジックは、Revisionメソッドで実行されます。
revisionがまだ実行されていなければ(if(!revision))、シンプレックスの頂点を並び替え、セントロイドを構築し、エージェントの最悪のシンプレックス頂点に対して反射を実行する必要があります。これは頂点に対する最初の操作で、常に「反射」です。この時点で各頂点の適応度がわかるので、ここで大域解を更新することができます。変数revisionの値をtrueに設定し、メソッドを終了します。
次の反復からは、シンプレックス頂点の適応度関数ではなく、特定の最適化エージェントの適応度関数を知ることになります。シンプレックス頂点の操作に従って、エージェントの既知の解を対応するシンプレックス頂点の頂点に割り当てます。エージェントの最適な適応度で大域解を更新する必要があります。
次に、前回の反復で頂点に対してどのような操作がおこなわれたかを確認します。これはシンプレックス頂点がベクトルであることを意味します。
「反射」がおこなわれた場合:
- エージェントの適応度をXrベクトルに割り当てる
- Xrベクトルの適応度をXwと比較し、その値の方が良ければ、Xwベクトルを更新する
- Xrベクトルの適応度をXbと比較し、その値の方が良ければ、「膨張」を実行する
- Xrベクトルの適応度をXg(最悪頂点から2番目)と比較し、その値の方が悪ければ「圧縮」を実行する
- 頂点が改善された場合は、事前の並び替えとセントロイドの構築を伴う「反射」を実行し、そうでなければ、Flyingを実行する
「膨張」をおこなった場合:
- エージェントの適応度をXeベクトルに割り当てる
- Xeベクトルの適応度をXwと比較し、その値の方が良ければ、Xwベクトルを更新する
- 予備的な並び替えとセントロイドの構築を伴って「反射」を実行する
「収縮」をおこなった場合:
- エージェントの適応度をXсベクトルに割り当てる
- Xсベクトルの適応度をXwと比較し、その値の方が良ければ、Xwベクトルを更新し、そうでなければ、Flyingを実行する
- 頂点が改善された場合は、事前の並び替えとセントロイドの構築を伴う「反射」を実行する
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Revision () { //---------------------------------------------------------------------------- if (!revision) { //sort agent simplex points by FF value and assign BGW for (int i = 0; i < popSize; i++) { Sorting (a [i].s.p); } //calculate the simplex centroid for (int i = 0; i < popSize; i++) { a [i].s.indG = simplexPoints - 2; a [i].s.indW = simplexPoints - 1; CalcCentroid (a [i].s, a [i].s.indW); } //assign the next type of operation - reflection for (int i = 0; i < popSize; i++) { Reflection (a [i], a [i].s.indW); a [i].s.operation = reflection; } //save the best point of the agents’ simplexes as a global solution for (int i = 0; i < popSize; i++) { if (a [i].s.p [0].f > fB) { fB = a [i].s.p [0].f; ArrayCopy (cB, a [i].s.p [0].c, 0, 0, WHOLE_ARRAY); } } revision = true; return; } //---------------------------------------------------------------------------- if (revision) { int pos = -1; for (int i = 0; i < popSize; i++) { if (a [i].f > fB) pos = i; } if (pos != -1) { fB = a [pos].f; ArrayCopy (cB, a [pos].c, 0, 0, WHOLE_ARRAY); } } //---------------------------------------------------------------------------- for (int i = 0; i < popSize; i++) { //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //if there was reflection ++++++++++++++++++++++++++++++++++++++++++++ if (a [i].s.operation == reflection) { a [i].s.Xr.f = a [i].f; bool needUpd = false; //------------------------------------------------------------------------ if (a [i].s.Xr.f > a [i].s.p [a [i].s.indW].f) //Xr > Xw |---w--Xr--g----------b---| { a [i].s.p [a [i].s.indW] = a [i].s.Xr; //replace Xw with Xr needUpd = true; } //------------------------------------------------------------------------ if (a [i].s.Xr.f > a [i].s.p [0].f) //if Xr > Xb |---w----g----------b---Xr| { Expansion (a [i]); //perform expansion continue; } //------------------------------------------------------------------------ if (a [i].s.Xr.f <= a [i].s.p [a [i].s.indG].f) //if Xr < Xg |---w---Xr--g----------b--| { Contraction (a [i], a [i].s.indG); //perform contraction continue; } if (needUpd) { Sorting (a [i].s.p); a [i].s.indG = simplexPoints - 2; a [i].s.indW = simplexPoints - 1; CalcCentroid (a [i].s, a [i].s.indW); Reflection (a [i], a [i].s.indW); } else Flying (a [i], a [i].s.indW); continue; } //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //if there was expansion +++++++++++++++++++++++++++++++++++++++++++++ if (a [i].s.operation == expansion) { a [i].s.Xe.f = a [i].f; if (a [i].s.Xe.f > a [i].s.Xr.f) a [i].s.p [a [i].s.indW] = a [i].s.Xe; else a [i].s.p [a [i].s.indW] = a [i].s.Xr; Sorting (a [i].s.p); a [i].s.indG = simplexPoints - 2; a [i].s.indW = simplexPoints - 1; CalcCentroid (a [i].s, a [i].s.indW); Reflection (a [i], a [i].s.indW); continue; } //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ //if there was contraction +++++++++++++++++++++++++++++++++++++++++++ if (a [i].s.operation == contraction) { a [i].s.Xc.f = a [i].f; if (a [i].s.Xc.f > a [i].s.p [a [i].s.indW].f) { a [i].s.p [a [i].s.indW] = a [i].s.Xc; Sorting (a [i].s.p); a [i].s.indG = simplexPoints - 2; a [i].s.indW = simplexPoints - 1; CalcCentroid (a [i].s, a [i].s.indW); Reflection (a [i], a [i].s.indW); } else Flying (a [i], a [i].s.indW); continue; } } } //——————————————————————————————————————————————————————————————————————————————
セントロイドを計算するために、CalcCentroidメソッドを書きます。これは、指定されたシンプレックスs内の座標の指定されたindWインデックスまでの平均値を計算します。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::CalcCentroid (S_Simplex &s, int indW) { double summ = 0.0; for (int c = 0; c < coords; c++) { summ = 0.0; for (int p = 0; p < simplexPoints; p++) { if (p != indW) summ += s.p [p].c [c]; } s.c [c] = summ / coords; } } //——————————————————————————————————————————————————————————————————————————————
「反射」操作は、「エージェント」シンプレックス頂点のindWインデックス頂点に対してReflectionメソッドによって実行されます。
ループ内では、Xr = Xo + reflectionCoeff * (Xo - Xw)という式を使用して反射値Xrを計算します。ここでreflectionCoeffは反射率であり、反射された頂点の初期頂点からのずれの度合いを決定します。
次に、SeInDiSpをXrに適用して、有効なrangeMin[c]とrangeMax[c]の範囲内にあることを、rangeStep[c]のステップで確認します。結果はagent.s.Xr.c[c]に保存されます。
agent.s.operationエージェントにreflection操作を設定することは、このエージェントに対してreflection操作が実行されたことを示します。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Reflection (S_Agent &agent, int indW) { double Xo; double Xr; double Xw; for (int c = 0; c < coords; c++) { Xo = agent.s.c [c]; Xw = agent.s.p [indW].c [c]; //Xr = Xo + RNDfromCI (0.0, reflectionCoeff) * (Xo - Xw); Xr = Xo + reflectionCoeff * (Xo - Xw); agent.s.Xr.c [c] = SeInDiSp (Xr, rangeMin [c], rangeMax [c], rangeStep [c]); agent.c [c] = agent.s.Xr.c [c]; } agent.s.operation = reflection; } //——————————————————————————————————————————————————————————————————————————————
「膨張」操作は、「エージェント」にたいしてExpansionメソッドによって実行されます。
ループの中で、Xe = Xo + expansionCoeff * (Xr - Xo)という式を使用してXeの膨張値を計算します。ここで、expansionCoeffは、拡大された頂点の初期頂点からのずれの増加の度合いを決定する拡大率です。
次に、SeInDiSpをXeに適用して、有効なrangeMin[c]とrangeMax[c]の範囲内にあることを、rangeStep[c]のステップで確認します。結果はagent.s.Xe.c[c]に保存されます。
agent.s.operationエージェントにexpansion操作を設定することは、このエージェントに対してexpansion操作が実行されたことを示します。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Expansion (S_Agent &agent) { double Xo; double Xr; double Xe; for (int c = 0; c < coords; c++) { Xo = agent.s.c [c]; Xr = agent.s.Xr.c [c]; //Xe = Xo + RNDfromCI (0.0, expansionCoeff) * (Xr - Xo); Xe = Xo + expansionCoeff * (Xr - Xo); agent.s.Xe.c [c] = SeInDiSp (Xe, rangeMin [c], rangeMax [c], rangeStep [c]); agent.c [c] = agent.s.Xe.c [c]; } agent.s.operation = expansion; } //——————————————————————————————————————————————————————————————————————————————
「縮約」操作は、「エージェント」シンプレックス頂点のindWインデックス頂点に対して、Contractionメソッドによって実行されます。
ループの中で、Xc = Xo + contractionCoeff * (Xw - Xo)という式を使用してXcの圧縮値を計算します。ここでcontractionCoeffは収縮率であり、圧縮された頂点の初期頂点からのずれが小さくなる度合いを決定します。
次に、SeInDiSpをXcに適用して、それがrangeStep[c]のステップで有効なrangeMin[c]とrangeMax[c]の範囲内にあることを確認します。結果はagent.s.Xc.c[c]に保存されます。
agent.s.operationエージェントのcontraction操作を設定することは、このエージェントに対してcontraction操作が実行されたことを示します。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Contraction (S_Agent &agent, int indW) { double Xo; double Xw; double Xc; for (int c = 0; c < coords; c++) { Xo = agent.s.c [c]; Xw = agent.s.p [indW].c [c]; //Xc = Xo + RNDfromCI (0.0, contractionCoeff) * (Xw - Xo); Xc = Xo + contractionCoeff * (Xw - Xo); agent.s.Xc.c [c] = SeInDiSp (Xc, rangeMin [c], rangeMax [c], rangeStep [c]); agent.c [c] = agent.s.Xc.c [c]; } agent.s.operation = contraction; } //——————————————————————————————————————————————————————————————————————————————
レヴィフライトの場合、「エージェント」に対して「フライング」操作を実行するFlyingメソッドを適用します。その目標は、調査対象空間の境界に近づくにつれて確率が0まで低下し、最適な大域的解に近づく確率を集中させるランダムな分布で、シンプレックスの頂点を検索空間内の新しい場所にプッシュすることです。
Flyingメソッドを適用した後、agent.s.operationにreflection操作を設定し、その後の反復において、reflection操作が実行されたかのようにアルゴリズムのロジックが実行されるようにします。
XrにSeInDiSpを適用し、それがrangeStep[c]のステップで有効なrangeMin[c]とrangeMax[c]の範囲内にあることを確認します。結果はagent.s.Xr.c[c]に保存されます。
したがって、Flyingメソッドはエージェントに対して「飛行」操作をおこない、大域的の最適解の現在の座標とランダムに生成された値に基づいてその座標を変更します。これにより、エージェントは解空間の新しい領域を探索し、より多くの最適解を見つける可能性があります。
//—————————————————————————————————————————————————————————————————————————————— void C_AO_NMm::Flying (S_Agent &agent, int indW) { for (int c = 0; c < coords; c++) { double r1 = RNDfromCI (-1.0, 1.0); r1 = r1 > 0.5 ? 1.0 : -1.0; double r2 = RNDfromCI (1.0, 20.0); r2 = pow (r2, -2.0); if (r1 > 0.0) Xr = cB [c] + Scale (r2, 0.0, 1.0, 0.0, rangeMax [c] - cB [c], false); else Xr = cB [c] - Scale (r2, 0.0, 1.0, 0.0, cB [c] - rangeMin [c], false); //Xr = RNDfromCI (rangeMin [c], rangeMax [c]); agent.s.Xr.c [c] = SeInDiSp (Xr, rangeMin [c], rangeMax [c], rangeStep [c]); agent.c [c] = agent.s.Xr.c [c]; } agent.s.operation = reflection; } //——————————————————————————————————————————————————————————————————————————————
3.テスト結果
NMテストスタンドの結果:
C_AO_NMm:5;1.0;2.0;0.5
=============================
5 Rastrigin's; Func runs 10000 result:77.96849465506082
Score:0.96607
25 Rastrigin's; Func runs 10000 result:61.68192953373487
Score:0.76427
500 Rastrigin's; Func runs 10000 result:50.62713783555849
Score:0.62730
=============================
5 Forest's; Func runs 10000 result:0.9552378700292385
Score:0.54033
25 Forest's; Func runs 10000 result:0.45877475613538293
Score:0.25951
500 Forest's; Func runs 10000 result:0.09902427974784639
Score:0.05601
=============================
5 Megacity's; Func runs 10000 result:5.6
Score:0.46667
25 Megacity's; Func runs 10000 result:2.304
Score:0.19200
500 Megacity's; Func runs 10000 result:0.40280000000000005
Score:0.03357
=============================
All score:3.90573
アルゴリズムの結果に基づいて、滑らかなRastrigin関数のかなりまともな結果に気づくことができます。
視覚化することで、個体数の減少やエージェントの一極集中の兆候が示されます。レヴィフライトを利用することで状況はいくらか改善されます。これは、対応する座標の個々の点で顕著です。収束グラフはあまり安定していません。各関数の最初の2つのテストだけを可視化しました。1000個の変数を使った3番目のテストは実行に時間がかかりすぎて、GIFで表示することができないからです。
Rastriginテスト関数のNMm
Forestテスト関数のNMm
Megacityテスト関数のNMm
多くの欠点があるにもかかわらず、ネルダー-ミード法のアルゴリズムは9位に終わりました。
# | AO | 詳細 | Rastrigin | Rastrigin最終 | Forest | Forest最終 | Megacity(離散) | Megacity最終 | 最終結果 | ||||||
10p(5F) | 50p(25F) | 1000p(500F) | 10p(5F) | 50p(25F) | 1000p(500F) | 10p(5F) | 50p(25F) | 1000p(500F) | |||||||
1 | SDSm | 確率的拡散探索M | 0.99809 | 1.00000 | 0.69149 | 2.68958 | 0.99265 | 1.00000 | 1.00000 | 2.99265 | 1.00000 | 1.00000 | 1.00000 | 3.00000 | 100.000 |
2 | SSG | 苗木の播種と育成 | 1.00000 | 0.92761 | 0.51630 | 2.44391 | 0.72120 | 0.65201 | 0.83760 | 2.21081 | 0.54782 | 0.61841 | 0.99522 | 2.16146 | 77.683 |
3 | DE | 差分進化 | 0.98295 | 0.89236 | 0.51375 | 2.38906 | 1.00000 | 0.84602 | 0.65510 | 2.50112 | 0.90000 | 0.52237 | 0.12031 | 1.54268 | 73.099 |
4 | HS | ハーモニー検索 | 0.99676 | 0.88385 | 0.44686 | 2.32747 | 0.99148 | 0.68242 | 0.37529 | 2.04919 | 0.71739 | 0.71842 | 0.41338 | 1.84919 | 70.623 |
5 | IWO | 侵入雑草最適化 | 0.95828 | 0.62227 | 0.27647 | 1.85703 | 0.70170 | 0.31972 | 0.26613 | 1.28755 | 0.57391 | 0.30527 | 0.33130 | 1.21048 | 48.250 |
6 | ACOm | 蟻コロニー最適化M | 0.34611 | 0.16683 | 0.15808 | 0.67103 | 0.86147 | 0.68980 | 0.64798 | 2.19925 | 0.71739 | 0.63947 | 0.05579 | 1.41265 | 47.387 |
7 | MEC | MindEvolutionaryComputation | 0.99270 | 0.47345 | 0.21148 | 1.67763 | 0.60244 | 0.28046 | 0.21324 | 1.09615 | 0.66957 | 0.30000 | 0.26045 | 1.23002 | 44.049 |
8 | SDOm | 螺旋ダイナミクス最適化M | 0.69840 | 0.52988 | 0.33168 | 1.55996 | 0.59818 | 0.38766 | 0.37600 | 1.36183 | 0.35653 | 0.15262 | 0.25842 | 0.76758 | 40.289 |
9 | NMm | ネルダー-ミード法M | 0.88433 | 0.48306 | 0.45945 | 1.82685 | 0.46155 | 0.24379 | 0.21903 | 0.92437 | 0.46088 | 0.25658 | 0.16810 | 0.88555 | 39.660 |
10 | COAm | カッコウ最適化アルゴリズムM | 0.92400 | 0.43407 | 0.24120 | 1.59927 | 0.57881 | 0.23477 | 0.13842 | 0.95200 | 0.52174 | 0.24079 | 0.17001 | 0.93254 | 37.830 |
11 | FAm | ホタルアルゴリズムM | 0.59825 | 0.31520 | 0.15893 | 1.07239 | 0.50637 | 0.29178 | 0.41704 | 1.21519 | 0.24783 | 0.20526 | 0.35090 | 0.80398 | 33.139 |
12 | ABC | 人工蜂コロニー | 0.78170 | 0.30347 | 0.19313 | 1.27829 | 0.53379 | 0.14799 | 0.11177 | 0.79355 | 0.40435 | 0.19474 | 0.13859 | 0.73768 | 29.766 |
13 | BA | コウモリアルゴリズム | 0.40526 | 0.59145 | 0.78330 | 1.78002 | 0.20664 | 0.12056 | 0.21769 | 0.54488 | 0.21305 | 0.07632 | 0.17288 | 0.46225 | 29.499 |
14 | CSS | 荷電系探索 | 0.56605 | 0.68683 | 1.00000 | 2.25289 | 0.13961 | 0.01853 | 0.13638 | 0.29452 | 0.07392 | 0.00000 | 0.03465 | 0.10856 | 27.930 |
15 | GSA | 重力探索法 | 0.70167 | 0.41944 | 0.00000 | 1.12111 | 0.31390 | 0.25120 | 0.27812 | 0.84322 | 0.42609 | 0.25525 | 0.00000 | 0.68134 | 27.807 |
16 | BFO | 細菌採餌の最適化 | 0.67203 | 0.28721 | 0.10957 | 1.06881 | 0.39364 | 0.18364 | 0.17298 | 0.75026 | 0.37392 | 0.24211 | 0.18841 | 0.80444 | 27.542 |
17 | EM | 電磁気学的アルゴリズム | 0.12235 | 0.42928 | 0.92752 | 1.47915 | 0.00000 | 0.02413 | 0.29215 | 0.31628 | 0.00000 | 0.00527 | 0.10872 | 0.11399 | 19.002 |
18 | SFL | ShuffledFrog-Leaping | 0.40072 | 0.22021 | 0.24624 | 0.86717 | 0.19981 | 0.02861 | 0.02221 | 0.25063 | 0.19565 | 0.04474 | 0.06607 | 0.30646 | 13.200 |
19 | MA | モンキーアルゴリズム | 0.33192 | 0.31029 | 0.13582 | 0.77804 | 0.09927 | 0.05443 | 0.07482 | 0.22852 | 0.15652 | 0.03553 | 0.10669 | 0.29874 | 11.777 |
20 | FSS | 魚群検索 | 0.46812 | 0.23502 | 0.10483 | 0.80798 | 0.12730 | 0.03458 | 0.05458 | 0.21647 | 0.12175 | 0.03947 | 0.08244 | 0.24366 | 11.332 |
21 | IWDm | インテリジェント水滴M | 0.26459 | 0.13013 | 0.07500 | 0.46972 | 0.28358 | 0.05445 | 0.05112 | 0.38916 | 0.22609 | 0.05659 | 0.05054 | 0.33322 | 10.423 |
22 | PSO | 粒子群最適化 | 0.20449 | 0.07607 | 0.06641 | 0.34696 | 0.18734 | 0.07233 | 0.18207 | 0.44175 | 0.16956 | 0.04737 | 0.01947 | 0.23641 | 8.426 |
23 | RND | ランダム | 0.16826 | 0.09038 | 0.07438 | 0.33302 | 0.13381 | 0.03318 | 0.03949 | 0.20648 | 0.12175 | 0.03290 | 0.04898 | 0.20363 | 5.054 |
24 | GWO | 灰色オオカミオプティマイザ | 0.00000 | 0.00000 | 0.02093 | 0.02093 | 0.06514 | 0.00000 | 0.00000 | 0.06514 | 0.23478 | 0.05789 | 0.02545 | 0.31812 | 1.000 |
まとめ
ネルダー-ミードアルゴリズムは、初期段階でシンプレックスの各頂点について適合度関数を計算する必要があるため、多次元探索空間に使用する能力を否定することになり、現代の重い最適化問題ではほとんど役に立ちません。そのため、このアルゴリズムは比較的良い結果を出し、評価表にも掲載されているにもかかわらず、使用を推奨することはできません。
SDSアルゴリズムは表から削除され、その修正版だけが残されました。
図2:関連テストによるアルゴリズムのカラーグラデーション
図3:アルゴリズムのテスト結果のヒストグラム(0から100までのスケールで、多ければ多いほど良いです。
アーカイブには、記事で適用されている方法で格付けを計算するためのスクリプトが含まれています)。
修正ネルダー-ミード(NMm)アルゴリズムの長所と短所:
長所
- 外部パラメータが少ない
- 変数の数が少ない滑らかな関数に関する結果が良い
短所
- 実装が複雑
- すべてのシンプレックス頂点の適応度を事前に計算する必要があるため、サードパーティのアプリケーションで使用するのは難しい
- 操作速度が低い
- 拡張性が非常に低い
- 結果のばらつきが大きい
- 極値で立ち往生する傾向
長所と短所は、特に修正版のシンプレックス検索に関するものです。結果が非常に乏しいので、通常のNMバージョンを表に含める意味はありません。
この記事には、過去の記事で説明したアルゴリズムコードの最新版を更新したアーカイブが添付されています。記事の著者は、正規アルゴリズムの説明の絶対的な正確さについて責任を負いません。検索機能を向上させるために、それらの多くに変更が加えられています。記事に示された結論と判断は、実験結果に基づいています。
MetaQuotes Ltdによってロシア語から翻訳されました。
元の記事: https://www.mql5.com/ru/articles/13805
- 無料取引アプリ
- 8千を超えるシグナルをコピー
- 金融ニュースで金融マーケットを探索