在MATLAB中，应用软间隔SVM算法（非线性）对鸢尾花数据进行分类，并画图
%% 数据预处理和导入
close;clear;clc
[train_data,train_label,test_data,test_label,m1,n1,m2,n2] = data_set(0.6,100);
%% 模型训练
Kernel = 'linear';% Kernel 核技巧备选：gaussian linear sigmoid mullinear
svm = train_svm(train_data',train_label',Kernel,10); % svm = train_svm(X,Y,kertype,C) C为变量上界（惩罚因子） svm为结构体
%% 模型测试
result = test_svm(svm,test_data',test_label',Kernel);
fprintf('训练完成！\n应用模型：SVM 支持向量机\n优化算法:interior-point-convex\n核函数：%s\n测试集识别率为：%f\n',Kernel,result.accuracy);
%% 作图显示数据以及训练结果；中间为支持向量[三维]
draw_svm(train_data,train_label,svm,3,Kernel);
function [train_data,train_label,test_data,test_label,m1,n1,m2,n2] = data_set(k,data_num)
% 功能说明：完成数据的预处理，setosa：1 versicolor：2 virginica：3
% 调用语法及参数说明：[data_iris,data_label] = data_set();
%
load('data_iris.mat');load('data_label.mat');
data_label = zeros(data_num,1);
for i = 1:data_num
switch species(i)
case 'setosa'
data_label(i) = 1;
case 'versicolor'
data_label(i) = -1;
% case 'virginica'
% data_label(i) = 3;
end
end
data_iris = iris(1:data_num,:);
% 乱序排列
randIndex = randperm(data_num);
data_new=data_iris(randIndex,:);
label_new=data_label(randIndex,:);
% 分为两组，比例k用于训练，剩余用于测试
k = k*data_num;
train_data=data_new(1:k,:);
train_label=label_new(1:k,:);
test_data=data_new(k+1:end,:);
test_label=label_new(k+1:end,:);
[m1,n1] = size(train_data);
[m2,n2] = size(test_data);
end
function draw_svm(train_data,train_label,svm,data_features,Kernel)
% 功能说明：根据数据特征的维数判断，进而分别绘图
% 函数语法及参数列表：draw_svm(inputArg1,data_features)
% input:
% train_data: 训练数据集
% train_label：训练集数据的类别
% svm：svm结构体(详见train_svm，help train_svm)
% data_features特征维数
switch data_features
case 2
plot(train_data(train_label==1,1),train_data(train_label==1,2),'ro',train_data(train_label==-1,1),train_data(train_label==-1,2),'go');hold on;
plot(svm.data(1,:),svm.data(2,:),'mo');hold on;title(['样本分布',Kernel]); % 显示支持向量 'mo'品红色的圈
[x1,x2] = meshgrid(4:0.01:7,2:0.01:5);
[rows,cols] = size(x1);
nt = rows*cols;
Xt = [reshape(x1,1,nt);reshape(x2,1,nt)];
Yt = ones(1,nt);
result = test_svm(svm, Xt, Yt,Kernel);
Yd = reshape(result.Y,rows,cols);
contour(x1,x2,Yd,'m');
otherwise
plot3(train_data(train_label==1,1),train_data(train_label==1,2),train_data(train_label==1,3),'r.');hold on;
plot3(train_data(train_label==-1,1),train_data(train_label==-1,2),train_data(train_label==-1,3),'gx');hold on;
plot3(svm.data(1,:),svm.data(2,:),svm.data(3,:),'mo');hold on;
title(['样本分布',Kernel]);
end
end
function K = kernel(X,Y,kerneltype)
% 功能：支持多种核运算;
% 语法结构：K = kernel(X,Y,kerneltype)，kerneltype选择核技巧
% 'linear'：线性内积
% K(v1,v2) = <v1,v2>
% 'gaussian'：高斯核 %
% K(v1,v2)=exp(-gama||v1-v2||^2)
% 'sigmoid':sigmoid核；双曲正切函数
% K(v1,v2)=tanh(gama<v1,v2>+c)
% 'mullinear':多项式核
% K（v1,v2）=<v1,v2>^d;d为多项式的次数
% 'triangle':三角核
% K（v1,v2）=-||v1-v2||^d
% 在svm中运用线性，高斯或者sigmoid效果比较好
switch kerneltype
case 'linear' % 线性内积
K = X'*Y;
case 'sigmoid'
belta = 0.01;
theta = 0.001;
K = tanh(belta*X*Y+theta);
case 'gaussian'% k(v1,v2) = exp(-||v1-v2||^2/(2sigma^2))
delta = 2*1.414;
delta = delta*delta;
XX = sum(X'.*X',2);
YY = sum(Y'.*Y',2);
XY = X'*Y;
K = abs(repmat(XX,[1 size(YY,1)]) + repmat(YY',[size(XX,1) 1]) - 2*XY);
K = exp(-K./delta);
case 'mullinear'
K = (X'*Y).^2;
% case'triangle'
% K = -norm(X-Y,1)^2;
end
end
function result = test_svm(svm, test_data, test_label, kerneltype)
% 功能说明：
% 完成测试集的预测以及准确率的输出
% 语法习惯核参数列表：result = test(svm, test_data, test_label, kerneltype)
% input:
% svm: train_svm函数返回的结构体（详见help train_svm）
% test_data: 测试数据
% test_label：测试集标签
% kerneltype：核技巧种类，形式参数，可选：linear gaussian sigmoid mullinear triangle
% output:
% result：结构体，属性如下
% result.Y:测试集中数据的预测类别 result.Y ∈{+1，-1}
% result.accuracy:测试集的准确率
% 教材非线性支持向量机学习算法的策略为选择a的一个正分量0< a <C进行计算
% 此处选择了对所有满足0< ai <C求得bi，并对b进行取平均运算
sum_b = svm.label - (svm.a'.* svm.label)*kernel(svm.data,svm.data,kerneltype);
b = mean(sum_b);
w = (svm.a'.* svm.label)*kernel(svm.data,test_data,kerneltype);% 统一起见，令 w = sigma(ai*yi*K(x,xi)
result.Y = sign(w+b);% 加外壳符号函数进行分类
result.plotx = min(test_data(1,:)):0.001: max(test_data(1,:));
result.accuracy = size(find(result.Y==test_label))/size(test_label);% 预测正确的数据数目/总测试集数目
end
function svm = train_svm(train_data,train_label,kertype,C)
% 功能说明：完成SVM训练
% 语法习惯与参数列表：svm = train_svm(train_data,train_label,kertype,B)
% input:
% train_data:训练数据
% train_label：训练数据的类别
% kertype：核函数的类别
% C 惩罚参数
% B 为变量约束中的上界
% output:
% svm：是一个结构体，包含属性如下：
% svm.a :得到的凸二次规划的解
% svm.data ： 支持向量
% svm.label ：支持向量的类别
% ------------*************************···········
% ------------关键函数quadprog的一些说明···········
% 函数quadprog：用于解二次规划问题
% 问题描述：
% min（x）: 0.5·x'·H·x + f'·x
%
% A·x <= b,
% s.t.: Aeq·x = beq;
% lb <=x <= ub;
%
% 全参数语法结构：x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options);
% 变量说明：
% H,A,Qeq是矩阵，f,b,beq,lb,ub,x是向量
% options:选择优化算法并进行设置
% 优化选项设置，对属性进行设置：
% 使用 optimoptions 或 optimset 创建选项(属性)；
% 指定为 optimoptions 的输出或 optimset 等返回的结构体。
% 变量初始化以及超参设置
n = length(train_label); % 对变量的自由约束，上下界
H = (train_label'*train_label).*kernel(train_data,train_data,kertype);% H为yi*yj*K(xi,xj)
f = -ones(n,1); % 保证f为列向量，原式中包含转置操作
A = [];% 不含不等约束
b = [];% 不含不等约束
Aeq = train_label; % s.t.: aY = 0;
beq = 0; % s.t.: aY = 0;
lb = zeros(n,1); % 解：a 的范围
ub = C*ones(n,1); % 0 <= a <= C
a0 = zeros(n,1); % a0是解的初始近似值
options = optimset; % 'interior-point-convex'（默认优化算法）
options.Display = 'iter-detailed';% 显示优化步骤
% x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options) 使用 options 中指定的优化选项求解上述问题。
% 使用 optimset 创建 options。如果不提供初始点，设置 x0 = []。
a = quadprog(H,f,A,b,Aeq,beq,lb,ub,a0,options);
% 寻找支持向量；a>e 则认为train_data为支持向量 函数find查找逻辑为真的索引
e = 1e-4;
sv_index = find(abs(a)>e);
svm.a = a(sv_index);
svm.data = train_data(:,sv_index);% 作图显示支持向量位置
svm.label = train_label(sv_index);
end