最急降下法!最急降下法!


GA大好きな私ですが、古典的な方法で最適化も試しておかなくてはなりません。
機械学習でも最急降下法がよく使われているとかで、この際便利なコードを書いておきましょう。とはいえ過去の資産を使わない手はありません。ググりましょう。
非線形最適化テスト
最急降下法は初期値の微分値を元に次の解のポイントを見つけます。f(x,y)に関する最小化問題を考えると、ある点x1,y1で、微分値x1′,y1’はその点での傾きを示します。それが正の数であれば、次に探索すべきポイントはx1,y1からλx1,λy1を引いた数となります。λは収束のための係数です。微係数x1′,y1’は解に近づけば小さくなるはずなので、0.0001とか小さな値を指定しておけばいいです。
しかし小さな値をかけると改善量は小さくなります。大きな値を指定すると最適解を通り越して振動してしまいます。

よって微係数が正負逆転した時にλを小さくすれば振動は収まるはずです。

f(x,y)=(x-3.5)*(x-3.5)+(y-4.8)*(y-4.8)+500

という関数についてx=50000,y=40000,λ=0.0000001で実験したところ、9333万回繰り返してなんとか近似解がでました。これに対してλを改善する方法ではx=50000,y=40000,λ=1で実験すると58回で最適解となりました。初期値に依存しないところがまあいい感じですね。最適化手法としてはアニーリングに近いのかもしれません。


#include <stdio.h>
#include <stdlib.h>
#include <math.h>

//評価する関数f
double f(double x, double y)
{
    return (x-3.5)*(x-3.5)+(y-4.8)*(y-4.8)+500;
}
// ∇f(x0, y0)を求める      
double* getGradF(double x, double y) {
    // 微小数
    static double result[2];  
    double h = 0.000001;  
    // xで偏微分
    double base=f(x,y);
    result[0] = (f(x+h,y  ) - base) / h;
    result[1] = (f(x  ,y+h) - base) / h;  
    return result;
}
//最急降下法
double fastdown(double &x0, double &y0){
    // 更新用  
    double e      = 0.0000000001;   // 収束条件  
    double lambda = 1;
    double distance;
    double gradF_old[2]={0};
    double x1;
    double y1;
    int cnt=0;
    while(true) {  
        double* gradF = getGradF(x0,y0);  
        // 更新  
        x1 = x0 - lambda*gradF[0];
        y1 = y0 - lambda*gradF[1];
        if( gradF[0]*gradF_old[0]<0 || gradF[1]*gradF_old[1]<0){
            lambda/= 2;
        }
        gradF_old[0] = gradF[0];
        gradF_old[1] = gradF[1];
        // 収束したら終了  
        distance = sqrt((x1-x0)*(x1-x0)+(y1-y0)*(y1-y0));
        if ( distance < e ){
            break;
        }
        //値を再設定
        x0 = x1;
        y0 = y1;
        cnt++;
    }
    printf( "count=%d\n", cnt);  
    return f(x0,y0);

}
//メイン
int main(int argc, char* argv[])
{
   double x = 50000;
   double y = 40000;
   double result;
   result = fastdown(x,y);
   printf( "x:%lf y:%lf diff:%lf\n", x,y,result);
   system("pause");
   return 0;
}

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です