SVM非线性分类原理实验

一、总结

一句话总结：

到更高维度去找可以分类的超平面（无限维度的平面中必然可分）

二、解密SVM系列（四）：SVM非线性分类原理实验

转自或参考：解密SVM系列（四）：SVM非线性分类原理实验 - gccbuaa - 博客园
https://www.cnblogs.com/gccbuaa/p/7299934.html

前面几节我们讨论了SVM原理、求解线性分类下SVM的SMO方法。本节将分析SVM处理非线性分类的相关问题。

一般的非线性分类例如以下左所看到的（后面我们将实战以下这种情况）：
这里写图片描写叙述
能够看到在原始空间中你想用一个直线分类面划分开来是不可能了，除非圆。而当你把数据点映射一下成右图所看到的的情况后，如今数据点明显看上去是线性可分的，那么在这个空间上的数据点我们再用前面的SVM算法去处理，就能够得到每一个数据点的分类情况了，而这个分类情况也是我们在低维空间的情况。也就是说，单纯的SVM是不能处理非线性问题的，说白了仅仅能处理线性问题，可是来了非线性样本怎么办呢？我们是在样本上做的文章，我把非线性样本变成线性样本，再去把变化后的线性样本拿去分类，经过这么一圈。就达到了把非线性样本分开的目的。所以仅仅看开头和结尾的话发现，SVM居然能够分非线性问题，事实上呢还是分的线性问题。

如今的问题是怎样找到这个映射关系对吧，就比方上面那个情况，我们能够人为计算出这种映射，比方一个样本点是用坐标表示的(x1,x2),它有个类标签，如果为1。那么把这个点映射到三维中变成(x21,2√x1x2,x22)。对每一个点我都这么去映射，如果一个原始点样本集是这种：
这里写图片描写叙述
然后依照上面那个公式去把每一个点映射成3维坐标点后，画出来是这种：

能够看到是线性可分的吧。如果还看不清把视角换个角度（右视图）：

如今能看清楚了吧。

那这是二维的点到三维。映射的关系就是上面的那个关系，那如果是三维到四维，四维到N维呢？这个关系你还想去找吗？理论上是找的到的。可是实际上人工去找你怎么去找？你怎么知道数据的映射关系是这种是那样的？不可能知道。然而我们真的须要找到这种关系吗？答案是不须要的。返回去看看前三节的关于SVM的理论部分能够看到，不管是计算ααα呀。还是b呀等等，仅仅要涉及到原始数据点的。都是以内积的形式出来的，也就是说是一个点的向量与还有一个点的向量相乘的。向量内积出来是一个值。

K(X1,X2)=XT1∗X2=(x1,y1,z1)T∗(x2,y2,z2)=x1y1+x2y2+x3y3=C2

最后也是得到一个值比方C2。既然SVM里面全部涉及到原始数据的地方都是以向量的形式出现的，那么我们还须要管它的映射关系吗？由于它也不须要你去计算说详细到比方说三维以后。三维里面的三个坐标值到底是多少，他须要的是内积以后的一个结果值。那么好办了，我就如果有一个黑匣子，输入原始数据维度下的两个坐标向量。然后经过黑匣子这么一圈，出来一个值，这个值我们就觉得是高维度下的值。

而黑匣子的潜在意义就相当于一个高维映射器一样。

更重要的是我们并不须要知道黑匣子到底是怎么映射的，仅仅须要知道它的低纬度下的形式就能够了。经常使用的黑匣子就是径向基函数。而这个黑匣子在数学上就叫做核函数，比如径向基函数的外在形式例如以下所看到的：

好了既然黑匣子是藏着的。那也就仅仅能说这么多了。

有趣的是上帝给的这个黑匣子不止一个，有好几个，仅仅是上面的那个普遍效果更好而已。

基于此，那么对于上节的SMO算法，如果拿来求解非线性数据的话，我们仅仅须要将当中相应的内积部分改成核函数的形式就可以。一个数据核函数程序例如以下：

function result = Kernel(data1,data2,sigma)
% data里面每一行数据是一个样本(的行向量)
[m1,~] = size(data1);
[m2,~] = size(data2);
result = zeros(m1,m2);
for i = 1:m1
    for j = 1:m2
        result(i,j) = exp(-norm(data1(i,:)-data2(j,:))/(2*sigma^2));
    end
end

有了此核函数。我们用上节的随机遍历α的方式（这个函数代码少一点）来实验一下非线性样本，非线性样本例如以下：
这里写图片描写叙述
然后把主程序相应的部分用上述核函数取代：

%%
% * svm 简单算法设计
%
%% 载入数据
% * 终于data格式：m*n，m样本数。n维度
% * label:m*1  标签必须为-1与1这两类
clc
clear
% close all
data = load('data_test1.mat');
data = data.data;
train_data = data(1:end-1,:)';
label = data(end,:)';
[num_data,d] = size(train_data);
data = train_data;
%% 定义向量机參数
alphas = zeros(num_data,1);
% 系数
b = 0;
% 松弛变量影响因子
C = 0.6;
iter = 0;
max_iter = 80;
% 核函数的參数
sigma = 4;
%%
while iter < max_iter
    alpha_change = 0;
    for i = 1:num_data
        %输出目标值
        pre_Li = (alphas.*label)'*Kernel(data,data(i,:),sigma) + b;
        %样本i误差
        Ei = pre_Li - label(i);
        % 满足KKT条件
        if (label(i)*Ei<-0.001 && alphas(i)<C)||(label(i)*Ei>0.001 && alphas(i)>0)
           % 选择一个和 i 不同样的待改变的alpha(2)--alpha(j)
            j = randi(num_data,1);  
            if j == i
                temp = 1;
                while temp
                    j = randi(num_data,1);
                    if j ~= i
                        temp = 0;
                    end
                end
            end
            % 样本j的输出值
            pre_Lj = (alphas.*label)'*Kernel(data,data(j,:),sigma) + b;
            %样本j误差
            Ej = pre_Lj - label(j);
            %更新上下限
            if label(i) ~= label(j) %类标签同样
                L = max(0,alphas(j) - alphas(i));
                H = min(C,C + alphas(j) - alphas(i));
            else
                L = max(0,alphas(j) + alphas(i) -C);
                H = min(C,alphas(j) + alphas(i));
            end
            if L==H  %上下限一样结束本次循环
                continue;end
            %计算eta
            eta = 2*Kernel(data(i,:),data(j,:),sigma)- ...
                Kernel(data(i,:),data(i,:),sigma)...
                - Kernel(data(j,:),data(j,:),sigma);
            %保存旧值
            alphasI_old = alphas(i);
            alphasJ_old = alphas(j);
            %更新alpha(2)，也就是alpha(j)
            alphas(j) = alphas(j) - label(j)*(Ei-Ej)/eta;
            %限制范围
            if alphas(j) > H
                alphas(j) = H;
            elseif alphas(j) < L
                    alphas(j) = L;
            end
            %如果alpha(j)没怎么改变，结束本次循环
            if abs(alphas(j) - alphasJ_old)<1e-4
                continue;end
            %更新alpha(1)。也就是alpha(i)
            alphas(i) = alphas(i) + label(i)*label(j)*(alphasJ_old-alphas(j));
            %更新系数b
            b1 = b - Ei - label(i)*(alphas(i)-alphasI_old)*...
                Kernel(data(i,:),data(i,:),sigma) - label(j)*...
                (alphas(j)-alphasJ_old)*Kernel(data(i,:),data(j,:),sigma);
            b2 = b - Ej - label(i)*(alphas(i)-alphasI_old)*...
                Kernel(data(i,:),data(j,:),sigma)- label(j)*...
                (alphas(j)-alphasJ_old)*Kernel(data(j,:),data(j,:),sigma);
            %b的几种选择机制
            if alphas(i)>0 && alphas(i)<C
                b = b1;
            elseif alphas(j)>0 && alphas(j)<C
                b = b2;
            else
                b = (b1+b2)/2;
            end
            %确定更新了，记录一次
            alpha_change = alpha_change + 1;
        end
    end
    % 没有实行alpha交换，迭代加1
    if alpha_change == 0
        iter = iter + 1;
    %实行了交换，迭代清0
    else
        iter = 0;
    end
    disp(['iter ================== ',num2str(iter)]);
end
%% 计算权值W
% W = (alphas.*label)'*data;
%记录支持向量位置
index_sup = find(alphas ~= 0);
%计算预測结果
predict = (alphas.*label)'*Kernel(data,data,sigma) + b;
predict = sign(predict);
%% 显示结果
figure;
index1 = find(predict==-1);
data1 = (data(index1,:))';
plot(data1(1,:),data1(2,:),'+r');
hold on
index2 = find(predict==1);
data2 = (data(index2,:))';
plot(data2(1,:),data2(2,:),'*');
hold on
dataw = (data(index_sup,:))';
plot(dataw(1,:),dataw(2,:),'og','LineWidth',2);
title(['核函数參数sigma = ',num2str(sigma)]);