有监督学习-----支持向量机(含python实现代码)
支持向量机（support vector machines, SVM）是一种二分类模型，它的基本模型是定义在特征空间上的间隔最大的线性分类器；SVM的的学习策略就是间隔最大化，SVM的的学习算法就是求解凸二次规划的最优化算法。SVM的实现中，最流行的一种实现是序列最小优化（Sequential Minimal Optimization, SMO）。SVM可以使用核函数（kernel）的方式，可以对
支持向量机（support vector machines, SVM）是一种二分类模型，它的基本模型是定义在特征空间上的间隔最大的线性分类器；SVM的的学习策略就是间隔最大化，SVM的的学习算法就是求解凸二次规划的最优化算法。
SVM的实现中，最流行的一种实现是序列最小优化（Sequential Minimal Optimization, SMO）。SVM可以使用核函数（kernel）的方式，可以对非线性可分的数据进行分类。
原理
希望找到离分隔超平面最近的点，确保它们离分隔面的距离尽可能远。这里点到分隔面的距离被称为间隔。我
们希望间隔尽可能地大，这是因为如果我们犯错或者在有限数据上训练分类器的话，我们希望分类器尽可能健壮。
支持向量（support vector)就是离分隔超平面最近的那些点。接下来要试着最大化支持向量到分隔面的距离，需要找到此问题的优化求解方法。
寻找最大间隔
分隔超平面的形式可以写成 w T x + b w^Tx+b wTx+b。要计算点A到分隔超平面的距离，就必须给出点到分隔面的法线或垂线的长度，该值为 ∣ w T A + b ∣ / ∣ ∣ w ∣ ∣ |w^TA+b|/ ||w|| ∣wTA+b∣/∣∣w∣∣。这里的常数b类似于逻辑回归中的结局 w 0 w_0 w0。这里的向量w和常数b一起描述了所给数据的分隔线或超平面。
分类器求解的优化问题
输入数据给分类器会输出一个类别标签。下面将使用单位阶跃函数对 w T x + b w^Tx+b wTx+b作用得到 f ( w T x + b ) f(w^Tx+b) f(wTx+b)，其中，当u<0时f(u)输出-1，反之则输出+1。这里的标签采用-1和+1，仅仅相差一个符号，方便数学上处理。可以通过一个统一的公式来表示间隔或者数据点到分隔超平面的距离，同时不必担心数据到底是属于-1还是+1。
当计算数据点到分隔面的距离并确定分隔面的放置位置时，间隔通过 l a b e l ∗ ( w T x + b ) label*(w^Tx+b) label∗(wTx+b)来计算，这时就能体现出-1和+1的好处了。如果数据点处于正方向（即+1类）并且离分隔超平面很远的位置时， w T x + b w^Tx+b wTx+b会是一个很大的正数，同时 l a b e l ∗ ( w T x + b ) label*(w^Tx+b) label∗(wTx+b)也会是一个很大的正数。如果数据点处于负方向（即-1类）并且离分隔超平面很远的位置时，此时由于类别标签为-1， l a b e l ∗ ( w T x + b ) label*(w^Tx+b) label∗(wTx+b)仍然是一个很大的正数。
现在的目标就是找出分类器定义中的w和b。为此，必须找到具有最小间隔的数据点，这些数据点也就是前面提到的支持向量。一旦找到具有最小间隔的数据点，我们就需要对该间隔最大化。可以写作：
a r g max w , b { min n ( l a b e l ∗ ( w T x + b ) ⋅ 1 ∥ w ∥ } arg \ \max_{w,b} \ \left\{ \ \min_{n}(label*(w^Tx+b)·\frac{1}{\Vert w \Vert} \right\} arg w,bmax { nmin(label∗(wTx+b)⋅∥w∥1}
直接求解相当困难，可以将它转换为另一种更容易求解的形式 。首先考虑上式中大括号内的部分。我们要做的是固定其中一个因子而最大化其他因子。如果令所有支持向量的 l a b e l ∗ ( w T x + b ) label*(w^Tx+b) label∗(wTx+b)都等于1，那么久可以通过求 ∥ w ∥ \Vert w \Vert ∥w∥的最大化来得到最终解。但是，并非所有数据点的 l a b e l ∗ ( w T x + b ) label*(w^Tx+b) label∗(wTx+b)都等于1，只有那些离分隔超平面最近的点得到值才为1.而离超平面越远的数据点，其 l a b e l ∗ ( w T x + b ) label*(w^Tx+b) label∗(wTx+b)值也就越大。
在上述优化问题中，给定了一些约束条件然后求最优值，因此该问题是一个带约束条件的优化问题。这里的约束条件就是 l a b e l ∗ ( w T x + b ) > = 1.0 label*(w^Tx+b)>=1.0 label∗(wTx+b)>=1.0.对于这类优化问题，有一个非常著名的求解方法，即拉格朗日乘子法。通过引入拉个朗日乘子，就可以基于约束条件来表述原来的问题。由于这里的约束条件都是基于数据点的，因此可以将超平面携程数据点的形式。于是优化目标函数最后可以写成：
max a [ ∑ i = 1 m α − 1 2 ∑ i , j = 1 m l a b e l ( i ) ⋅ l a b e l ( j ) ⋅ α i ⋅ α j ⟨ x ( i ) , x ( j ) ⟩ ] \max_{a}{\left [ \sum_{i=1}^{m}{\alpha} -\frac{1}{2}\sum_{i,j=1}^{m}{label^{(i)}·label^{(j)}·\alpha_i·\alpha_j\langle x^{(i)},x^{(j)}\rangle} \right ]} amax[i=1∑mα−21i,j=1∑mlabel(i)⋅label(j)⋅αi⋅αj⟨x(i),x(j)⟩]
其约束条件为
α ≥ 0 , 和 ∑ i − 1 m α i ⋅ l a b e l ( i ) = 0 \alpha \geq0, \ 和\sum_{i-1}^{m}{\alpha_i·label^{(i)}}=0 α≥0, 和i−1∑mαi⋅label(i)=0
但是这里有个假设：数据必须100%线性可分。但是，几乎所有数据都不那么“干净”。这时就可以通过引入松弛变量来允许有些数据点可以处于分隔面错误一侧。这里的优化目标就能保持仍然不变，但是此时新的约束条件则变为：
C ≥ α ≥ 0 , 和 ∑ i − 1 m α i ⋅ l a b e l ( i ) = 0 C \geq \alpha \geq0, \ 和\sum_{i-1}^{m}{\alpha_i·label^{(i)}}=0 C≥α≥0, 和i−1∑mαi⋅label(i)=0
这里的常数C用于控制“最大化间隔”和“保证大部分点的函数间隔小于1.0”这两个目标的权重。在优化算的实现代码中，常数C是一个参数，因此就可以通过调节该参数得到不同的结果。一旦求出了所有的alpha，那么分隔超平面就可以通过这些alpha来表达。SVM的主要工作就是求解这些alpha。
SMO高效优化算法
序列最小优化（Sequential Minimal Optimization, SMO），Platt的SMO算法是将大优化问题分解为多个小优化问题求解的。这些小优化问题往往很容易求解，并且对它们进行顺序求解的结果与将它们作为整体来求解的结果是完全一致的。在结果完全相同的同时，SMO算法的求解时间短很多。
SMO算法的目标是求出一系列alpha和b，一旦求出这些alpha，就很容易计算出权重向量，并得到分隔超平面。
SMO算法的工作原理是：每次循环中选择两个alpha进行优化处理。一旦找到一堆合适的alpha，那么就增大其中一个同时减小另一个。这里所谓的“合适”就是指两个alpha必须要符合一定的条件，条件之一就是这两个alpha必须要在间隔边界之外，而其中第二个条件则是这两个alpha还没有进行过区间化处理或者不在边界上。
应用简化版SMO算法处理小规模数据集
简化版代码虽然量少，但执行速度慢。简化版首先在数据集上遍历每个alpha,然后在剩下的alpha集合中随机选择另一个alpha，从而构建alpha对。即就是要同时改变两个alpha。
为此，将构建一个辅助函数，用于在某个区间范围内随机选择一个整数。同时需要另一个辅助函数，用于在数值太大时对其进行调整。
python实现
import numpy as np
def loadDataSet(filename):
"""
加载数据
:param filename: 文件名
:return: 数据集 和 类别标签
"""
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):
"""
用于在[0, m)区间范围内随机选择一个整数，这个整数不等于i
:param i: 第一个α的下标
:param m: 所有α的数目
:return: 第二个α的下标
"""
j = i
while (j == i):
j = int(np.random.uniform(0, m)) # 在[0，m)区间随机生成一个实数，并将该数转换为整数
return j
def clipAlpha(aj, H, L):
"""
用于调整大于H或小于L的α值
:param aj: 当前α值
:param H:最大值
:param L:最小值
:return: 调整后的α值
"""
if aj > H:
aj = H
if aj < L:
aj = L
return aj
def smoSimple(dataIn, classLabels, C, toler, maxIter):
"""
简化版的SMO（最小序列优化）
:param dataIn: 数据集
:param classLabels: 类别标签
:param C: 松弛变量，用于控制“最大化间隔”和“保证大部分点的函数间隔小于1.0” 这两个目标的权重
:param toler: 容错率
:param maxIter: 最大迭代次数
:return: 常数b和矩阵alphas
"""
dataMat = np.mat(dataIn) # 将数据集转换为numpy矩阵
labelMat = np.mat(classLabels).transpose() # 转换为列向量
b = 0
m, n = np.shape(dataMat)
alphas = np.mat(np.zeros((m, 1))) # 构建alphas列矩阵，初始化元素都为0
iter = 0 # 在没有任何α改变的情况下遍历数据集的次数
while (iter < maxIter):
alphaPairsChanged = 0 # 记录alpha是否巳经进行优化
for i in range(m):
# fXi为预测的类别
fXi = float(np.multiply(alphas, labelMat).T * \
(dataMat * dataMat[i, :].T)) + b
Ei = fXi - float(labelMat[i]) # 预测值和真实值的误差
# 如果误差很大，可以对该数据实例所对应的α值进行优化
if ((labelMat[i] * Ei < -toler) and (alphas[i] < C)) or \
((labelMat[i] * Ei > toler) and (alphas[i] > 0)):
j = selectJrand(i, m) # 随机选择第二个α
fXj = float(np.multiply(alphas, labelMat).T * \
(dataMat * dataMat[j, :].T)) + b
Ej = fXj - float(labelMat[j])
# 方便后面将新的α和老的α进行比较
alphaIold = alphas[i].copy()
alphaJold = alphas[j].copy()
# 保证α在0与C之间
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])
# 如果L和H相等，则不作任何改变
if L == H:
print("L==H")
continue
# eta 是alphas[j]的最优修改量
eta = 2.0 * dataMat[i, :] * dataMat[j, :].T - \
dataMat[i, :] * dataMat[i, :].T - \
dataMat[j, :] * dataMat[j, :].T
# 如果eta>= 0，则需要退出for循环的当前迭代过程
if eta >= 0:
print("eta>=0")
continue
# 计算出一个新的alphas[j]
alphas[j] -= labelMat[j] * (Ei - Ej) / eta
alphas[j] = clipAlpha(alphas[j], H, L)
# 检查alphas[j]是否有轻微改变
if (abs(alphas[j] - alphaJold) < 0.00001):
print("j not moving enough")
continue
# 对i进行修改，修改量与j相同，但方向相反
alphas[i] += labelMat[j] * labelMat[i] * (alphaJold - alphas[j])
# 给这alphas[i]和alphas[j]值设置一个常数项b
b1 = b - Ei - labelMat[i] * (alphas[i] - alphaIold) * \
dataMat[i, :] * dataMat[i, :].T - \
labelMat[j] * (alphas[j] - alphaJold) * \
dataMat[i, :] * dataMat[j, :].T
b2 = b - Ej - labelMat[i] * (alphas[i] - alphaIold) * \
dataMat[i, :] * dataMat[j, :].T - \
labelMat[j] * (alphas[j] - alphaJold) * \
dataMat[j, :] * dataMat[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 # 表示已经成功的改变了一对alpha
print("iter: %d i: %d, pairs changed %d" % ( \
iter, i, alphaPairsChanged))
# 检查alpha是否做了更新
if (alphaPairsChanged == 0):
iter += 1
else: # 如果有更新则将iter为0后继续运行程序
iter = 0
print("iteration number: %d" % iter)
return b, alphas
def calcWs(alphas,dataArr,classLabels):
"""
计算权重w
:param alphas: alphas矩阵
:param dataArr: 数据集
:param classLabels: 类别标签
:return: 权重w
"""
X = np.mat(dataArr)
labelMat = np.mat(classLabels).transpose()
m,n = np.shape(X)
w = np.zeros((n,1))
for i in range(m):
# 多个数的乘积，由于大部分alphas0，非零alpha所对应的也就是支持向量，所以最终起作用只有支持向量
w += np.multiply(alphas[i]*labelMat[i],X[i,:].T)
return w
运行下面代码
dataIn, classLabels = loadDataSet('testSet.txt')
b, alphas = smoSimple(dataIn, classLabels, 0.6, 0.001, 40)
print("b:\n", b)
# 由于alpha矩阵本身零元素太多，所以只观察大于0的元素
print("alphas大于0：\n", alphas[alphas > 0])
# 为了了解哪些是支持向量
print("*"*50)
print("支持向量如下：")
for i in range(100):
if alphas[i]>0.0:
print(dataIn[i],classLabels[i])
# 计算权重
ws = calcWs(alphas, dataIn, classLabels)
print("ws:\n",ws)
运行结果如下：
在原始数据集上对这些支持向量画圈之后的结果如下
对于只有一百个点组成的小规模数据集，简化的SMO算法是没有问题的，但是，在更大的数据集上运行速度就会很慢。
利用完整Platt SMO 算法加速优化
Platt SMO 算法是通过一个外循环来选择第一个alpha值的，并且其选择过程会在梁中锋哪个是之间进行交替：一种是在所有数据集上进行单遍扫描，另一种是在非边界alpha中实现单遍扫描。而所谓的非边界alpha指的就是那些不等于边界0或C的alpha值。对整个数据集的扫描相当容易，而实现非边界alpha值的扫描时，首先需要建立这些alpha值的列表，然后再对这个表进行遍历。同时，该步骤会跳过那些已知的不会改变的alpha值。
在选择第一个alpha值后，算法会通过一个内循环来选择第二个alpha值。在优化过程中，会通过最大步长的方式来获得第二个alpha值。在简化版SMO算法中，我们会在选择j之后计算错误率Ej。但是在这里，我们会建立一个全局的缓存用于保存误差值，并从中选择使得步长或者Ei-Ej最大的alpha值。
python实现
class optStruct:
"""构建一个类，以实现其成员变量的填充"""
def __init__(self, dataIn, classLabels, C, toler):
self.X = dataIn
self.labelMat = classLabels
self.C = C
self.tol = toler
self.m = np.shape(dataIn)[0]
self.alphas = np.mat(np.zeros((self.m, 1)))
self.b = 0
self.eCache = np.mat(np.zeros((self.m, 2))) # 误差缓存 [eCache是否有效，实际的E值]
def calcEk(oS, k):
"""
计算误差E值
:param oS:类对象
:param k:k
:return:误差Ek
"""
fXk = float(np.multiply(oS.alphas, oS.labelMat).T * \
(oS.X * oS.X[k, :].T)) + oS.b
Ek = fXk - float(oS.labelMat[k])
return Ek
def selectJ(i, oS, Ei):
"""
用于选择第二个alpha值（即内循环的alpha值）
目标是选择合适的第二个alpha值以保证在每次优化中采用最大步长
:param i:第一个alpha的下标
:param oS:对象
:param Ei:误差
:return:第二个alpha的下标，误差Ej
"""
maxK = -1
maxDeltaE = 0
Ej = 0
oS.eCache[i] = [1, Ei] # 将输人值Ei在缓存中设置成为有效的,这里的有效是它已经计算好了
validEcacheList = np.nonzero(oS.eCache[:, 0].A)[0] # 非零E值所对应的alpha值
if (len(validEcacheList)) > 1:
for k in validEcacheList:
if k == i:
continue
Ek = calcEk(oS, k)
deltaE = abs(Ei - Ek)
if (deltaE > maxDeltaE): # 选择具有最大步长的j
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):
"""
计算误差值并存人缓存当
:param oS:对象
:param k:k
"""
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]
def innerL(i, oS):
"""
完整Platt SMO算法中的用于寻找决策边界的优化例程
:param i: 第一个alpha的下标
:param oS: 类对象
:return: 如果任意一对alpha发生改变，返回1，否则，返回0
"""
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) # 在alpha值改变时更新误差缓存
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
def smoP(dataIn,classLabels,C,toler,maxIter,kTup=('lin',0)):
"""
完整版Platt SMO 的外循环代码
:param dataIn: 数据集
:param classLabels: 类别标签
:param C: 松弛变量
:param toler: 容错率
:param maxIter: 最大迭代次数
:param kTup: 包含核函数信息的元组，主要在核函数处用到
:return:返回常数b 和 alphas值
"""
oS = optStruct(np.mat(dataIn),np.mat(classLabels).transpose(),C,toler)
iter = 0
entireSet = True
alphaPairsChanged = 0
# 当迭代次数超过指定的最大值 或者遍历整个集合都未对任意alpha对进行修改时，就退出循环
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((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
elif (alphaPairsChanged == 0):
entireSet = True
print("iteration number: %d" % iter)
return oS.b, oS.alphas
def plotfig_SVM(dataIn, classLabels, ws, b, alphas):
"""
绘制支持向量
:param dataIn: 数据集
:param classLabels: 类别标签
:param ws: 权重
:param b: 常数
:param alphas: 矩阵alpha
"""
import matplotlib.pyplot as plt
from matplotlib.patches import Circle
xMat = np.mat(dataIn)
yMat = np.mat(classLabels)
# b原来是矩阵，先转为数组类型后其数组大小为（1,1），所以后面加[0]，变为(1,)
b = np.array(b)[0]
fig = plt.figure()
ax = fig.add_subplot(111)
# 注意flatten的用法
ax.scatter(xMat[:, 0].flatten().A[0], xMat[:, 1].flatten().A[0])
x = np.arange(-2.0, 12.0, 0.1)
# 根据wT*x + b = 0 得到w0.x1 + w1.x2 + b = 0,转换形式得到x2(x2就是y值)
y = (-b - ws[0, 0] * x) / ws[1, 0]
ax.plot(x, y) # 画出分隔面
# 画出两类数据的散点图，并用颜色、形状区分
for i in range(np.shape(yMat[0, :])[1]):
if yMat[0, i] > 0:
ax.scatter(xMat[i, 0], xMat[i, 1], marker='s', s=50, c='blue')
else:
ax.scatter(xMat[i, 0], xMat[i, 1], marker='o', s=50, c='red')
# 找到支持向量，并在图中用圈标记
for i in range(100):
if alphas[i] > 0.0:
circle = Circle((xMat[i, 0], xMat[i, 1]), 0.5, facecolor='none',
edgecolor=(0, 0.8, 0.8), linewidth=3, alpha=0.5)
ax.add_patch(circle)
# 限定X Y轴的区间范围
ax.axis([-2, 12, -8, 6])
plt.show()
运行下面代码
b,alphas = smoP(dataIn, classLabels,0.6,0.001,40)
print("b:", b)
print("alphas大于0", alphas[alphas > 0])
ws = calcWs(alphas,dataIn, classLabels)
print("ws:\n",ws)
dataMat = np.mat(dataIn)
# 对数据进行分类
# 该值大于0，那么其属于1类；如果该值小于0，那么则属于-1类
print("预测值：",dataIn[0]*np.mat(ws) + b)
print("真实值：",dataIn[0])
print("预测值：", dataIn[1] * np.mat(ws) + b)
print("真实值：", dataIn[1])
print("预测值：", dataIn[2] * np.mat(ws) + b)
print("真实值：", dataIn[2])
# 绘制支持向量
plotfig_SVM(dataIn, classLabels, ws, b, alphas)
运行结果如下
在原始数据集上对这些支持向量画圈之后的结果如下
相对于简化版SMO算法，完整版Platt SMO算法运行速度更快。这里的数据点都是分布在一条直线的两边。但是，倘若两类数据蒂娜分别分布在一个圆的内部和外部，那么会得到怎样的分类面呢？下一节介绍一种方法对分类器进行修改，以说明类别区域形状不同情况下的数据集分隔问题。
在复杂数据上应用核函数
使用核函数的工具将数据转换成易于分类器理解的形式。首先解释核函数的概念，并介绍他们在支持向量机中的使用方法。然后，介绍一种称为径向基函数的最流行的核函数。最后，将该核函数应用与前面得到的分类器。
利用核函数将数据映射到高维空间
对于下图，如果只在x和y轴构成的坐标系中插入直线进行分类的话，并不会得到想要的结果。或许可以对圆中的数据进行某种形式的转换，从而得到某些新的变量来表示数据。在这种表示情况下，就更容易得到大于0或者小于0的测试结果。我们可以将数据从一个特征空间转换到另一个特征空间。在新的空间下，可以很容易利用已有的工具对数据进行处理。数学家们喜欢将这个过程称为从一个特征空间到另一个特征空间的映射。在通常情况下，这种映射会将低维特征空间映射到高维空间。
这种从某个特征空间到另一个特征空间的映射是通过核函数来实现的。核函数具有多种类型，经过空间转换后，可以在高维空间解决线性问题，这也就等价于在低维空间解决非线性问题。
SVM优化中一个特别好的地方及时，所有的运算都可以写成内积的形式。向量的内积指的是两个向量相乘，之后得到单个标量或者数值。可以把内积运算替换成核函数，而不必做建华处理。将内积替换成核函数的方式称为核技巧。
很多机器学习算法都用到核函数。其中，最流行的核函数是径向基核函数
径向基核函数
径向基函数是SVM中常用的一个核函数。径向基函数是一个采用向量作为自变量的函数，能够基于向量距离运算赎回粗一个标量。径向基函数的搞死版本，具体公式如下
k ( x , y ) = e x p ( − ∥ x − y ∥ 2 2 σ 2 ) k(x,y)=exp\left( \frac{-\Vert x-y \Vert^2}{2\sigma^2} \right) k(x,y)=exp(2σ2−∥x−y∥2)
其中， σ \sigma σ是用户定义的用于确定到达率或函数值跌落到0的速度参数。
上述高斯核函数将数据从其特征空间映射到更高维的空间，它会得到一个理想的结果，也能得到低错误率。
python实现
def kernelTrans(X, A, kTup):
"""
利用核函数将数据从低维特征空间因社会到高维特征空间
:param X:数据集
:param A:某条数据
:param kTup:包含核函数信息的元组，第一个参数是描述所有核函数的类型，第二个参数是可能需要的可选参数
:return:转换后的数据
"""
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循环中对于矩阵的每个元素计算高斯函数的值
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, dataIn, classLabels, C, toler, kTup):
"""
初始化
:param dataIn:数据集
:param classLabels: 类别标签
:param C: 松弛变量
:param toler: 容错率
:param kTup: 包含核函数信息的元组
"""
self.X = dataIn
self.labelMat = classLabels
self.C = C
self.tol = toler
self.m = np.shape(dataIn)[0]
self.alphas = np.mat(np.zeros((self.m, 1)))
self.b = 0
self.eCache = np.mat(np.zeros((self.m, 2)))
# 构建矩阵K，并通过kernelTrans()方法进行填充，全局的K值只需计算一次
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()函数进行修改，直接用下面的即可
def innerL(i, oS):
"""
完整Platt SMO算法中的用于寻找决策边界的优化例程
:param i: 第一个alpha的下标
:param oS: 类对象
:return: 如果任意一对alpha发生改变，返回1，否则，返回0
"""
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.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.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
def calcEk(oS, k):
"""
计算误差E值
:param oS:类对象
:param k: k
:return:误差Ek
"""
fXk = float(np.multiply(oS.alphas, oS.labelMat).T * oS.K[:, k] + oS.b)
Ek = fXk - float(oS.labelMat[k])
return Ek
def plotSupportVectors_RBF(dataIn, classLabels, alphas):
"""
绘制径向基核函数的支持向量
:param dataIn: 数据集
:param classLabels: 类别标签
:param alphas: 矩阵alpha
"""
import matplotlib.pyplot as plt
from matplotlib.patches import Circle
xMat = np.mat(dataIn)
yMat = np.mat(classLabels)
fig = plt.figure()
ax = fig.add_subplot(111)
# 注意flatten的用法
ax.scatter(xMat[:, 0].flatten().A[0], xMat[:, 1].flatten().A[0])
# 画出两类数据的散点图，并用颜色、形状区分
for i in range(np.shape(yMat[0, :])[1]):
if yMat[0, i] > 0:
ax.scatter(xMat[i, 0], xMat[i, 1], marker='o', s=50, c='blue')
else:
ax.scatter(xMat[i, 0], xMat[i, 1], marker='s', s=50, c='red')
# 找到支持向量，并在图中用圈标记
for i in range(100):
if alphas[i] > 0.0:
circle = Circle((xMat[i, 0], xMat[i, 1]), 0.06, facecolor='none',
edgecolor=(0, 0.8, 0.8), linewidth=3, alpha=0.5)
ax.add_patch(circle)
# 限定X Y轴的区间范围
ax.axis([-1.5, 1.5, -1.0, 1.5])
plt.show()
# 利用核函数进行分类的径向基测试函数
def testTbf(k1=1.3):
"""
利用核函数进行分类的径向基测试函数
:param k1:输入参数
"""
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()
# 获取大于0的alpha值对应的索引，以及对应的数据集和标签
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):
# 利用核函数进行分类
kenelEval = kernelTrans(sVs,dataMat[i,:],('rbf',k1))
predict = kenelEval.T * np.multiply(labelSV,alphas[svInd]) + b
# 如果预测的值和测试的值不一样，则错误数加1
if np.sign(predict) != np.sign(labelMat[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):
kenelEval = kernelTrans(sVs,dataMat[i,:],('rbf',k1))
predict = kenelEval.T * np.multiply(labelSV,alphas[svInd]) + b
if np.sign(predict) != np.sign(labelMat[i]):
errorCount += 1
print("the test error rate is : %f" %(float(errorCount)/m))
plotSupportVectors_RBF(dataArr, labelArr, alphas)
运行下面语句
testTbf(k1=1.3)
得到结果
可以尝试不同的k1参数来观察测试集的错误率、训练集出错误率、支持向量个数随k1的变化情况。
上图是k1=0.1时的结果，该参数减少了每个支持向量的影响程度。因此需要更多的支持向量。可以发现上图中共有100个数据点，其中有88个支持向量。
将k1=1.3会得到下图的结果，这里的支持向量少于上图，而且，此时的测试错误率也在下降。该数据集在这个设置的某处存在着最优值。如果降低k1,那么训练错误率会降低，但是测试错误率却会上升。
支持向量的数目存在一个最优值。SVM的优点在于它能对数据进行高效分类。如果支持向量太少，就可能会得到一个很差的决策边界；如果支持向量太多，就相当于每次都利用整个数据集进行分类。
示例：手写识别
主要包括以下方法：
- img2vector()：把一个32*32的二进制图像矩阵转换为1*1024的向量
- loadImage()：对图片进行加载并处理，获取数据和类别标签
- smoP()：完整版Platt SMO 的外循环代码，见利用完整Platt SMO 算法加速优化的python实现
- testDigits()：用核函数进行手写数字的识别
def img2vector(filename):
"""
把一个32*32的二进制图像矩阵转换为1*1024的向量
:param filename: 文件名，带后缀的
:return: 转换为向量的数据
"""
returnVect = np.zeros((1, 1024))
fr = open(filename)
for i in range(32):
lineStr = fr.readline()
for j in range(32):
returnVect[0, i * 32 + j] = int(lineStr[j])
return returnVect
def loadImage(dirName):
"""
对图片进行加载并处理，获取数据和类别标签
:param dirName: 文件夹名称
:return:数据和类别标签
"""
from os import listdir
hwLabels = []
trainingFileList = listdir(dirName)
m = len(trainingFileList)
traningMat = np.zeros((m, 1024))
for i in range(m):
fileNameStr = trainingFileList[i] # 文件名，带后缀的
fileStr = fileNameStr.split('.')[0] # 文件名，不带后缀
classNumStr = int(fileStr.split('_')[0])
if classNumStr == 9: # 一旦碰到数字9，则标签为-1
hwLabels.append(-1)
else: # 否则标签为1
hwLabels.append(1)
traningMat[i, :] = img2vector('%s/%s' % (dirName, fileNameStr))
return traningMat, hwLabels
def testDigits(kTup=('rbf', 10)):
"""
用核函数进行手写数字的识别
:param kTup:包含核函数信息的元组
"""
# 获得数据和类别标签
dataArr, labelArr = loadImage('trainingDigits')
b, alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, kTup)
dataMat = np.mat(dataArr)
labelMat = np.mat(labelArr).transpose()
svInd = np.nonzero(alphas > 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, :], kTup)
predict = kernelEval.T * np.multiply(labelSV, alphas[svInd]) + b
if np.sign(predict) != np.sign(labelMat[i]):
errorCount += 1
print("the training error rate is : %f" % (float(errorCount) / m))
dataArr, labelArr = loadImage('testDigits')
dataMat = np.mat(dataArr)
labelMat = np.mat(labelArr).transpose()
m, n = np.shape(dataMat)
errorCount = 0
for i in range(m):
kernelEval = kernelTrans(sVs, dataMat[i, :], kTup)
predict = kernelEval.T * np.multiply(labelSV, alphas[svInd]) + b
if np.sign(predict) != np.sign(labelMat[i]):
errorCount += 1
print("the test error rate is : %f" % (float(errorCount) / m))
运行代码如下：
testDigits(('rbf', 20))
运行结果如下：
可以尝试不同的k值和线性核函数，来对比结果。
优缺点
优点：
- 有严格的数学理论支持，可解释性强，不依靠统计方法，从而简化了通常的分类和回归问题；
- 能找出对任务至关重要的关键样本，即：支持向量；
- 采用核方法后，可以处理非线性分类/回归任务；
- 最终决策函数只由少数的支持向量所确定，计算的复杂性取决于支持向量的数目，而不是样本空间的维数，这在某种意义上避免了“维数灾难”。
缺点：
-
对参数调节和核函数的选择敏感，原始分类器不加修改仅适用于处理二分类问题
-
训练时间长。当采用 SMO 算法时，由于每次都需要挑选一对参数，因此时间复杂度为 O ( N 2 ) O(N^2) O(N2) ，其中 N 为训练样本的数量；
-
采用核方法时，如果需要存储核矩阵，则空间复杂度为 O ( N 2 ) O(N^2) O(N2) ；
-
模型预测时，预测时间与支持向量的个数成正比。当支持向量的数量较大时，预测计算复杂度较高。
因此，支持向量机目前只适合小批量样本的任务，无法适应百万甚至上亿样本的任务。
总结
支持向量机是一种分类器。它具有良好的学习能力，且学习的结果具有很好的推广性。支持向量机中最流行的算法是SMO算法，该算法可以通过每次只优化2个apha值来加快SVM的训练速度。相对于简化版的SMO，完整版算法不仅大大提高了优化的速度，还使其存在一些进一步提高运行速度的空间。
核方法，可以将一个在低维空间中的非线性问题转换成高维空间下的线性问题来求解。其中，径向基函数是一个常用的度量两个向量距离的核函数。
支持向量机是一个二类分类器。当用支持向量机解决多类问题时，需要额外的方法对其进行扩展。
参考链接
书籍：周志华的《机器学习实战》
视频：吴恩达的《机器学习》
强烈推荐推导看下面几个
魔乐社区（Modelers.cn) 是一个中立、公益的人工智能社区，提供人工智能工具、模型、数据的托管、展示与应用协同服务，为人工智能开发及爱好者搭建开放的学习交流平台。社区通过理事会方式运作，由全产业链共同建设、共同运营、共同享有，推动国产AI生态繁荣发展。
更多推荐
2
0- 0
已为社区贡献3条内容
所有评论(0)