Monday, March 4, 2013

Support Vector Classifier in MATLAB

% the following code is for the book The Elements of Statistical Learning, section 12.2
% choose right directories before and after this step
cvx_setup

% cvx can be downloaded @ http://cvxr.com/cvx/
%%%%%%%%%%%
%separable case
clf
var = 1.25;
n = 100; mu1 = [3; 8]; mu2 = [8; 3];
x1 = repmat(mu1',n/2,1) + randn(n/2,2)*sqrt(var);
x2 = repmat(mu2',n/2,1) + randn(n/2,2)*sqrt(var);
x = [x1; x2];
Y = [repmat(1,50,1); repmat(-1,50,1)];
hold on
scatter(x1(:,1),x1(:,2),'or','filled')
scatter(x2(:,1),x2(:,2),'ob','filled')
axis([min(min(x(:,1))) max(max(x(:,1))) min(min(x(:,2))) max(max(x(:,2)))]);
hold off
 


cvx_begin
cvx_precision('high')
variable beta0;
variable beta(2);
minimize( .5*square_pos(norm(beta)));
subject to
Y.*(x*beta + beta0) >= 1;
cvx_end

hold on
scatter(x1(:,1),x1(:,2),'or','filled')
scatter(x2(:,1),x2(:,2),'ob','filled')
axis([min(min(x(:,1))) max(max(x(:,1))) min(min(x(:,2))) max(max(x(:,2)))]);
xs = linspace(min(min(x)),max(max(x)),1000);
plot(xs,(-beta0 - xs*beta(1))/beta(2),'k')
M = 1/norm(beta);
plot(xs,(-beta0 - xs*beta(1)+1)/beta(2) ,'--k')
plot(xs,(-beta0 - xs*beta(1)-1)/beta(2),'--k')
hold off


%%%%%%%%%%%%%%
%non-seperable case
clf
var = 2.5;
n = 100; mu1 = [3; 8]; mu2 = [8; 3];
x1 = repmat(mu1',n/2,1) + randn(n/2,2)*sqrt(var);
x2 = repmat(mu2',n/2,1) + randn(n/2,2)*sqrt(var);
x = [x1; x2];
Y = [repmat(1,50,1); repmat(-1,50,1)];
hold on
scatter(x1(:,1),x1(:,2),'or','filled')
scatter(x2(:,1),x2(:,2),'ob','filled')
axis([min(min(x(:,1))) max(max(x(:,1))) min(min(x(:,2))) max(max(x(:,2)))]);
hold off

gam = .1;
cvx_begin
cvx_precision('high')
variable beta0;
variable beta(2);
variable epsilon(n);
minimize( .5*square_pos(norm(beta)) + gam*sum(epsilon));
subject to
for i=1:n
    Y(i)*(x(i,:)*beta + beta0) >= 1 - epsilon(i);
end
epsilon >= 0;
cvx_end
clf
hold on
scatter(x1(:,1),x1(:,2),'or','filled')
scatter(x2(:,1),x2(:,2),'ob','filled')
axis([min(min(x(:,1))) max(max(x(:,1))) min(min(x(:,2))) max(max(x(:,2)))]);
xs = linspace(min(min(x)),max(max(x)),1000);
plot(xs,(-beta0 - xs*beta(1))/beta(2),'k')
M = 1/norm(beta);
plot(xs,(-beta0 - xs*beta(1)+1)/beta(2) ,'--k')
plot(xs,(-beta0 - xs*beta(1)-1)/beta(2),'--k')
ind = epsilon>1e-6;
scatter(x(ind,1),x(ind,2),'+k','SizeData',150)
ind = abs(x*beta + beta0 + 1)<1e-4;
scatter(x(ind,1),x(ind,2),'Xk','SizeData',150)
ind = abs(x*beta + beta0 - 1)<1e-4;
scatter(x(ind,1),x(ind,2),'Xk','SizeData',150)
hold off

No comments:

Post a Comment