您的位置:首页 > 其它

[完]机器学习实战 第六章 支持向量机(Support Vector Machine)

2016-09-06 23:03 579 查看
[参考] 机器学习实战(Machine Learning in Action)

本章内容

支持向量机(Support Vector Machine)是最好的现成的分类器,“现成”指的是分类器不加修改即可直接使用。基本形式的SVM分类器就可得到低错误率的结果。SVM有很多实现,文中采用最流行的一种实现,即序列最小优化Sequential Minimal Optimization,一种求解支持向量机二次规划的算法)算法。还会介绍如何使用一种称为核函数kernel)的方式将SVM扩展到更多的数据集上。

优点:泛化错误率低,计算开销不大,结果易理解。

这些优点是其十分流行,有人认为他是监督学习中最好的定式算法。与下章中AdaBoost最好的监督学习的方法相对应。

缺点:对参数调节和核函数的选择敏感,原始分类器不加修改仅适用于二类问题。

使用数据类型:数值型和标称型数据。

希望找到离分隔超平面最近的点,确保它们离分隔面的距离尽可能远。这里点到分隔面的距离被称为间隔margin)。支持向量(support vector)就是离分隔超平面最近的那些点

分隔超平面的形式可以写成wT+b,要计算点A到分隔超平面的距离,就必须给出点到分隔面的法线或者垂线的长度,该值为|wTA+b|/∥w∥。常数b类似于Logistic回归中的截距w0。向量w和常数b一起描述所给数据的分割线或超平面。

其他概念:线性可分(linearly separable),分隔超平面(separating hyperplane),间隔(margin)

SVM使用类似海维赛德阶跃函数(即单位阶跃函数)的函数对wT+b作用得到f(wT+b),其中当u<0时f(u)=−1,反之则输出+1。不像Logistic回归,采用0和1,是因为-1和+1仅仅相差一个符号,方便数学上的处理。可以通过一个统一公式来表示间隔或者数据点到分隔超平面的距离,同时不必担心数据到底是属于-1还是+1类。

间隔通过label∗(wTx+b)来表示计算目标就是找出分类器定义中的w和b,需要找到具有最小间隔的点,这些点就是支持向量。然后,对间隔进行最大化。可以写作:

arg maxw,b{minn(label⋅(wTx+b))⋅1∥w∥}

上面是对乘积进行优化的事情,需要固定其中一个因子而最大化其他因子。若令所有支持向量的label∗(wTx+b)都为1,那么可以求∥w∥−1的最大值来得到最终解。但是并非所有数据点的label∗(wTx+b)都等1,只有那些离分隔超平面最近的点得到的值才为1。而离超平面越远的数据点,其label∗(wTx+b)的值也就越大。

上述优化问题给定了一些约束条件然后求最优值,该问题是一个带约束条件的优化问题。这里的优化条件是label∗(wTx+b)≥1.0。此类问题的一个著名的求解方法是拉格朗日乘子法。通过引入拉格朗日乘子,可基于约束条件表达原来的问题。由于这里的约束条件是基于数据点的,因此就可以将超平面写成数据点的形式。优化目标函数最后可以写成:

maxa⎡⎣∑i=1mα−12∑i,j=1mlabel(i)⋅label(i)⋅αi⋅αj⟨xi,xj⟩⎤⎦(1)

其中尖括号表示xi,xj两个向量的内积,整个目标函数的约束条件为:

α≥0,和∑i=1mαj⋅label(i)=0

因为并非所有数据都100%线性可分,所以需要引入松弛变量(slack variable),来允许有些数据点可处于分隔面的错误一侧。这样优化目标可保持不变,但此时约束条件则变为:

C≥α≥0,和∑i=1mαi⋅label(i)=0(2)

常数C用于控制“最大化间隔”和“保证大部分点的间隔小于1.0”这两个目标的权重。在优化算法的实现代码中,C是一个参数,可通过调节该参数得到不同的结果。一旦求出所有的alpha,那么分隔超平面可通过这些alpha来表达。SVM的主要工作就是求解这些alpha。

其中(1)是最小化目标函数,(2)是优化过程中必须遵循的约束条件。以前,采用二次规划求解工具(quadratic solver)求解上述最优化问题。

SMO(Sequential Minimal Optimization)是一个强大的算法,表示序列最小优化。SMO将大优化问题分解成多个小优化问题来求解。此算法目标是求出一系列alpha和b,一旦求出这些alpha,就容易计算出权重向量w。

SMO算法的工作原理:每次循环中选择两个alpha进行优化处理。一旦找到一对“合适”的alpha,那么就增大其中一个同时减小另一个。“合适”指两个alpha必须要符合一定条件,条件之一就是这两个alpha必须要在间隔边界之外,而其第二个条件则是这两个alpha还没有进行过区间化处理或者不在边界上。

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

特征空间转换,称为从一个特征空间到另一个特征空间的映射。一般,这种空间映射会将低维特征空间映射到高维空间。这种映射可通过核函数(kernel)来实现。可以将核函数想象成一个包装器(wrapper)或者接口(interface)。经过空间转换后,可在高维空间中解决线性问题,这等价于在低维空间中解决非线性问题。

SVM优化中,所有运算都可以写成内积(inner product,也成为点积)的形式。向量的内积指的是两个向量相乘,之后得到单个标量或者数值。可将内积运算替换成核函数,而不必做简化处理。将内积替换成核函数的方式成为核技巧(kernel trick)或者核“变电”(kernel substation)

径向基核函数(radial basis function),是一个采用向量作为自变量的函数,能够基于向量距离运算输出一个标量。这个距离可以是从<0,0>向量或者其他向量开始计算的距离。可使用径向基核函数的高斯版本,具体公式:

k(x,y)=exp(−∥x−y∥22σ2)

其中,σ是用户定义的用于确定到达率(reach)或者函数值跌落到0的速度参数。高斯核函数可将数据映射到一个无穷维的空间。

使用函数

函数功能
str.strip(rm)删除str字符串中开头、结尾处,位于rm删除序列的字符,当rm为空时,默认删除空白符(含
'\n','\r','\t',' '
str.lstrip(rm)删除str字符串中开头处,位于rm删除序列的字符
str.rstrip(rm)删除str字符串中结尾处,位于rm删除序列的字符
random.uniform(x, y)在[x,y]范围内,随机生成一个实数
np.multiply(a, b)对应下标元素相乘
np.multiply.outer(a, b)列表a中元素分别和b相乘,得到一个二维数组
list.copy()列表的浅复制
alphas[alphas>0]数组过滤,返回alphas中大于0的元素
mat.A将matrix数据转换成array
mat[ind]获取ind列表中数字对应的mat矩阵的行
listdir(dirName)from os import listdir,获取给定文件夹下的文件名列表,不含文件路径

程序代码

# coding=utf-8

from numpy import *

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 selectJrand(i,m) :
j=i
while(j==i) :
j = int(random.uniform(0,m))
return j

def clipAlpha(aj, H, L) :
if aj > H :
aj = H
if L > aj :
aj = L
return aj

# dataMatIn:        数据集
# classLabels:      类别标签
# C:                常数C
# toler:            容错率
# maxIter:          退出前最的循环次数
def smoSimple(dataMatIn, classLabels, C, toler, maxIter) :
dataMatrix = mat(dataMatIn)
labelMat = mat(classLabels).transpose()
b = 0;
m,n = shape(dataMatrix)
alphas = mat(zeros((m,1)))
iter = 0
while (iter < maxIter) :
# alphaPairsChanged记录alpha是否已经进行优化
alphaPairsChanged = 0
for i in range(m) :
# fXi预测的类别
fXi = float(multiply(alphas, labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + b
# 预测结果与真实结果比对,计算误差Ei
Ei = fXi - float(labelMat[i])
# 如果误差很大,那么可对该数据实例所对应的alpha值进行优化,分别对正间隔和负间隔做了测试,
# 并且检查了alpha值,保证其不能等于0或者C,由于后面alpha小于0或者大于C时将被调整为0或C,
# 所以一旦该if语句中它们等于这两个值得话,那么它们就已经在“边界”上了,因而不再能够减小或增大,
# 因此也就不值得对它们进行优化
if((labelMat[i]*Ei < -toler) and (alphas[i] < C)) or ((labelMat[i]*Ei > toler) and (alphas[i] > 0)) :
# 利用辅助函数,随机选择第二个alpha值
j = selectJrand(i,m)
fXj = float(multiply(alphas, labelMat).T*(dataMatrix*dataMatrix[j,:].T)) + b
Ej = fXj - float(labelMat[j])
alphaIold = alphas[i].copy()
alphaJold = alphas[j].copy()
# L和H用于将alphas[j]调整到0-C之间。如果L==H,就不做任何改变,直接执行continue语句
if labelMat[i] != labelMat[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])
if L==H : print 'L==H'; continue
# eta是alphas[j]的最优修改量,如果eta==0,需要退出for循环的当前迭代过程
eta = 2.0 * dataMatrix[i, :] * dataMatrix[j, :].T - dataMatrix[i, :]*dataMatrix[i, :].T \
- dataMatrix[j, :]*dataMatrix[j, :].T
if eta >= 0 : print 'eta>=0' ; continue
# 计算出一个新的alphas[j]值
alphas[j] -= labelMat[j]*(Ei - Ej)/eta
# 并使用辅助函数,以及L和H对其进行调整
alphas[j] = clipAlpha(alphas[j], H, L)
# 检查alphas[j]是否有轻微改变,如果是的话,退出for循环
if( abs(alphas[j] - alphaJold) < 0.00001) : print 'j not moving enough'; continue
# 然后alphas[i]和alphas[j]同样进行改变,虽然改变的大小一样,但是改变的方向正好相反
alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])
# 在对alphas[i]和alphas[j]进行优化后,给这两个alpha值设置一个常数项b
b1 = b - Ei - labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i, :]*dataMatrix[i, :].T\
- labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[i, :]*dataMatrix[j, :].T
b2 = b - Ej - labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i, :]*dataMatrix[j, :].T\
- labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[j, :]*dataMatrix[j, :].T
if (0 < alphas[i]) and (C > alphas[i]) : b = b1
elif (0 < alphas[j]) and (C > alphas[j]) : b = b2
else: b = (b1 + b2)/2.0
alphaPairsChanged += 1
print 'iter: %d i:%d, pairs changed %d' % (iter, i, alphaPairsChanged)
if alphaPairsChanged == 0 : iter += 1
else: iter = 0
print 'iteration number: %d' % iter
return b, alphas

# 建立一个数据结构来保存所有重要值
'''
class optStruct:
def __init__(self, dataMatIn, classLabels, C, toler) :
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)))
'''

# 使用Kernel函数的class optStruct
class optStruct :
def __init__(self, dataMatIn, classLabels, C, toler, kTup) :
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)))
self.K = mat(zeros((self.m, self.m)))
for i in range(self.m) :
self.K[:,i] = kernelTrans(self.X, self.X[i,:], kTup)

'''
# 对于给定的alpha值,calcEk()能够计算E值并返回
def calcEk(oS, k) :
fXk = float(multiply(oS.alphas, oS.labelMat).T*(oS.X*oS.X[k,:].T)) + oS.b
Ek = fXk - float(oS.labelMat[k])
return Ek
'''
# 使用Kernel函数的calcEk
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

# 内循环中的启发式方法,用于选择第二个alpha或者内循环的alpha值
# 目标:选择合适的第二个alpha值以保证在每次优化中采用最大步长
# 该函数的误差值与第一个alpha值Ei和下标i有关
def selectJ(i, oS, Ei) :
maxK = -1; maxDeltaE = 0; Ej = 0
# 将输入值Ei在缓存中设置成为有效的
oS.eCache[i] = [1, Ei]
# nonzero()返回一个列表,列表包含以输入列表为目录的列表值,这里的值并非零
# 返回的非零E值对应的alpha值,而不是E值本身
validEcacheList = nonzero(oS.eCache[:,0].A)[0]
if len(validEcacheList) > 1 :
for k in validEcacheList :
if k == i : continue
Ek = calcEk(oS, k)
deltaE = abs(Ei - Ek)
# 选择具有最大步长的j
if deltaE > maxDeltaE :
maxK = k; maxDeltaE = deltaE; Ej = Ek
return maxK, Ej
else :# 第一次循环,随机选择一个alpha值
j = selectJrand(i, oS.m)
Ej = calcEk(oS, j)
return j, Ej

# 计算误差值并存入缓存当中
def updateEk(oS, k) :
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]

# 完整的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) ) :
j, Ej = selectJ(i, oS, Ei)
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.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].T
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] - alphaJold) < 0.00001 :
print "j not moving enough"; return 0
oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])
updateEk(oS, i)
b1 = oS.b - Ei - oS.labelMat[i]*(oS.alphas[i] - alphaIold)*oS.X[i,:]*oS.X[i,:].T -\
oS.labelMat[j]*(oS.alphas[j] - alphaJold)*oS.X[i,:]*oS.X[j,:].T
b2 = oS.b - Ej - oS.labelMat[i]*(oS.alphas[i] - alphaIold)*oS.X[i,:]*oS.X[j,:].T -\
oS.labelMat[j]*(oS.alphas[j] - alphaJold)*oS.X[j,:]*oS.X[j,:].T
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
'''

# 使用Kernel函数的innerL
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)
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.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].T
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] - alphaJold) < 0.00001 :
print "j not moving enough"; return 0
oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])
updateEk(oS, i)
'''
b1 = oS.b - Ei - oS.labelMat[i]*(oS.alphas[i] - alphaIold)*oS.X[i,:]*oS.X[i,:].T -\
oS.labelMat[j]*(oS.alphas[j] - alphaJold)*oS.X[i,:]*oS.X[j,:].T
b2 = oS.b - Ej - oS.labelMat[i]*(oS.alphas[i] - alphaIold)*oS.X[i,:]*oS.X[j,:].T -\
oS.labelMat[j]*(oS.alphas[j] - alphaJold)*oS.X[j,:]*oS.X[j,:].T
'''
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

# 完整版Platt SMO算法的外循环代码
def smoP(dataMatIn, classLabels, C, toler, maxIter, kTup=('lin', 0)) :
oS = optStruct(mat(dataMatIn), mat(classLabels).transpose(), C, toler, kTup)
iter = 0
entireSet = True; alphaPairsChanged = 0
# 退出循环条件:1、迭代次数超过指定最大值;2、遍历整个集合都未对任意alpha对进行修改。
while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)) :
alphaPairsChanged = 0
if entireSet :
# 在数据集上遍历任意可能的alpha,使用innerL()来选择第二个alpha,并在可能时对其进行优化
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 = nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
# 遍历所有非边界alpha值,也就是不在边界0或C上的值
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

# 利用alpha值,进行分类
def calcWs(alphas, dataArr, classLabels) :
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

# 核转换函数
# kTup: 元组,核函数的信息,元组的一个参数是描述所用核函数类型的一个字符串,另外2个参数则都是核函数可能需要的可选参数
# 一个调用实例:kernelTrans(sVs, dataMat[i,:], ('rbf', k1))
# 其中k1是径向基核函数高斯版本中的sigma
def kernelTrans(X, A, kTup) :
m,n = shape(X)
# 构建一个列向量
K = mat(zeros((m, 1)))
# 检查元组以确定核函数的类型
if kTup[0] == 'lin' : K = X * A.T
# 在径向基核函数的情况下
elif kTup[0] == 'rbf' :
# for循环中对于矩阵的每个元素计算高斯函数的值
for j in range(m) :
deltaRow = X[j, :] - A
K[j] = deltaRow*deltaRow.T
# 将计算过程应用到整个向量,元素间的除法
K = exp(K/(-1*kTup[1]**2))
else : raise NameError('Houston We Have a Problem -- That Kernel is not recognized')
return K

# 利用核函数进行分类的径向基测试函数
# k1: 高斯径向基函数中一个用户定义变量
# 此函数从文件中读取数据集,然后在该数据集上运行Platt SMO算法,其中核函数的类型是'rbf'
def testRbf(k1 = 1.3) :
dataArr, labelArr = loadDataSet('C:\python27\ml\\testSetRBF.txt')
b, alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, ('rbf', k1))
dataMat = mat(dataArr); labelMat = mat(labelArr).transpose()
svInd = nonzero(alphas.A>0)[0]
sVs = dataMat[svInd]
labelSV = labelMat[svInd]
print "there are %d Support Vectors" % shape(sVs)[0]
m,n = shape(dataMat)
errorCount = 0
for i in range(m) :
# for循环中前两行,给出了如何利用核函数进行分类
kernelEval = kernelTrans(sVs, dataMat[i,:], ('rbf', k1))
predict = kernelEval.T * multiply(labelSV, alphas[svInd]) + b
if sign(predict) != sign(labelArr[i]) : errorCount += 1
print "the training error rate is: %f" % (float(errorCount)/m)
dataArr, labelArr = loadDataSet('C:\python27\ml\\testSetRBF2.txt')
errorCount = 0
dataMat = mat(dataArr); labelMat = mat(labelArr).transpose()
m,n = shape(dataMat)
for i in range(m) :
kernelEval = kernelTrans(sVs, dataMat[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)

# 基于SVM的手写数字识别
def loadImages(dirName) :
from os import listdir
hwLabels = []
trainingFileList = listdir(dirName)
m = len(trainingFileList)
trainingMat = zeros((m, 1024))
for i in range(m) :
fileNameStr = trainingFileList[i]
fileStr = fileNameStr.split('.')[0]
classNumStr = int(fileStr.split('_')[0])
if classNumStr == 9 : hwLabels.append(-1)
else : hwLabels.append(1)
trainingMat[i, :] = img2vector('%s/%s' % (dirName, fileNameStr))
return trainingMat, hwLabels

def testDigits(kTup = ('rbf', 10)) :
dataArr, labelArr = loadImages('c:\python27\ml\\trainingDigits')
b, alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, kTup)
dataMat = mat(dataArr)
labelMat = mat(labelArr).transpose()
svInd = nonzero(alphas.A > 0)[0]
sVs = dataMat[svInd]
labelSV = labelMat[svInd]
print "there are %d Support Vectors" % shape(sVs)[0]
m,n = shape(dataMat)
errorCount = 0
for i in range(m) :
kernelEval = kernelTrans(sVs, dataMat[i, :], kTup)
predict = kernelEval.T * multiply(labelSV, alphas[svInd]) + b
if sign(predict) != sign(labelArr[i]) : errorCount += 1
print "the training error rate is: %f" % (float(errorCount)/m)
dataArr, labelArr = loadImages('c:\python27\ml\\testDigits')
errorCount = 0
dataMat = mat(dataArr)
labelArr = mat(labelArr).transpose()
m,n = shape(dataMat)
for i in range(m) :
kernelEval = kernelTrans(sVs, dataMat[i, :], kTup)
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)


在命令行中执行:

# 加载数据
>>> import ml.svmML as svmML
>>> dataArr, labelArr = svmML.loadDataSet('c:\python27\ml\\testSet.txt')
>>> labelArr
[-1.0, -1.0, 1.0, -1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, 1.0,
-1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, 1.0, -1.0, -1.0, 1.0, 1.0, -1.0, -1
.0, -1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 1.0, -1.0, -1.0, -1.0
, -1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0, -1.0, -1.0, 1.0, 1.0, -1.0, 1.0
, -1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 1.0, 1.0, 1.0, -1.0, 1.0,
1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0,
-1.0, 1.0, -1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0]

# 执行简化版SMO算法
>>> b, alphas = svmML.smoSimple(dataArr, labelArr, 0.6, 0.001, 40)
......
iteration number: 23
j not moving enough
j not moving enough
j not moving enough
iteration number: 24
iter: 24 i:29, pairs changed 1
j not moving enough
j not moving enough
iteration number: 0
j not moving enough
j not moving enough
j not moving enough
iteration number: 1
......
>>> b               # 因为SMO算法的随机性,每次运行的结果不一样
matrix([[-3.84267434]])
>>> alphas[alphas>0]
matrix([[ 0.08891067,  0.27233877,  0.03016793,  0.33108151]])
>>> import numpy as np
>>> np.shape(alphas[alphas>0])              # 得到支持向量的个数
(1, 4)
>>> for i in range(100):                    # 得到哪些点是支持向量
...     if alphas[i]>0.0: print dataArr[i], labelArr[i]
...
[4.658191, 3.507396] -1.0
[3.457096, -0.082216] -1.0
[5.286862, -2.358286] 1.0
[6.080573, 0.418886] 1.0

# 完整版Platt SMO
>>> reload(svmML)
<module 'ml.svmML' from 'C:\Python27\ml\svmML.py'>
>>> dataArr, labelArr = svmML.loadDataSet('c:\python27\ml\\testSet.txt')
>>> b, alphas = svmML.smoP(dataArr, labelArr, 0.6, 0.001, 40)
......
non-bound, iter: 1 i:55, pairs changed 0
non-bound, iter: 1 i:94, pairs changed 0
iteration number: 2
j not moving enough
j not moving enough
j not moving enough
j not moving enough
j not moving enough
j not moving enough
j not moving enough
L==H
j not moving enough
j not moving enough
L==H
j not moving enough
j not moving enough
L==H
L==H
j not moving enough
fullSet, iter: 2 i: 99, pairs changed 0
iteration number: 3
>>> alphas[alphas>0]
matrix([[ 0.06961952,  0.0169055 ,  0.0169055 ,  0.0272699 ,  0.04522972,
0.0272699 ,  0.0243898 ,  0.06140181,  0.06140181]])

# 进行数据点分类
>>> ws = svmML.calcWs(alphas, dataArr, labelArr)
>>> ws
array([[ 0.65307162],
[-0.17196128]])
>>> import numpy as np
>>> dataMat = np.mat(dataArr)
>>> dataMat[0]
matrix([[ 3.542485,  1.977398]])
>>> b
matrix([[-2.89901748]])
>>> dataMat[0]*np.mat(ws) + b
matrix([[-0.92555695]])
>>> labelArr[0]
-1.0
>>> dataMat[2]*np.mat(ws) + b
matrix([[ 2.30436336]])
>>> labelArr[2]
1.0
>>> dataMat[1]*np.mat(ws) + b
matrix([[-1.36706674]])
>>> labelArr[1]
-1.0

# 在测试中使用核函数
>>> reload(svmML)
<module 'ml.svmML' from 'C:\Python27\ml\svmML.pyc'>
>>> svmML.testRbf()
......
j not moving enough
L==H
L==H
L==H
L==H
fullSet, iter: 6 i: 99, pairs changed 0
iteration number: 7
there are 27 Support Vectors
the training error rate is: 0.030000
the test error rate is: 0.040000

# 基于SVM的手写数字识别
>>> reload(svmML)
<module 'ml.svmML' from 'C:\Python27\ml\svmML.py'>
>>> svmML.testDigits(('rbf', 20))
......
j not moving enough
j not moving enough
j not moving enough
fullSet, iter: 7 i: 1933, pairs changed 0
iteration number: 8
there are 157 Support Vectors
the training error rate is: 0.000000
the test error rate is: 0.012685
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: