SVM-Light interface for MATLAB/Octave

SVM-LightをMatlab/Octaveで使うためのインターフェイスを書いた．
正直言って，データ行列とクラスラベルをファイルとして保存するため，HDDの読み書きがあるので実行が遅くなるかもしれませんが，晒しておきます．

SVM-Lightのインストール

とりあえずソースからコンパイルします．
コンパイル済みのバイナリを使いたい人はそれでも良いです．

以下コンパイルしたい人向け．
SVM-Light Support Vector Machineよりsvm-light.tar.gzをダウンロード．

コンパイルします．

$ cd {Your downloaded directory}
$ mkdir svm_light
$ mv svm_light.tar.gz svm_light/.
$ gunzip -c svm_light.tar.gz | tar xvf -
$ cd svm_ligth
$ make

以上の手順を踏むと，svm_learnとsvm_classifyという実行ファイルができる．

バイナリでダウンロードした人もコンパイルした人も実行ファイル（svm_learn，svm_classify）にパスを通す．
例えば実行ファイルを/usr/local/binに置く．

$ sudo mv svm_learn /usr/local/bin/.
$ sudo mv svm_classify /usr/local/bin/.

rootじゃない人は，好きなディレクトリ（例えば，~/Myfunc/svm_light）において，パスを通しておく．

$ export PATH=~/Myfunc/svm_light/:$PATH

Windowsな人は良く分かりませんが，MATLAB/Octaveで使うディレクトリにsvm_learn.exeとsvm_classify.exeを置けばよいと思います．

以上をやってから，Terminalエミュレータ上でsvm_learnが使えるか確認する．

$ svm_learn

Not enough input parameters!

SVM-light V6.XX: Support Vector Machine, learning module
...

SVM-Lightの使い方

ググれ．

SVM-Light wrapper for MATLAB and Octave

ダウンロードは以下からどうぞ．
for svm_learn
for svm_classify
Sample code for test

svm_learn.m

function [out svm] = svm_learn(X, y, opt, samplefile, outfile)

% SVM-Light interface for MATLAB/Octave (svm_learn).
%
% INPUTS
% 	- X: d by N data matrix as [x_1 x_2, ..., x_N],
% 		where x_n is a feature vector in larning dataset.
% 	- y: N-dimensional label row vector, y_i in {-1, +1}.
% 	- opt: options for learning.
% 		Learning options
% 		c [float]: soft margin (def. (mean(x*x)^{-1})
% 		j [float]: Cost: cost-factor, by which training errors on
% 				postive examples out weighte erros on negative
% 				example (def. = 1)
% 		b {0, 1}: if 1, uses biased hyperplane (i.e. w'*x + b)
% 				if 0, uses unbiased hyperplane (i.e. w'*x)
% 				(def. = 1)
% 		i {0, 1}: remove inconsistent training examples and retrain
% 				(def. 0)
% 		Transduction options
% 		p [0, 1]: fraction of unlabeled examples to be classified
% 				into the positive class (default is the ratio of
% 				positive and negative examples in the training
% 				data)
% 		Kernel options
% 		t {0, 1, 2, 3, 4}: type of kernel function
% 			0: linear
% 			1: polynomial
% 			2: radial basis function (RBF)
% 			3: sigmoid tanh
% 			4: user defined kernel from kernel.h
% 		d [int]: parameter d in polynimial kernel
% 		g [float]: parameter gamma in RBF kernel
% 		s [float]: parameter c in sigmoid/polynomial kernel
% 		u [string]: parameter of user dfined kernel kernel
% 		Optimization options
% 		q [2, inf]: maximum size of QP-subprblems (default 10)
% 		n [2, q]: number of new variables entering the working set
% 				in each iteration (def. n = q).
% 				Set n < q to prevent zig-zagging
% 		m [5, inf]: size of cache for kernel evaluations in MB
% 				(def. = 40). The larger the faster calculating
% 		e [float]: Allow that error for termination criterion
% 				y(w'*x + b) - 1 = e (def. = 0.001)
% 		h [5, inf]: number of iterations a variable needs to be
% 				optimal before considered for shrinking
% 				(def. = 100)
% 		f {0, 1}: do final optimality check for variables removed
% 				by shrinking. Akthough this test is usually positive,
% 				there is no guarantee that the optimum was found if
% 				the test is omitted (def. = 1)
% 		y [string]: if option is given, reads alphas from file with
% 				given and uses them as starting point.
% 				(def. 'disabled')
% 		# [int]: terminate optimization, if no prgress after
% 				this number of iterations (def = 100000)

%
%% temporary files
if exist('samplefile', 'var') == 0
	samplefile = tempname(pwd);
end

if exist('outfile', 'var') == 0
	outfile = tempname(pwd);
end

%
%% options
opt_all = ['c' 'w' 'j' 'b' 'i' 'p' 't' 'd' 'g' 's' 'r' 'u',...
	'q' 'n' 'm' 'e' 'h' 'f' 'y' '#']';
opt_type = ['f' 'f' 'f' 'd' 'd' 'f' 'd' 'd' 'f' 'f' 'f' 's',...
	'f' 'd' 'd' 'f' 'd' 'd' 's' 'd']';
input_opt = cell2mat(fieldnames(opt));
Nopt = length(input_opt);
opt_str = '';
for ii = 1:Nopt
	input_val = eval(sprintf('opt.%s', input_opt(ii)));
	opt_str = make_opt(input_opt(ii), input_val,...
		opt_all, opt_type, opt_str);
end

%
% write sample file
[d N] = size(X);

fid = fopen(samplefile, 'w');
if fid < 0
	error(sprintf('cannot open %s', samplefile));
end
for ii = 1:N
	fprintf(fid, '%d ', y(ii));
	fprintf(fid, '%d:%e ', [[1:d]; X(:, ii)']);
	fprintf(fid, '\n');
end
fclose(fid);
%system(sprintf('cat %s', samplefile))

%
% excute svm_learn
sysc = sprintf('svm_learn -v 0 -z c %s %s %s',...
	opt_str, samplefile, outfile);
[a b] = system(sysc);
delete(samplefile)

%
% read output file of svm_learn to cell array
%system(sprintf('cat %s', outfile))
fid = fopen(outfile, 'r');
if fid < 0
	error(sprintf('cannot open %s', outfile));
end
cnt = 1;
outl1 = 0;
while outl1 ~= -1;
	outl1 = fgetl(fid);
	if outl1 ~= -1
		out{cnt} = outl1;
		cnt = cnt + 1;
	end
end
fclose(fid);
delete(outfile);

%
% read w and b
svm.version = out{1};
svm.kernel = str2num(strtok(out{2}, ' '));
svm.d = str2num(strtok(out{3}, ' '));
svm.g = str2num(strtok(out{4}, ' '));
svm.s = str2num(strtok(out{5}, ' '));
svm.r = str2num(strtok(out{6}, ' '));
svm.u = strtok(out{7}, '#');
svm.d = str2num(strtok(out{8}, ' '));
svm.N = str2num(strtok(out{9}, ' '));
svm.Nsv = str2num(strtok(out{10}, ' ')) - 1;
svm.b = str2num(strtok(out{11}, ' '));

alpha = zeros(1, svm.Nsv);
sv = zeros(svm.d, svm.Nsv);

for ii=1:svm.Nsv
	val = out{ii + 11}; % SVs are written from line 12;
	val(val==':')= ' ';
	val(val=='#')= ' ';
	val = str2num(val);
	alpha(ii) = val(1);
	sv(:, ii) = val(3:2:end)';
end
svm.w = alpha*sv';


%% ==================================
%% making string for options
function str = make_opt(name, val, opt_all, opt_type, str)
egs = findstr(opt_all, name);
if length(egs) == 1
	type_str = sprintf('%%%s', opt_type(egs));
	str = sprintf(['%s -%s ' type_str], str, name, val);
else
	error('Invalid option is inputted')
end

簡単に解説すると，
サンプルデータはXという行列．
サンプルデータのラベルはベクトルy．
入力変数optは構造体になってて，そこにSVM-Lightに使うオプションを並べれば良い．
一時ファイル名は入れても入れなくても良い．Linuxだったらtmpfsを使えば早く読み書きできるかもしれない．
出力をそのままsvm_classify.mに入れる．
svm_classifyを使いたくないときは，出力svm.wとsvm.bなどが使えます．

svm_classify.m

function [labs z] = svm_classify(X, y, model, samplefile, outfile, modelfile)

% SVM-Light interface for MATLAB/Octave (svm_classify).
%
% INPUTS
% 	- X: d by N data matrix as [x_1 x_2, ..., x_N],
% 		where x_n is a feature vector in test samples.
% 	- y: N-dimensional label row vector, y_i in {-1, +1}.
% 	- model: Model of SVM.

%
%% temporary files
if exist('samplefile', 'var') == 0
	samplefile = tempname(pwd);
end

if exist('outfile', 'var') == 0
	outfile = tempname(pwd);
end

if exist('modelfile', 'var') == 0
	modelfile = tempname(pwd);
end

%
%% correct labels
[d N] = size(X);
if length(y) == 0
	y = zeros(1, N);
end

%
% write sample file
fid = fopen(samplefile, 'w');
if fid < 0
	error(sprintf('cannot open %s', samplefile));
end
for ii = 1:N
	fprintf(fid, '%d ', y(ii));
	fprintf(fid, '%d:%e ', [[1:d]; X(:, ii)']);
	fprintf(fid, '\n');
end
fclose(fid);

%
% write model file
fid = fopen(modelfile, 'w');
if fid < 0
	error(sprintf('cannot open %s', modelfile));
end
Nml = size(model, 2);
for ii = 1:Nml
	fprintf(fid, '%s', model{ii});
	fprintf(fid, '\n');
end
fclose(fid);
%system(sprintf('cat %s', modelfile))

%
% excute svm_classify
sysc = sprintf('svm_classify -v 0 %s %s %s',...
	samplefile, modelfile, outfile);
[a b] = system(sysc);
delete(samplefile)
delete(modelfile)
%system(sprintf('cat %s', outfile))

%
% read out file
fid = fopen(outfile, 'r');
if fid < 0
	error(sprintf('cannot open %s', outfile));
end
cnt = 1;
outl1 = 0;
while outl1 ~= -1;
	outl1 = fgetl(fid);
	if outl1 ~= -1
		z(1, cnt) = str2num(outl1);
		cnt = cnt + 1;
	end
end
fclose(fid);
delete(outfile);
labs(find(z < 0)) = -1;
labs(find(z >= 0)) = 1;

サンプルデータはXという行列．
サンプルデータのラベルはベクトルy．

test_svm_light.m

clear all
randn('state',1);

N = 20; % number of samples for a class
X1 = randn(2, N) + [2 2]'*ones(1, N);
X2 = randn(2, N) + [0 0]'*ones(1, N);

%X1 = randn(2, N);
%X1 = X1./repmat(sqrt(diag(X1'*X1))', 2, 1) + .5*randn(2, N);
%X2 = randn(2, N);
%X2 = X2./repmat(sqrt(diag(X2'*X2))'*.5, 2, 1) + .5*randn(2, N);

X = [X1 X2]; % input learning samples
y = [ones(1, N) -1*ones(1, N)]; % input learning labels

% options for SVM
opt.t = 0; % 0: linear, 1: polynomial, 2: radial basis function (RBF),
		   % 3: sigmoid tanh, 4: user defined kernel from kernel.h
opt.b = 1; % use bias
opt.g = 1; % gamma for RBF kernel
opt.d = 2; % order for polynomial kernel
opt.c = 10; % soft margin parameter

% svm learn
[out svm] = svm_learn(X, y, opt);

% svm classify for training accuracy
[out_label z] = svm_classify(X, [], out);
acc = length(find((out_label - y) == 0))/(2*N);
fprintf('training accuracy: %.2f%%\n', acc*100)

% plot hyperplane
x = -4:0.1:4;
[XX YY] = meshgrid(x);
Xmesh(1, :) = reshape(XX, 1, prod(size(XX)));
Xmesh(2, :) = reshape(YY, 1, prod(size(XX)));
[out_label z] = svm_classify(Xmesh, [], out);
Xmesh1 = Xmesh(:, find(out_label == 1));
Xmesh2 = Xmesh(:, find(out_label == -1));

%% plot results
plot(Xmesh1(1, :),  Xmesh1(2, :),  'y.',...
	Xmesh2(1, :),  Xmesh2(2, :),  'r.',...
	X1(1, :),  X1(2, :),  'o', 'linewidth', 5,...
	X2(1, :),  X2(2, :),  's', 'linewidth', 5)
grid on
axis([-4 4 -4 4])
xlabel('x(1)')
ylabel('x(2)')

サンプルデータの学習をして，識別境界をプロットする．

実行例

サンプルデータを

X1 = randn(2, N) + [2 2]'*ones(1, N);
X2 = randn(2, N) + [0 0]'*ones(1, N);

という感じに作り，線形カーネルによって境界面を作ると下図のように識別できる．
（青と緑がそれぞれのトレーニングサンプル，赤と黄の境界が識別境界）

サンプルデータを

X1 = randn(2, N);
X1 = X1./repmat(sqrt(diag(X1'*X1))', 2, 1) + .5*randn(2, N);
X2 = randn(2, N);
X2 = X2./repmat(sqrt(diag(X2'*X2))'*.5, 2, 1) + .5*randn(2, N);

という感じに作り，2次の多項式カーネルで境界面作ると，

となる．
同じデータで，ガウスカーネルによって識別すると，

となる．過学習ぎみデスネ．

以上です．
LIBSVMにもMATLAB interfaceがあるようなので，それも試してみる予定．
あと，SVM-KMというのもある．