目的

Octave でちゃちゃっと SVM を試してみます．

データの準備

関数と分類するためのデータをちゃちゃっと用意します．A と B はそれぞれ適当な位置を中心にして適当な幅を持った分布からサンプリングします．rot() は分布を回転させるための行列です．特徴量は X で分類のためのラベルは Y です．

# 回転行列
rot = @(t) [cos(t) -sin(t); sin(t) cos(t)];

# 適当にデータを用意します
A = bsxfun(@plus, ...
     randn(300,2)*eye(5*rand(1,2))*rot(pi*rand()), 10*rand(2,1));
B = bsxfun(@plus, ...
     randn(300,2)*eye(5*rand(1,2))*rot(pi*rand()), 10*rand(2,1));

# SVM 計算用にデータを整理します
X = [[A;B], ones(600,1)];
Y = [ones(300,1); -ones(300,1)];

最適化

SVM によって xy 平面にデータを 2 分するための直線を引きます．

w_1x + w_2y + w_3 = 0

ベクトル w を求めます．ソフトマージン SVM では以下の式を w について最小化することで分類するためのパラメタ w を求めることができます．

\mathcal{L}({\bf w}) = \sum_n \max
\bigl(
1 - Y_n \left({\bf X}_n \cdot {\bf w} \right) 
\bigr)
+ C ||{\bf w}||_{2}

Octave の表式に合わせて目的関数 L を定義します． C は正則化項の係数です．

L = @(w,C) sum(max(1.0-Y.*(X*w(:)),0)) + C*sumsq(w);

L の最小値を求めるために sqp を利用します．関数を線形近似して iterative に二次計画法を解くことによって，非線形関数の最小値を求めているようです．初期値は適当に入れます．

# 正則化項の係数を C = 0.01 にした無名関数を引数に渡します
w = sqp([-1;1;1], @(w) L(w,1e-2))

# とあるデータに対する結果
# w =
# 
#   -0.82591
#    0.97024
#    0.82691
#

せっかくなので可視化してみます．

# 可視化のために dirty な関数を用意
Lplot = @(x,y) log10((reshape(cellfun(@(w) L(w,1e-5), ...
          mat2cell([x(:),y(:),zeros(size(x(:)))], ...
          ones(numel(x),1),3)),size(x,1),size(x,2))));

# プロットしてみる
subplot(1,2,2);
  # w(3) = 0 平面での目的関数の値をカラーマップで表示
  ezsurf(Lp,[-8,2,-4,6],64); view(2);
  title("Objective Function (at w(3)=0)");
  xlabel("w(1)"); ylabel("w(2)");
subplot(1,2,1);
  # A, B と分離平面(直線)を表示
  plot(A(:,1),A(:,2),'1o', B(:,1), B(:,2),'3p');
  hold('on');
  ezplot(@(x,y) w(1)*x+w(2)*y+w(3),[-8,14,-8,15]);
  hold('off');

それっぽいところに直線が引けました．最適化の初期値やグラフの軸は適当に調整してください．

Octave でソフトマージン SVM を計算してみる

目的

データの準備

最適化