最急降下法

理論はクソ簡単なので, 特に説明しません. 坂道転がってくだけです.
実装では自分で微分した式を使っちゃってるので, すごいダサく仕上がってます. 将来的にはn次元空間でも使えるように, 微分の計算部分も実装したい.

zの空間は, こちらの記事をそのまま使っちゃってて,

$z = 3x^2 + 5y^2 - 6xy$

となっています. 一応グラフにしとくとこんな感じ.
f:id:ket-30:20170202231513p:plain:w400:h400
まだgnuplotに慣れてなくて, 汚い図です. いずれgnuplotの使い方をまとめるかも.

そして当たり前ですが

$z \geq 0, \,\,z(0,0) = 0$

なので, 目的となる点は(0,0)です.
何度か勾配計算 & 更新を繰り返して, 0に辿りつければ成功です.

結果 & コード

初期値は適当に(-7, 4)に落として, そのあと100回更新しました. 更新幅でbreakさせた方が良かったですね. 他の手法を実装する時そうします. f:id:ket-30:20170202233824p:plain:w400:h100
まぁとりあえずオッケー!

以下全コード

extern crate gnuplot;

use gnuplot::*;
use std::vec::Vec;

fn main() {
    // --- 最適化したい空間の設定 & 作図 ---
    let x = linspace(-3.0, 3.0, 1000);
    let y = linspace(-3.0, 3.0, 1000);

    let z = make_zspace(&x, &y);

    let mut fg = Figure::new();
    fg.axes3d()
        .surface(z.iter(), x.len(), y.len(), Some((-4.0, -4.0, 4.0, 4.0)), &[])
        .show_contours(true, true, Linear, Fix(""), Auto)
        .set_title("z space", &[]);
    fg.set_terminal("png", "zspace.png");
    fg.show();

    let mut init_p: Vec<f32> = vec![-7.0, 4.0];
    let mut gd = GradientDescent::new(init_p, 0.1);

    for i in 0..100 {
        gd.update();
        println!("{:?}", gd.p);
    }
}

// 連続なベクトルの生成
fn linspace(start: f32, limit: f32, n: usize) -> Vec<f32> {
    if n <= 1 {
        vec![limit]
    } else {
        let mut v = vec![];
            let s = start;
            let d = (limit - start) / (n - 1) as f32;
            for i in (0..n) {
                v.push(s + i as f32 * d);
            }
        v
    }
}

// 最適化したい空間の生成
fn make_zspace (x: &Vec<f32>, y: &Vec<f32>) -> Vec<f32> {
    let mut z = Vec::with_capacity((x.len()*y.len()) as usize);
    for xi in x.iter() {
        for yi in y.iter() {
            z.push( 3f32*xi.powi(2) + 5f32*yi.powi(2) - 6f32*xi*yi );
        }
    }
    z
}

struct GradientDescent {
    p: Vec<f32>,
    eta: f32,
}

impl GradientDescent {

    fn new(p: Vec<f32>, eta: f32) -> GradientDescent {
        GradientDescent{p: p, eta: eta}
    }

    fn update (&mut self, ) {
        self.p[0] -= self.eta*self.calc_grads()[0];
        self.p[1] -= self.eta*self.calc_grads()[1];
    }

    // 勾配も&で受けれるようにメンバーに加えるべきだったかも.
    // 毎回grad_vを作成しているのはセンスない.
    // 次から気をつける.
    fn calc_grads(&mut self) -> Vec<f32> {
        let mut grad_v = vec![];
        grad_v.push( 6f32 * self.p[0] - 6f32 * self.p[1]);
        grad_v.push(10f32 * self.p[0] - 6f32 * self.p[1]);

        grad_v
    }
}