支持向量机
1.基本思想
2.数学知识补充一
1、超平面的表示方法:
其中,w表示超平面的法向量;b为相对于过原点超平面的偏移。
下面两个公式分别表示超平面两侧的区域
2、点到超平面的距离公式
可以理解,超平面上的点都满足超平面方程,那么超平面外的点到超平面的距离与超平面方程左侧绝对值成正比,同时可知法向量w可以任意倍数缩放为不改变超平面,所以对w要归一化处理。
3、凸优化问题
凸函数定义域是凸集,且函数本身满足Jensen不等式。
凸优化问题可以从两部分描述:优化目标 f(x) 和 约束条件 如g(x) h(x)
凸优化问题的标准数学形式:
3.基本原型和理论推导
将要找的超平面记作 wTx+b =0,因为要将两类数据完美分开,所以对于标签y =+1的和y=-1的样本,
这里只考察可以完美分开并且流出间隔的情况,称为线性可分,线性不可分的情况稍后讨论。
超平面方程的右侧是一个常数,而w和b可以放缩,所以常数大小可以进行设置,为了方便计算姑且设置为1.此时。
可以理解如果存在一个间隔,那么处于间隔边界的点,正好支撑着整个类别的数据集不约过边界,进入被预留出来的距离为y的间隔中,这样的点(也即距离超平面最近的点)就称为支撑向量。而支撑向量使得上面不等式取等号。
d+表示使不等式取+1的支撑向量,间隔width =d- + d+
支撑向量还有一个重要的特点:它是真正决定分界面位置的样本点。
整理前面公式可以把SVM模型写成一个优化问题,即
4.拉格朗日乘子法
等价于
以下为解释。
对于无约束优化问题可以利用求偏导取0,来寻找凸函数的唯一的全局最优值,但是对于有约束优化问题,情况并非如此,此时拉格朗日乘子法是一个处理一般有约束问题的强大工具。
依照拉格朗日函数,可导出关于拉格朗日乘子的对偶优化问题。通常,对偶问题是原问题最优解的一个下界,当满足KKT条件时,对偶问题与原问题的解时等价的。通常,支持向量机所对应的二次优化问题实在对偶形式下求解的。
支持向量机的对偶形式的优化问题传达两点重要信息。
1、表明寻求最大间隔决策面等价于寻找支持向量(即那些非零拉格朗日乘子所对应的训练样本)
2、表明该优化问题由训练实例的逐对内积(即Gram矩阵)所确定。
上面还只是讨论了线性可分的情况,如果数据不是线性可分的呢?
5.数学知识补充二
这就是核方法的基本出发点,即通过讲线性不可分的数据映射到高位空间,从而使其线性可分,甚至不需要显式的给出这个映射表达。
7.核方法
回到问题
首先考虑样本点从低维到高维映射,
此时原问题就变为
而Φ(x)是一个列向量,Φ(x) T Φ(x)就是向量的内积。
预设了Φ(x)是高维度的,甚至无穷维,如果每一个x都做一个到Φ(x)的映射,在两两计算,计算两会非常大,直接凉凉。能不能有一个二元函数值刚好是某个Φ(x)的内积呢。
该想法是可行的,函数的具体形式为,并称为核函数 (kernel function)。
Mercer定理:只要函数对定义域内的x和y形成的二维矩阵总是半正定的,那么它就能被当作一个核函数来使用,也就意味着它能表示一个低维映射到高维的高维向量之间的内积。
8.软间隔
无论是线性可分的SVM还是核方法的SVM都有一个共同点,都是严格的将两种类别映射到分界面的两边,并且留出巨大的间隔。这种使两类样本点必须都要分类正确的SVM称为硬间隔,会带来一些弊端。介绍另一种策略,即软间隔。
软间隔,仍然采取间隔最大化的思想来寻找SVM的超平面,但是不完全严格的要求所有点都被正确无误的分到被支持向量确定出来的边界以内。允许错分,但是错分的点要受到惩罚。该惩罚写入目标函数,与间隔最大化一起优化。
引入松弛变量,使得约束条件修正为以下形式。
实验
SMO
为更加容易的求解SVM的最优解,SMO可以将大优化问题分解为多个小优化问题来求解。
SMO的工作原理是:每次循环中选择两个alpha进行优化处理。一旦找到一对合适的alpha,那么就增大其中一个同时减小另外一个。合适的alpha满足两个条件,1、两个alpha必须要在分隔边界之外,2、两个alpha没有进行过区间优化处理或不在边界上。
只改比那一个alpha又可能导致约束条件alpha与对应标签的积的累加=0失效,因此同时改变两个alpha。
完整版SMO函数的伪代码:
通过一个外循环来选择第一个alpha值,并且其选择过程会在两种方式之间进行交替:
1、在所有数据集上进行单遍扫描,2、在非边界alpha中实现单遍扫描(不等于0或C的alpha值)
再选择第一个alpha值后,会通过一个内循环选择第二个alpha值。会通过最大化步长的方式获得第二个alpha。在这里我们会建立一个全局的缓存用于保存误差值,并从中选择使得步长或Ei - Ej最大的alpha值。
添加用于清理代码的数据结构和3个用于对E进行缓存的辅助函数。
1、函数说明:维护所有需要操作的值
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]
# 根据矩阵行数初始化alphas矩阵,一个m行1列的全零列向量
self.alphas = np.mat(np.zeros((self.m, 1)))
# 初始化b参数为0
self.b = 0
# 根据矩阵行数初始化误差缓存矩阵,第一列为是否有效标志位,第二列为实际的误差E的值
self.eCache = np.mat(np.zeros((self.m, 2)))
# 初始化核K
self.K = np.mat(np.zeros((self.m, self.m)))
# 计算所有数据的核K
for i in range(self.m):
self.K[:, i] = kernelTrans(self.X, self.X[i, :], kTup)
2、函数说明:计算误差
def calcEk(oS, k):
# multiply(a,b)就是个乘法,如果a,b是两个数组,那么对应元素相乘
# .T为转置
fXk = float(np.multiply(oS.alphas, oS.labelMat).T * oS.K[:, k] + oS.b)
# 计算误差项
Ek = fXk - float(oS.labelMat[k])
# 返回误差项
return Ek
3、函数说明:内循环启发方式
def selectJ(i, oS, Ei):
# 初始化
maxK = -1
maxDeltaE = 0
Ej = 0
# 根据Ei更新误差缓存
oS.eCache[i] = [1, Ei]
# 对一个矩阵.A转换为Array类型
# 返回误差不为0的数据的索引值
validEcacheList = np.nonzero(oS.eCache[:, 0].A)[0]
# 有不为0的误差
if (len(validEcacheList) > 1):
# 遍历,找到最大的Ek
for k in validEcacheList:
# 不计算k==i节省时间
if k == i:
continue
# 计算Ek
Ek = calcEk(oS, k)
# 计算|Ei - Ek|
deltaE = abs(Ei - Ek)
# 找到maxDeltaE
if (deltaE > maxDeltaE):
maxK = k
maxDeltaE = deltaE
Ej = Ek
# 返回maxK,Ej
return maxK, Ej
# 没有不为0的误差
else:
# 随机选择alpha_j的索引值
j = selectJrand(i, oS.m)
# 计算Ej
Ej = calcEk(oS, j)
# 返回j,Ej
return j, Ej
4 函数说明:计算Ek,并更新误差缓存
def updateEk(oS, k):
# 计算Ek
Ek = calcEk(oS, k)
# 更新误差缓存
oS.eCache[k] = [1, Ek]
5
def kernelTrans(X, A, kTup): #calc the kernel or transform data to a higher dimensional space
m,n = np.shape(X)
K = np.mat(np.zeros((m,1)))
if kTup[0]=='lin': K = X * A.T #linear kernel
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
6 函数说明:优化的SMO算法
def innerL(i, oS):
# 步骤1:计算误差Ei
Ei = calcEk(oS, i)
# 优化alpha,设定一定的容错率
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)):
# 使用内循环启发方式2选择alpha_j,并计算Ej
j, Ej = selectJ(i, oS, Ei)
# 保存更新前的alpha值,使用深层拷贝
alphaIold = oS.alphas[i].copy()
alphaJold = oS.alphas[j].copy()
# 步骤2:计算上界H和下界L
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
# 步骤3:计算eta
eta = 2.0 * oS.K[i, j] - oS.K[i, i] - oS.K[j, j]
if eta >= 0:
print("eta >= 0")
return 0
# 步骤4:更新alpha_j
oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej) / eta
# 步骤5:修剪alpha_j
oS.alphas[j] = clipAlpha(oS.alphas[j], H, L)
# 更新Ej至误差缓存
updateEk(oS, j)
if (abs(oS.alphas[j] - alphaJold) < 0.00001):
print("alpha_j变化太小")
return 0
# 步骤6:更新alpha_i
oS.alphas[i] += oS.labelMat[i] * oS.labelMat[j] * (alphaJold - oS.alphas[j])
# 更新Ei至误差缓存
updateEk(oS, i)
# 步骤7:更新b_1和b_2:
b1 = oS.b - Ei - oS.labelMat[i] * (oS.alphas[i] - alphaIold) * oS.K[i, i] - oS.labelMat[j] * (
oS.alphas[j] - alphaJold) * oS.K[j, i]
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]
# 步骤8:根据b_1和b_2更新b
if (0 < oS.alphas[i] < oS.C):
oS.b = b1
elif (0 < oS.alphas[j] < oS.C):
oS.b = b2
else:
oS.b = (b1 + b2) / 2.0
return 1
else:
return 0
7 函数说明:完整的线性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
entrieSet = True
alphaPairsChanged = 0
# 遍历整个数据集alpha都没有更新或者超过最大迭代次数,则退出循环
while (iter < maxIter) and ((alphaPairsChanged > 0) or (entrieSet)):
alphaPairsChanged = 0
# 遍历整个数据集
if entrieSet:
for i in range(oS.m):
# 使用优化的SMO算法
alphaPairsChanged += innerL(i, oS)
print("全样本遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter, i, alphaPairsChanged))
iter += 1
# 遍历非边界值
else:
# 遍历不在边界0和C的alpha
nonBoundIs = np.nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
for i in nonBoundIs:
alphaPairsChanged += innerL(i, oS)
print("非边界遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter, i, alphaPairsChanged))
iter += 1
# 遍历一次后改为非边界遍历
if entrieSet:
entrieSet = False
# 如果alpha没有更新,计算全样本遍历
elif (alphaPairsChanged == 0):
entrieSet = True
print("迭代次数:%d" % iter)
# 返回SMO算法计算的b和alphas
return oS.b, oS.alphas
8 函数说明:计算w
def calcWs(alphas,dataArr,classLabels):
X = np.mat(dataArr); labelMat = np.mat(classLabels).transpose()
m,n = np.shape(X)
w = np.zeros((n,1))
for i in range(m):
w += np.multiply(alphas[i]*labelMat[i],X[i,:].T)
return w
9 函数说明:测试
def test(dataArr,labelArr):
datMat = np.mat(dataArr)
m,n = np.shape(datMat)
error = 0
w = datMat * np.mat(ws)+b
for i in range(m):
if w[i]*labelArr[i] < 0:
error=error+1;
return 1-error/m
10 函数说明:读取数据
def loadDataSet(fileName):
# 数据矩阵
dataMat = []
# 标签向量
labelMat = []
# 打开文件
fr = open(fileName)
# 逐行读取
for line in fr.readlines():
# 去掉每一行首尾的空白符,例如'\n','\r','\t',' '
# 将每一行内容根据'\t'符进行切片
lineArr = line.strip().split('\t')
# 添加数据(100个元素排成一行)
dataMat.append([float(lineArr[0]), float(lineArr[1])])
# 添加标签(100个元素排成一行)
labelMat.append(float(lineArr[2]))
return dataMat, labelMat
11 运行模型测试
dataArr, labelArr = loadDataSet('testSet.txt')
b, alphas = smoP(dataArr, labelArr, 0.6, 0.001, 40)
ws = calcWs(alphas,dataArr,labelArr)
test(dataArr,labelArr)
//正确率 1
模型还只是针对线性可分的情况,对于非线性可分的情况,还可以如前面分析的采取软间隔和核函数的方法进行划分。