Svm算法原理及实现

      Svm（support Vector Mac）又称为支持向量机，是一种二分类的模型。固然若是进行修改以后也是能够用于多类别问题的分类。支持向量机能够分为线性核非线性两大类。其主要思想为找到空间中的一个更够将全部数据样本划开的超平面，而且使得本本集中全部数据到这个超平面的距离最短。

1、基于最大间隔分隔数据

1.1支持向量与超平面

    在了解svm算法以前，咱们首先须要了解一下线性分类器这个概念。好比给定一系列的数据样本，每一个样本都有对应的一个标签。为了使得描述更加直观，咱们采用二维平面进行解释，高维空间原理也是同样。举个简单子：以下图所示是一个二维平面，平面上有两类不一样的数据，分别用圆圈和方块表示。咱们能够很简单地找到一条直线使得两类数据正好可以彻底分开。可是能将据点彻底划开直线不止一条，那么在如此众多的直线中咱们应该选择哪一条呢？从直观感受上看图中的几条直线，是否是要更好一些呢？是的咱们就是但愿寻找到这样的直线，使得距离这条直线最近的点到这条直线的距离最短。这读起来有些拗口，咱们从图三中直观来解释这一句话就是要求的两条外面的线之间的间隔最大。这是能够理解的，由于假如数据样本是随机出现的，那么这样分割以后数据点落入到其类别一侧的几率越高那么最终预测的准确率也会越高。在高维空间中这样的直线称之为超平面，由于当维数大于三的时候咱们已经没法想象出这个平面的具体样子。那些距离这个超平面最近的点就是所谓支持向量，实际上若是肯定了支持向量也就肯定了这个超平面，找到这些支持向量以后其余样本就不会起做用了。

                                  图 1                                                 图2

 

1.2寻找最大间隔

1.2.1点到超平面的距离公式

      既然这样的直线是存在的，那么咱们怎样寻找出这样的直线呢？与二维空间相似，超平面的方程也能够写成一下形式：

                                                                                             （1.1）

有了超平面的表达式以后以后，咱们就能够计算样本点到平面的距离了。假设为样本的中的一个点，其中表示为第个特征变量。那么该点到超平面的距离就能够用以下公式进行计算：

                                                                         （1.2）

 

其中||W||为超平面的范数，常数b相似于直线方程中的截距。

上面的公式能够利用解析几何或高中平面几何知识进行推导，这里不作进一步解释。

1.2.2最大间隔的优化模型

    如今咱们已经知道了如何去求数据点到超平面的距离，在超平面肯定的状况下，咱们就可以找出全部支持向量，而后计算出间隔margin。每个超平面都对应着一个margin，咱们的目标就是找出全部margin中最大的那个值对应的超平面。所以用数学语言描述就是肯定w、b使得margin最大。这是一个优化问题其目标函数能够写成：

                                         （1.3）

其中表示数据点的标签，且其为-1或1。距离用计算，这是就能体会出-1和1的好处了。若是数据点在平面的正方向(即+1类)那么是一个正数，而当数据点在平面的负方向时(即-1类)，依然是一个正数，这样就可以保证始终大于零了。注意到当w和b等比例放大时，d的结果是不会改变的。所以咱们能够令全部支持向量的u为1，而其余点的u大1这是能够办经过调节w和b求到的。所以上面的问题能够简化为：                     （1.4）

为了后面计算的方便，咱们将目标函数等价替换为：

                                                               （1.5）

这是一个有约束条件的优化问题，一般咱们能够用拉格朗日乘子法来求解。拉格朗日乘子法的介绍能够参考这篇博客。应用拉格朗日乘子法以下：

令                                （1.6）

求L关于求偏导数得：                          （1.7）

将(1.7)代入到(1.6)中化简得：

                                      （1.8）

原问题的对偶问题为：

                                              （1.9）

该对偶问题的KKT条件为

                                （1.10）

    到此，彷佛问题就可以完美地解决了。可是这里有个假设：数据必须是百分之百可分的。可是实际中的数据几乎都不那么“干净”，或多或少都会存在一些噪点。为此下面咱们将引入了松弛变量来解决这种问题。

1.2.3松弛变量

    由上一节的分析咱们知道实际中不少样本数据都不可以用一个超平面把数据彻底分开。若是数据集中存在噪点的话，那么在求超平的时候就会出现很大问题。从图三中课看出其中一个蓝点误差太大，若是把它做为支持向量的话所求出来的margin就会比不算入它时要小得多。更糟糕的状况是若是这个蓝点落在了红点之间那么就找不出超平面了。

                                   

                                                         图 3

所以引入一个松弛变量ξ来容许一些数据能够处于分隔面错误的一侧。这时新的约束条件变为:

  （1.11）

式中ξi的含义为容许第i个数据点容许偏离的间隔。若是让ξ任意大的话，那么任意的超平面都是符合条件的了。因此在原有目标的基础之上，咱们也尽量的让ξ的总量也尽量地小。因此新的目标函数变为：

（1.12）

（1.13）

其中的C是用于控制“最大化间隔”和“保证大部分的点的函数间隔都小于1”这两个目标的权重。将上述模型完整的写下来就是：

（1.14）

新的拉格朗日函数变为：

（1.15）

接下来将拉格朗日函数转化为其对偶函数，首先对分别求ξ的偏导，并令其为0,结果以下：

（1.16）

代入原式化简以后获得和原来同样的目标函数：

（1.17）

可是因为咱们获得而，所以有因此对偶问题写成：

（1.18）

通过添加松弛变量的方法，咱们如今可以解决数据更加混乱的问题。经过修改参数C，咱们能够获得不一样的结果而C的大小到底取多少比较合适，须要根据实际问题进行调节。

1.2.4核函数

    以上讨论的都是在线性可分状况进行讨论的，可是实际问题中给出的数据并非都是线性可分的，好比有些数据多是如图4样子。      

                                                          图4

那么这种非线性可分的数据是否就不能用svm算法来求解呢？答案是否认的。事实上，对于低维平面内不可分的数据，放在一个高维空间中去就有可能变得可分。以二维平面的数据为例，咱们能够经过找到一个映射将二维平面的点放到三维平面之中。理论上任意的数据样本都可以找到一个合适的映射使得这些在低维空间不能划分的样本到高维空间中以后可以线性可分。咱们再来看一下以前的目标函数：

（1.19）

定义一个映射使得将全部映射到更高维空间以后等价于求解上述问题的对偶问题：

（1.20）

这样对于线性不可分的问题就解决了，如今只须要找出一个合适的映射便可。当特征变量很是多的时候在，高维空间中计算内积的运算量是很是庞大的。考虑到咱们的目的并非为找到这样一个映射而是为了计算其在高维空间的内积，所以若是咱们可以找到计算高维空间下内积的公式，那么就可以避免这样庞大的计算量，咱们的问题也就解决了。实际上这就是咱们要找的核函数，即两个向量在隐式映射后的空间中的内积。下面的一个简单例子能够帮助咱们更好地理解核函数。

经过以上例子，咱们能够很明显地看到核函数是怎样运做的。上述问题的对偶问题能够写成以下形式：

（1.21）

那么怎样的函数才能够做为核函数呢？下面的一个定理能够帮助咱们判断。

Mercer定理：任何半正定的函数均可以做为核函数。其中所谓半正定函数是指拥有训练集数据集合，咱们定义一个矩阵的元素，这个矩阵是的矩阵，若是这个矩阵是半正定的，那么就称为半正定函数。

值得注意的是，上述定理中所给出的条件是充分条件而非充要条件。由于有些非正定函数也能够做为核函数。

下面是一些经常使用的核函数：

                                             表1 经常使用核函数表

核函数名称	核函数表达式	核函数名称	核函数表达式
线性核		指数核
多项式核		拉普拉斯核
高斯核		Sigmoid核

    如今咱们已经了解了一些支持向量机的理论基础，咱们经过对偶问题的的转化将最开始求的问题转化为求的对偶问题。只要找到全部的(即找出全部支持向量)，咱们就可以肯定。而后就能够经过计算数据点到这个超平面的距离从而判断出该数据点的类别。

2、Smo算法原理

2.1 约束条件

     根据以上问题的分析，咱们已经将原始问题转化为了求的值，即求下面优化模型的解：

（2.1）

求解的值的方法有不少，Smo算法就是其中一种比较经常使用的方法。该算法是由John Platt在1996年发布，他的思路是将大的优化问题转化为小的优化问题。而这些小的优化问题每每更容易求解，而且对它们进行顺序求解的结果和将它们做为总体求解的结果彻底一致可是Smo算法的时间要小得多。

     Smo算法的原理为：每次任意抽取两个乘子和，而后固定和之外的其它乘子，使得目标函数只是关于和的函数。而后增大其中一个乘子同时减小另一个。这样，不断的从一堆乘子中任意抽取两个求解，不断的迭代求解子问题，最终达到求解原问题的目的。

     而原对偶问题的子问题的目标函数能够表达成：

（2.2）

其中：

（2.3）

这里之因此算两个是由于的限制，若是只改变其中的一个量，那么这个约束条件可能就不成立了。要解决这个问题，咱们必须得选取这样的两个乘子。那么怎样肯定这样的和呢？这是咱们首先要考虑的问题，在《机器学习实战》这本书中，做者首先给出了一种简化版的方法，遍历每个而后在剩余的中随机选取一个进行优化。虽然样也可以解决问题，可是运算量太大，所以考虑找一种更好的方法来寻找对。

为了表述方便，定义一个特征到输出结果的输出函数：

（2.4）

该对偶问题中KKT条件为：

（2.5）

根据上述问题的KKT条件能够得出目标函数中的的含义以下：

一、 ，代表是正常分类，在边界外；

二、，代表是支持向量，在边界上；

三、，代表在两边界之间。

最优解须要知足KKT条件，所以须要知足以上的三个条件都知足。而不知足这三个条件的状况也有三种：

一、<=1可是<C则是不知足的，而本来=C;

二、>=1可是>0则是不知足的，而本来=0;

三、=1可是=0或者=C则代表不知足的，而本来应该是0<<C.

也就是说若是存在不知足这些KKT条件的，咱们就要更新它，这就是约束条件之一。其次，还受到约束条件的限制，所以假设选择的两个因子为和，他们在更新前分别为在更新后为，为了保证上述约束条件成立必需要保证下列等式的成立：

（2.6）

其中为常数。

2.2参数优化

    由于两个因子很差同时求解，因此可先求第二个乘子的解，而后再用的解表示的解。为了求解，得先肯定的取值范围。假设它的上下边界分别为H和L，那么有：（2.6）

接下来，综合和这两个约束条件，求取的取值范围。

当时，根据可得，因此有。

当时，一样根据可得：，因此有。

回顾第二个约束条件 ：（2.7）

令其两边同时乘y1，可得：

                          . （2.8）

其中：.

所以能够用表示，即：               （2.9）

令                    （2.10）

通过转化以后可得：

                   ,.        (2.11)

那么如何来选择乘子呢？对于第一个乘子，咱们能够按照3种不知足KTT的条件来寻找。对于第二个乘子，咱们能够寻找知足条件的乘子。

而b在知足如下条件时须要更新：

                       (2.12)

且更新后的和以下：

                      (2.13）

每次更新完两个乘子以后，都须要从新计算b以及对应的E。最后更新完全部的,y和b,这样模型也就出来了，从而能够计算出开始所说的分类函数：

（2.14）

3、Svm的python实现与应用

     第二节已经对Smo算法进行了充分的说明并进行了推导，如今一切都准备好了。接下来须要作的就是实现这些算法了，这里我参考了《机器学习实战》这本书中的代码，并利用该程序对一个问题进行了求解。因为代码数量过大，所以这里就再也不列出，而是放在附录中。有兴趣的朋友能够直接下载，也能够去官网下载源代码。笔者在读这些代码的过程当中，也遇到了许多困难，大部分都根据本身的状况进行了注释。

3.1问题描述

    这里我选取的一个数据集为声呐数据集，该问题为须要根据声呐从不一样角度返回的声音强度来预测被测物体是岩石仍是矿井。数据集中共有208个数据，60个输入变量和1个输出变量。每一个数据的前60个元素为不一样方向返回的声音强度，最后一个元素为本次用于声呐测试的物体类型。其中标签M表示矿井，标签R为岩石。

3.2数据预处理

    所给数据中没有缺失值和异常值，所以不须要对数据集进行清洗。注意到所给数据集的标签为字母类型，而svm算法的标准标签为“-1”和“1”的形式，因此须要对标签进行转化，用“-1”、“1”分别代替M和R。因为该数据集中所给标签相同的数据都放在了一块儿，所以先对数据顺序进行了打乱。代码以下：

def loadDataSet(fileName):    dataMat=[];labelMat=[];data=[]    fr=open(fileName)    for line in fr.readlines():        line=line.replace(',','\t')#去除逗号        line=line.replace('M','-1')#对标签进行替换        line=line.replace('R','1')        lineArr=line.strip('\n').split('\t')#分割数据        data.append([float(lineArr[i]) for i in range(len(lineArr))])    random.shuffle(data)  #随机打乱列表    for i in range(len(data)):        dataMat.append(data[i][0:len(data)-1])        labelMat.append(data[i][-1])    return dataMat,labelMat

3.3模型训练及测试

    首先测试了一下将全部数据即做为训练集又做为测试集，而后用Smo模型进行训练找到全部的支持向量。最后根据训练好的模型进行求解，最终测试的准确率平均为54%。若是把数据集分红测试集和训练集，发现测试的准确率也在这附近。而根据网上数据统计该数据集测试的准确率最高为84%，通常状况下为百分之六十几。数据集自己是形成测试准确率低的一个缘由，可是另一个更加剧要的缘由大概是参数的选择问题。如何选取合适的参数是一个值得思考的问题，在接下来的学习过程当中我也会注意一下参数选取这个问题。到这里，关于svm的算法就告一段落了。

#svm算法的实现
from numpy import*
import random
from time import*
def loadDataSet(fileName):#输出dataArr(m*n),labelArr(1*m)其中m为数据集的个数
    dataMat=[];labelMat=[]
    fr=open(fileName)
    for line in fr.readlines():
        lineArr=line.strip().split('\t')#去除制表符，将数据分开
        dataMat.append([float(lineArr[0]),float(lineArr[1])])#数组矩阵
        labelMat.append(float(lineArr[2]))#标签
    return dataMat,labelMat

def selectJrand(i,m):#随机找一个和i不一样的j
    j=i
    while(j==i):
        j=int(random.uniform(0,m))
    return j

def clipAlpha(aj,H,L):#调整大于H或小于L的alpha的值
    if aj>H:
        aj=H
    if aj<L:
        aj=L
    return aj

def smoSimple(dataMatIn,classLabels,C,toler,maxIter):
    dataMatrix=mat(dataMatIn);labelMat=mat(classLabels).transpose()#转置
    b=0;m,n=shape(dataMatrix)#m为输入数据的个数，n为输入向量的维数
    alpha=mat(zeros((m,1)))#初始化参数，肯定m个alpha
    iter=0#用于计算迭代次数
    while (iter<maxIter):#当迭代次数小于最大迭代次数时（外循环）
        alphaPairsChanged=0#初始化alpha的改变量为0
        for i in range(m):#内循环
            fXi=float(multiply(alpha,labelMat).T*\
                      (dataMatrix*dataMatrix[i,:].T))+b#计算f(xi)
            Ei=fXi-float(labelMat[i])#计算f(xi)与标签之间的偏差
            if ((labelMat[i]*Ei<-toler)and(alpha[i]<C))or\
                    ((labelMat[i]*Ei>toler)and(alpha[i]>0)):#若是能够进行优化
                j=selectJrand(i,m)#随机选择一个j与i配对
                fXj=float(multiply(alpha,labelMat).T*\
                          (dataMatrix*dataMatrix[j,:].T))+b#计算f(xj)
                Ej=fXj-float(labelMat[j])#计算j的偏差
                alphaIold=alpha[i].copy()#保存原来的alpha(i)
                alphaJold=alpha[j].copy()
                if(labelMat[i]!=labelMat[j]):#保证alpha在0到c之间
                    L=max(0,alpha[j]-alpha[i])
                    H=min(C,C+alpha[j]-alpha[i])
                else:
                    L=max(0,alpha[j]+alpha[i]-C)
                    H=min(C,alpha[j]+alpha[i])
                if L==H:print('L=H');continue
                eta=2*dataMatrix[i,:]*dataMatrix[j,:].T-\
                    dataMatrix[i,:]*dataMatrix[i,:].T-\
                    dataMatrix[j,:]*dataMatrix[j,:].T
                if eta>=0:print('eta=0');continue
                alpha[j]-=labelMat[j]*(Ei-Ej)/eta
                alpha[j]=clipAlpha(alpha[j],H,L)#调整大于H或小于L的alpha
                if (abs(alpha[j]-alphaJold)<0.0001):
                    print('j not move enough');continue
                alpha[i]+=labelMat[j]*labelMat[i]*(alphaJold-alpha[j])
                b1=b-Ei-labelMat[i]*(alpha[i]-alphaIold)*\
                    dataMatrix[i,:]*dataMatrix[i,:].T-\
                    labelMat[j]*(alpha[j]-alphaJold)*\
                    dataMatrix[i,:]*dataMatrix[j,:].T#设置b
                b2=b-Ej-labelMat[i]*(alpha[i]-alphaIold)*\
                    dataMatrix[i,:]*dataMatrix[i,:].T-\
                    labelMat[j]*(alpha[j]-alphaJold)*\
                    dataMatrix[j,:]*dataMatrix[j,:].T
                if (0<alpha[i])and(C>alpha[j]):b=b1
                elif(0<alpha[j])and(C>alpha[j]):b=b2
                else:b=(b1+b2)/2
                alphaPairsChanged+=1
                print('iter:%d i:%d,pairs changed%d'%(iter,i,alphaPairsChanged))
        if (alphaPairsChanged==0):iter+=1
        else:iter=0
        print('iteraction number:%d'%iter)
    return b,alpha
#定义径向基函数
def kernelTrans(X, A, kTup):#定义核转换函数（径向基函数）
    m,n = shape(X)
    K = mat(zeros((m,1)))
    if kTup[0]=='lin': K = X * A.T   #线性核K为m*1的矩阵
    elif kTup[0]=='rbf':
        for j in range(m):
            deltaRow = X[j,:] - A
            K[j] = deltaRow*deltaRow.T
        K = exp(K/(-1*kTup[1]**2)) #divide in NumPy is element-wise not matrix like Matlab
    else: raise NameError('Houston We Have a Problem -- \
    That Kernel is not recognized')
    return K

class optStruct:
    def __init__(self,dataMatIn, classLabels, C, toler, kTup):  # Initialize the structure with the parameters
        self.X = dataMatIn
        self.labelMat = classLabels
        self.C = C
        self.tol = toler
        self.m = shape(dataMatIn)[0]
        self.alphas = mat(zeros((self.m,1)))
        self.b = 0
        self.eCache = mat(zeros((self.m,2))) #first column is valid flag
        self.K = mat(zeros((self.m,self.m)))
        for i in range(self.m):
            self.K[:,i] = kernelTrans(self.X, self.X[i,:], kTup)

def calcEk(oS, k):
    fXk = float(multiply(oS.alphas,oS.labelMat).T*oS.K[:,k] + oS.b)
    Ek = fXk - float(oS.labelMat[k])
    return Ek

def selectJ(i, oS, Ei):
    maxK = -1; maxDeltaE = 0; Ej = 0
    oS.eCache[i] = [1,Ei]
    validEcacheList = nonzero(oS.eCache[:,0].A)[0]
    if (len(validEcacheList)) > 1:
        for k in validEcacheList:   #loop through valid Ecache values and find the one that maximizes delta E
            if k == i: continue #don't calc for i, waste of time
            Ek = calcEk(oS, k)
            deltaE = abs(Ei - Ek)
            if (deltaE > maxDeltaE):
                maxK = k; maxDeltaE = deltaE; Ej = Ek
        return maxK, Ej
    else:   #in this case (first time around) we don't have any valid eCache values
        j = selectJrand(i, oS.m)
        Ej = calcEk(oS, j)
    return j, Ej

def updateEk(oS, k):#after any alpha has changed update the new value in the cache
    Ek = calcEk(oS, k)
    oS.eCache[k] = [1,Ek]

def innerL(i, oS):
    Ei = calcEk(oS, i)
    if ((oS.labelMat[i]*Ei < -oS.tol) and (oS.alphas[i] < oS.C)) or ((oS.labelMat[i]*Ei > oS.tol) and (oS.alphas[i] > 0)):
        j,Ej = selectJ(i, oS, Ei) #this has been changed from selectJrand
        alphaIold = oS.alphas[i].copy(); alphaJold = oS.alphas[j].copy()
        if (oS.labelMat[i] != oS.labelMat[j]):
            L = max(0, oS.alphas[j] - oS.alphas[i])
            H = min(oS.C, oS.C + oS.alphas[j] - oS.alphas[i])
        else:
            L = max(0, oS.alphas[j] + oS.alphas[i] - oS.C)
            H = min(oS.C, oS.alphas[j] + oS.alphas[i])
        if L==H: print("L==H"); return 0
        eta = 2.0 * oS.K[i,j] - oS.K[i,i] - oS.K[j,j] #changed for kernel
        if eta >= 0: print("eta>=0"); return 0
        oS.alphas[j] -= oS.labelMat[j]*(Ei - Ej)/eta
        oS.alphas[j] = clipAlpha(oS.alphas[j],H,L)
        updateEk(oS, j) #added this for the Ecache
        if (abs(oS.alphas[j] - alphaJold) < 0.00001): print("j not moving enough"); return 0
        oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])#update i by the same amount as j
        updateEk(oS, i) #added this for the Ecache                    #the update is in the oppostie direction
        b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,i] - oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[i,j]
        b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,j]- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[j,j]
        if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1
        elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2
        else: oS.b = (b1 + b2)/2.0
        return 1
    else: return 0
#smoP函数用于计算超平的alpha,b
def smoP(dataMatIn, classLabels, C, toler, maxIter,kTup=('lin', 0)):    #完整的Platter SMO
    oS = optStruct(mat(dataMatIn),mat(classLabels).transpose(),C,toler, kTup)
    iter = 0#计算循环的次数
    entireSet = True; alphaPairsChanged = 0
    while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):
        alphaPairsChanged = 0
        if entireSet:   #go over all
            for i in range(oS.m):
                alphaPairsChanged += innerL(i,oS)
                print("fullSet, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged))
            iter += 1
        else:#go over non-bound (railed) alphas
            nonBoundIs = nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
            for i in nonBoundIs:
                alphaPairsChanged += innerL(i,oS)
                print("non-bound, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged))
            iter += 1
        if entireSet: entireSet = False #toggle entire set loop
        elif (alphaPairsChanged == 0): entireSet = True
        print("iteration number: %d" % iter)
    return oS.b,oS.alphas
#calcWs用于计算权重值w
def calcWs(alphas,dataArr,classLabels):#计算权重W
    X = mat(dataArr); labelMat = mat(classLabels).transpose()
    m,n = shape(X)
    w = zeros((n,1))
    for i in range(m):
        w += multiply(alphas[i]*labelMat[i],X[i,:].T)
    return w

#值得注意的是测试准确与k1和C的取值有关。
def testRbf(k1=1.3):#给定输入参数K1
    #测试训练集上的准确率
    dataArr,labelArr = loadDataSet('testSetRBF.txt')#导入数据做为训练集
    b,alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, ('rbf', k1)) #C=200 important
    datMat=mat(dataArr); labelMat = mat(labelArr).transpose()
    svInd=nonzero(alphas.A>0)[0]#找出alphas中大于0的元素的位置
    #此处须要说明一下alphas.A的含义
    sVs=datMat[svInd] #获取支持向量的矩阵，由于只要alpha中不等于0的元素都是支持向量
    labelSV = labelMat[svInd]#支持向量的标签
    print("there are %d Support Vectors" % shape(sVs)[0])#输出有多少个支持向量
    m,n = shape(datMat)#数据组的矩阵形状表示为有m个数据，数据维数为n
    errorCount = 0#计算错误的个数
    for i in range(m):#开始分类，是函数的核心
        kernelEval = kernelTrans(sVs,datMat[i,:],('rbf', k1))#计算原数据集中各元素的核值
        predict=kernelEval.T * multiply(labelSV,alphas[svInd]) + b#计算预测结果y的值
        if sign(predict)!=sign(labelArr[i]): errorCount += 1#利用符号判断类别
        ### sign（a）为符号函数：若a>0则输出1，若a<0则输出-1.###
    print("the training error rate is: %f" % (float(errorCount)/m))
    #二、测试测试集上的准确率
    dataArr,labelArr = loadDataSet('testSetRBF2.txt')
    errorCount = 0
    datMat=mat(dataArr)#labelMat = mat(labelArr).transpose()此处能够不用
    m,n = shape(datMat)
    for i in range(m):
        kernelEval = kernelTrans(sVs,datMat[i,:],('rbf', k1))
        predict=kernelEval.T * multiply(labelSV,alphas[svInd]) + b
        if sign(predict)!=sign(labelArr[i]): errorCount += 1
    print("the test error rate is: %f" % (float(errorCount)/m))
def main():
    t1=time()
    dataArr,labelArr=loadDataSet('testSet.txt')
    b,alphas=smoP(dataArr,labelArr,0.6,0.01,40)
    ws=calcWs(alphas,dataArr,labelArr)
    testRbf()
    t2=time()
    print("程序所用时间为%ss"%(t2-t1))

if __name__=='__main__':
    main()

 

                                                   后记

    这是第一次写博客，其中不免会出错，所以但愿你们可以批评指正。首先很是感谢网上的一些朋友，在理解svm这算法他们给了我不少启发，在公式推导中给了我不少参考的地方。本文主要参考的资料是《机器学习实战》和《惊呼！理解svm的三种境界》这篇博客。对于svm虽然学的时间不长，可是对它一直有种特别的感受。第一次据说svm是在作一个验证码识别问题的时候，但那时候我使用的是KNN算法，尽管效果还不错，可是我一直但愿可以用svm算法来完成这个题目。原本此次是打算把这个问题一块儿解决的，可是因为时间关系，没有来得及作。只能等下次有空闲的时候再来作这个问题了。