非线性分类问题

• 非线性分类问题是指通过利用非线性模型才能很好地进行分类的问题。对于下面左图，是用一条椭圆曲线（非线性模型）将正负样例正确分类的。然而非线性问题往往不好求解，所以希望能用解线性分类的方法求解原来的非线性问题。如下面右图所示，通过变换将左图中椭圆曲线变换成直线，正确将样例分类。
• 上面的例子说明，用线性分类方法求解非线性分类问题分为两步：首先使用一个变换将原区间的数据映射到新空间，然后在新空间里用线性分类方法从训练数据中学习分类模型。所以，非线性分类的关键就在于如何进行这样的映射。下面即将要介绍的核方法就可以完成这样的映射。

利用核函数将数据映射到高维空间

• 由上图可知，我们将数据由一个特征空间转换到另一个特征空间。在新空间下，我们可以很容易利用已有的工具对数据进行处理。数学家喜欢将这个过程称之为从一个特征空间到另一个特征空间的映射。在通常情况下，这种映射会将低维特征空间映射到高维空间。
• 这种从某个特征空间到另一个特征空间的映射是通过核函数来实现的。核函数可以把数据从某个很难处理的形式转换成为另一个较容易处理的形式。
• 经过空间转换之后，我们可以在高维空间中解决线性问题，这也就等价于在低维空间中解决非线性问题。
• SVM 优化中一个特别好的地方就是，所有的运算都可以写成内积(也称点积)的形式。向量的内积指的是两个向量相乘，之后得到单个标量或者数值。我们可以把内积运算替换成核函数，而不必做简化处理。将内积替换成核函数的方式被称为核技巧 (kemel trick) 或者核 “ 变电”(kemel substation)。

径向基核函数

• 径向基核函数是 SVM 中常用的一个核函数。
  

Code

数据集

• 在文本文件 testSetRBF.txt 中，共计 100 行 3 列。形式如下：
• 图像显示为:

import numpy as np
import matplotlib.pyplot as plt

# 加载数据
def loadDataSet(filename):
    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 plotFigure():
    x, y = loadDataSet('testSetRBF.txt')
    xarr = np.array(x)
    n = np.shape(x)[0]
    x1 = []; y1 = []
    x2 = []; y2 = []
    for i in np.arange(n):
        if int(y[i]) == 1:
            x1.append(xarr[i,0]); y1.append(xarr[i,1])
        else:
            x2.append(xarr[i,0]); y2.append(xarr[i,1])
    
    plt.scatter(x1, y1, s = 30, c = 'r', marker = 's')
    plt.scatter(x2, y2, s = 30, c = 'g')
    plt.show()

if __name__ == '__main__':
    plotFigure()

核转换函数

def kernelTrans(X, A, kTup):
    m, n = np.shape(X)
    K = np.mat(np.zeros((m, 1)))
    if kTup[0] == 'lin':
        K = X * A.T
    elif kTup[0] == 'rbf':
        for j in range(m):
            deltaRow = X[j,:] - A
            K[j] = deltaRow * deltaRow.T
        # 元素间的除法
        K = np.exp(K / (-1 * kTup[1] ** 2))
    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):
        self.X = dataMatIn
        self.labelMat = classLabels
        self.C = C
        self.tol = toler
        self.m = np.shape(dataMatIn)[0]
        self.alphas = np.mat(np.zeros((self.m, 1)))
        self.b = 0
        # 误差缓存,第一列为是否有效标志位，第二列为实际的Ｅ值
        self.eCache = np.mat(np.zeros((self.m, 2))) 
        self.K = np.mat(np.zeros((self.m, self.m)))
        for i in range(self.m):
            self.K[:, i] = kernelTrans(self.X, self.X[i,:], kTup)

使用核函数时需要对 innerL() 及 calcEk() 函数进行修改

# 完整 Platt SMO 算法中的优化例程
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) ):
            # 用于选择第二个 alpha 或者说内循环的　alpha 值
            j, Ej = selectJ(i, oS, Ei)
            alphaIoId = oS.alphas[i].copy()
            alphaJoId = 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]
            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) # 更新误差缓存
            
            if (abs(oS.alphas[j] - alphaJoId) < 0.00001):
                # print('j not moving enough')
                return 0
            oS.alphas[i] += oS.labelMat[j] * oS.labelMat[i] * (alphaJoId - oS.alphas[j])
            updateEk(oS, i) # 更新误差缓存
            
            b1 = oS.b - Ei - oS.labelMat[i] * (oS.alphas[i] - alphaIoId) * \
                oS.K[i, i] - oS.labelMat[j] * \
                (oS.alphas[j] - alphaJoId) * oS.K[i, j]
            b2 = oS.b - Ej - oS.labelMat[i] * (oS.alphas[i] - alphaIoId) * \
                oS.K[i, j] - oS.labelMat[j] * \
                (oS.alphas[j] - alphaJoId) * oS.K[j, j]
                
            if (0 < oS.alphas[i]) and (oS.alphas[i] < oS.C):
                oS.b = b1
            elif (0 < oS.alphas[j]) and (oS.alphas[j] < oS.C):
                oS.b = b2
            else:
                oS.b = (b1 + b2) / 2.0
            return 1
    else:
        return 0

# 完整版 SMO 算法中的外循环代码
def smoP(dataMatIn, classLabels, C, toler, maxIter, kTup = ('lin', 0)):
    oS = optStruct(np.mat(dataMatIn), np.mat(classLabels).transpose(), C, toler, kTup)
    iter = 0
    entireSet = True
    alphaPairsChanged = 0
    while(iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):
        alphaPairsChanged = 0
        if entireSet:
            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:
            # 遍历非边界值
            nonBoundIs = np.nonzero((0 < oS.alphas.A) * (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
        elif (alphaPairsChanged == 0):
            entireSet = True
        # print('iteration number: %d' % iter)
    return oS.b, oS.alphas

利用核函数进行分类的径向基测试函数

• 第一个 for 和第二个 for 仅仅数据集不同，前者训练集后者测试集

def testRbf(k1 = 1.3):
    dataArr, labelArr = loadDataSet('testSetRBF.txt')
    b, alphas =  smoP(dataArr, labelArr, 200, 0.0001, 10000, ('rbf', k1))
    dataMat = np.mat(dataArr)
    labelMat = np.mat(labelArr).transpose()
    svInd = np.nonzero(alphas.A > 0)[0]
    sVs = dataMat[svInd]
    labelSV = labelMat[svInd]
    print('there are %d Support Vectors' % np.shape(sVs)[0])
    m, n = np.shape(dataMat)
    errorCount = 0
    for i in range(m):
        kernelEval = kernelTrans(sVs, dataMat[i,:], ('rbf', k1))
        predict = kernelEval.T * np.multiply(labelSV, alphas[svInd]) + b
        if np.sign(predict) != np.sign(labelArr[i]):
            errorCount += 1
    print('the training error rate is: %f' %(float(errorCount) / m))
    
    dataArr, labelArr = loadDataSet('testSetRBF2.txt')
    errorCount = 0
    dataMat = np.mat(dataArr)
    labelMat = np.mat(labelArr).transpose()
    m, n = np.shape(dataMat)
    for i in range(m):
        kernelEval = kernelTrans(sVs, dataMat[i,:], ('rbf', k1))
        predict = kernelEval.T * np.multiply(labelSV, alphas[svInd]) + b
        if np.sign(predict) != np.sign(labelArr[i]):
            errorCount += 1
    print('the best error rate is: %f' %(float(errorCount) / m))

主函数

if __name__ == '__main__':
    testRbf()

参考

• GitHub源码