机器学习 支持向量机（SVM） 从理论到放弃，从代码到理解
基本概念
支持向量机（support vector machines，SVM）是一种二分类模型，它的基本模型是定义在特征空间上的间隔最大的线性分类器。支持向量机还包括核技巧，这使它成为实质上的非线性分类器。其学习策略就是间隔最大化，可形式化为一个求解凸二次规划（convex quadratic programming）的问题，也等价于正则化的合页损失函数的最小化问题，支持向量机的学习算法是求解凸二次规划的最优化算法。html
支持向量机学习方法包含构建由简至繁的模型：线性可分支持向量机（linear support vector machine in linearly separable case）、线性支持向量机（linear support vector machine）及非线性支持向量机（non-linear support vector machine）。当训练数据线性可分时，经过硬间隔最大化（hard margin maximization），学习一个线性的分类器；当训练数据近似线性可分时，经过软间隔最大化（soft margin maximization），也学习一个线性的分类器；当训练数据不可分时，经过核技巧（kernel trick）及软间隔最大化，训练非线性的分类器。python
最大边缘超平面
对于一个二分类问题，若数据集是线性可分的，那么咱们能够找到这样一个超平面，使得数据的两个label分别位于平面两侧。而且，能够注意到，这样的超平面咱们能够找到无数个。算法
然而，虽然有无数个这样的超平面，咱们并不能保证这些超平面在未知实例上的效果都同样的好。所以，咱们须要找到一个具备很好的泛化偏差的超平面。这里，咱们能够选取最大边缘超平面。下面将讲述缘由。数组
函数间隔
下图中有A，B，C三点，表示3个实例，对于这3个实例的预测，咱们预测的确信度会是A>B>C这样一个结果。通常来讲，一个点距离超平面的远近能够表示分类预测的确信程度。在超平面 $w \cdot x+b=0$ 肯定的状况下， $|w \cdot x +b|$ 可以相对地表示点 $x$ 距离超平面地远近。而 $w \cdot x + b$ 的符号与相似标记 $y$ 的符号是否一致可以表示分类是否正确。因此可用量 $y(w \cdot x + b)$ 来表示分类的正确性及确信度，这就是函数间隔（function margin）。缓存
几何间隔
函数间隔能够表示分类预测的正确性及确信度，可是选择分离超平面时，只有函数间隔还不够，由于只要成比例地该边$w和b$，超平面并无改变，但函数间隔却变为原来地两倍。所以，咱们须要对分离超平面的法向量 $w$ 加某些约束，如规范化。这时，函数间隔便成了几何间隔（geometric margin）。app
间隔最大化
对训练数据集找到几何间隔最大的超平面意味着以充分大的确信度对训练数据进行分类。即，不只将正负实例点分开，并且对最难分的实例点也有足够大的确信度将它们分开，这样的超平面对未知的新实例有很好的分类预测能力。这也是为何咱们要选取最大边缘超平面的缘由。dom
线性可分支持向量机
最大间隔法
最优化
运用拉格朗日对偶性，经过求解对偶问题（dual problem）获得原始问题的最优解。机器学习
优势：ide
- 对偶问题每每更易求解
- 天然引入核函数，进而推广到非线性分类问题
首先构建拉格朗日函数，对每一个不等式约束引进拉格朗日乘子 $\alpha_i \ge 0,i=1,2, \cdots,N$ ,定义拉格朗日函数： $$ L(w,b,\alpha)=\frac{1}{2}||w||^2-\sum_{i=1}^{N}\alpha_iy_i(w\cdot x_i+b)+\sum_{i=1}^{N}\alpha_i $$ 其中，$\alpha =(\alpha_1,\alpha_2, \cdots ,\alpha_N)^T$ 为拉格朗日乘子向量。函数
根据拉格朗日对偶性，原始问题的对偶问题就是极大极小问题： $$ \max_a\min_{w,b}L(w,b,a) $$ 附上推导过程
线性支持向量机
当数据集没法线性可分时，以下图所示
这时，咱们能够对每个样本点引进一个松弛变量 $\zeta_i \ge 0 $ ，使函数间隔加上松弛变量大于等于1，这样，约束条件就变成 $$ y_i(w \cdot x+b)\ge 1-\zeta_i $$ 同时，对每一个松弛变量 $ \zeta_i $ ，支付一个代价 $\zeta_i $ 。目标函数由原来的$\frac{1}{2}||w||^2$ 变成 $$ \frac{1}{2}||w||^2+C\sum_{i=1}^{N}\zeta_i $$ 这里，$ C > 0 $ 称为惩罚参数，通常由应用问题决定。最小化目标函数包含两层含义：使 $ \frac{1}{2}||w||^2 $ 尽可能小，即间隔尽可能大，同时使偏差分类点的个数尽可能小，C是调和两者的系数。
同上述推导方法，咱们可写出下面算法
非线性支持向量机
有时，咱们面对的是非线性问题，以下图左。非线性问题每每很差求解，因此但愿能用解线性分类问题的方法解决这个问题，采起的方法是进行一个非线性变换，将非线性问题变换为线性问题，经过解变换后的线性问题的方法求解原来的非线性问题。
原空间： $$ \mathcal{X} \subset \mathbf{R}^{2}, x=\left(x^{(1)}, x^{(2)}\right)^{\mathrm{T}} \in \mathcal{X} $$ 新空间： $$ \mathcal{Z} \subset \mathbf{R}^{2}, z=\left(z^{(1)}, z^{(2)}\right)^{\mathrm{T}} \in \mathcal{Z} \quad z=\phi(x)=\left(\left(x^{(1)}\right)^{2},\left(x^{(2)}\right)^{2}\right)^{\mathrm{T}} $$
核函数
设 $\mathcal{X}$ 是输入空间(欧氏空间 $\mathbf{R}^{n}$ 的子集或离散集合)，又设 $\mathcal{H}$ 为特征空间(希尔伯特空间)，若是存在一个从 $\mathcal{X}$ 到 $\mathcal{H}$ 的映射 $$ \phi(x) : \mathcal{X} \rightarrow \mathcal{H} $$ 使得对全部 $$ x, z \in \mathcal{X} $$ 函数K(x,z)知足条件 $$ K(x, z)=\phi(x) \cdot \phi(z) $$ 则称 K(x,z) 为核函数， $\phi(x) $ 为映射函数， 式中 $\phi(x) \cdot \phi(z)$ 为 $\phi(x) $ 和 $\phi(z)$ 的内积。
在学习与预测中只定义核函数K(x,z)，而不显式地定义映射函数，一般，直接计算K(x,z)比较容易，而经过 $\phi(x) $ 和 $\phi(z)$ 计算K(x, z)并不容易。
注：φ是输入空间 $\mathbf{R}^{n}$ 到特征空间 $\mathcal{H}$ 的映射，特征空间 $\mathcal{H}$ 通常是高维，映射能够不一样。
引入核函数后，目标函数改成 $$ W(\alpha)=\frac{1}{2} \sum_{i=1}^{N} \sum_{j=1}^{N} \alpha_{i} \alpha_{j} y_{i} y_{j} K\left(x_{i}, x_{j}\right)-\sum_{i=1}^{N} \alpha_{i} $$
正定核
水平有限，不会推导，就直接给结论好了。。orz
一般咱们所说的核函数就是正定核函数，设 $\mathrm{K} : \mathcal{X} \times \mathcal{X} \rightarrow \mathbf{R}$ 是对称函数，则K(x,z)为正定核函数 的充要条件是对任意 $x_{i} \in \mathcal{X}, \quad i=1,2, \cdots, m$ ，K(x,z)对应的Gram矩阵 $$ K=\left[K\left(x_{i}, x_{j}\right)\right]_{m \times m} $$ 是半正定的。
这必定义在构造核函数时颇有用。但对于一个具体函数 K(x,z) 来讲，检验它是否为正定核函数并不容易，由于 要求对任意有限输入集 $\left{x_{1}, x_{2}, \cdots, x_{m}\right}$ 验证K对应的 Gram矩阵是否为半正定的。在实际问题中每每应用己有的核函数。
经常使用核函数
- 多项式核函数（Polynomial kernel function)
$$ K(x, z)=(x \cdot z+1)^{p} $$
对应的支持向量机为P次多项式分类器，分类决策函数： $$ f(x)=\operatorname{sign}\left(\sum_{i=1}^{N_{s}} a_{i}^{} y_{i}\left(x_{i} \cdot x+1\right)^{p}+b^{}\right) $$
- 高斯核函数 （Gaussian Kernel Function)
$$ K(x, z)=\exp \left(-\frac{|x-z|^{2}}{2 \sigma^{2}}\right) $$
决策函数： $$ f(x)=\operatorname{sign}\left(\sum_{i=1}^{N_{s}} a_{i}^{} y_{i} \exp \left(-\frac{|x-z|^{2}}{2 \sigma^{2}}\right)+b^{}\right) $$
- 字符串核函数（string kernel function）：
核函数不只能够定义在欧式空间上，还能够定义在离散数据的集合上。好比，字符串核是定义在字符串集合上的核函数。字符串核函数在文本分类、信息检索、生物信息学等方面都有应用。 $$ k_n(s,t)=\sum_{u\in \sum^n}[\phi_n(s)]u[\phi_n(t)]u=\sum{u\in \sum^n}\sum{(i,j):s(i)=t(j)=u}\lambda^{l(i)}\lambda^{l(j)} $$
字符串核函数 $k_n(s,t)$ 给出了字符串s和t中长度等于n的全部字串组成的特征向量的余弦类似度。
算法
序列最小优化（SMO）算法
基本思路
若是全部变量的解都知足此最优化问题的KKT条件，那么获得解； 不然，选择两个变量，固定其它变量，针对这两个变量构建一个二次规划问题，称为子问题，可经过解析方法求解，提升计算速度。
子问题的两个变量：一个是违反KKT条件最严重的那个，另 一个由约束条件自动肯定。假设 $\alpha_1,\alpha_2$ 为两个变量， $\alpha_3,\alpha_4, \cdots , \alpha_ N$ 固定，那么由等式约束可知 $$ \alpha_1=-y_1\sum_{i=2}^{N}\alpha_iy_i $$
变量选择
SMO算法在每一个子问题中选择两个变量优化，其中至少一个变量是违反KKT条件的。
- 第一个变量的选择：外循环，违反KKT最严重的样本点，检验样本点是否知足KKT条件：
$$ \alpha_{i}=0 \Leftrightarrow y_{i} g\left(x_{i}\right) \geqslant 1 \ \begin{aligned} 0<\alpha_{i} &<C \Leftrightarrow y_{i} g\left(x_{i}\right)=1 \ \alpha_{i} &=C \Leftrightarrow y_{i} g\left(x_{i}\right) \leqslant 1 \ g\left(x_{i}\right) &=\sum_{j=1}^{N} \alpha_{j} y_{j} K\left(x_{i}, x_{j}\right)+b \end{aligned} $$
-
第二个变量的检查: 内循环， 选择的标准是但愿能使目标函数有足够大的变化，即对应 $|E_1-E_2|$ 最大 ,若是内循环经过上述方法找到的点不能使目标函数有足够的降低，则：遍历间隔边界上的样本点，测试目标函数降低 ，若是降低不大，则遍历全部样本点 ，若是依然降低不大，则丢弃外循环点，从新选择。
-
每次完成两个变量的优化后，从新计算 $b,E_i$
算法
代码
这是我用的数据集数据的分布
感受仍是看代码理解的快，先贴一份根据《机器学习实战》上写的代码
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import random
import math
from sklearn.model_selection import train_test_split
from sklearn.svm import SVC
"""
Author:
ITryagain
Modify:
2019-04-09
"""
def loadDataSet(fileName):
"""
读取数据
Parameters:
fileName - 文件名
Returns:
dataMat - 数据矩阵
labelMat - 数据标签
"""
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 np.array(dataMat), np.array(labelMat)
class optStruct:
def __init__(self, dataMatIn, classLabels, C, toler, kTup):
"""
:param dataMatIn: 输入数据 [X1, X2, ... , XN]
:param classLabels: 分类标签 [y]
:param C: 松弛变量
:param toler: 容错率
"""
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)
def calcEk(oS, k):
"""
:param oS:
:param k: 第k行
:return: 预测值与实际值的差
"""
fXk = float(np.multiply(oS.alphas, oS.labelMat).T * oS.k[:, k]) + oS.b
Ek = fXk - float(oS.labelMat[k])
return Ek
def clipAlpha(aj, H, L):
"""
调整aj的值，使aj处于 L<=aj<=H
:param aj: 目标值
:param H: 最大值
:param L: 最小值
:return:
aj 目标值
"""
if aj > H:
aj = H
elif L > aj:
aj = L
return aj
def selectJrand(i, m):
"""
随机选择一个整数
:param i: 第一个alpha的下标
:param m: 全部alpha的数目
:return:
j 返回一个不为i的随机数，在0~m之间的整数值
"""
j = i
while j == i:
j = int(random.uniform(0, m))
return j
def selectJ(i, oS, Ei):
"""
内循环的启发式方法，选择第二个alpha的值
:param i: 第一个alpha的下标
:param oS:
:param Ei: 预测结果与真实结果比对，计算偏差Ei
:return:
j 随机选出的第j一行
Ej 预测结果与真实结果比对，计算偏差Ej
"""
maxK = -1
maxDeltaE = 0
Ej = 0
oS.eCache[i] = [1, Ei] # 将输入值Ei在缓存中设置成为有效的。这里的有效意味着它已经计算好了。
validEcacheList = np.nonzero(oS.eCache[:, 0].A)[0] # 非零E值的行的list列表，所对应的alpha值
if(len(validEcacheList)) > 1:
for k in validEcacheList: # 在全部的值上进行循环，并选择其中使得改变最大的那个值
if k ==1:
continue
Ek = calcEk(oS, k)
deltaE = abs(Ei - Ek)
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):
"""
计算偏差值并存入缓存中
:param oS:
:param k: 某一列的行号
"""
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]
def kernelTrans(X, A, kTup):
"""
核转换函数
:param X: dataMatIn数据集
:param A: dataMatIn数据集的第i行的数据
:param kTup:核函数的信息
:return:
"""
m, n = np.shape(X)
K = np.mat(np.zeros((m, 1)))
if kTup[0] == 'lin':
# linear kernel: m*n * n*1 = m*1
K = X * A.T
elif kTup[0] == 'rbf':
for j in range(m):
deltaRow = X[j, :] - A
K[j] = deltaRow * deltaRow.T
# 径向基函数的高斯版本
# print(len(K))
# K = math.exp(K / (-1 * kTup[1] ** 2))
for i in range(m):
K[i] = math.exp(K[i] / (-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
def innerL(i, oS):
"""
内循环代码
:param i:具体的某一行
:param oS:
:returns:
0 找不到最优的值
1 找到了最优的值，而且oS.Cache到缓存中
"""
Ei = calcEk(oS, i)
# 约束条件 (KKT条件是解决最优化问题的时用到的一种方法。咱们这里提到的最优化问题一般是指对于给定的某一函数，求其在指定做用域上的全局最小值)
# 0<=alphas[i]<=C，但因为0和C是边界值，咱们没法进行优化，由于须要增长一个alphas和下降一个alphas。
# 表示发生错误的几率：labelMat[i]*Ei 若是超出了 toler， 才须要优化。至于正负号，咱们考虑绝对值就对了。
'''
# 检验训练样本(xi, yi)是否知足KKT条件
yi*f(i) >= 1 and alpha = 0 (outside the boundary)
yi*f(i) == 1 and 0<alpha< C (on the boundary)
yi*f(i) <= 1 and alpha = C (between the boundary)
'''
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进行优化。效果更明显
j ,Ej = selectJ(i, oS, Ei)
alphaIold = oS.alphas[i].copy()
alphaJold = oS.alphas[j].copy()
# L和H用于将alphas[j]调整到0-C之间。若是L==H，就不作任何改变，直接return 0
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是alphas[j]的最优修改量，若是eta==0，须要退出for循环的当前迭代过程
# 参考《统计学习方法》李航-P125~P128<序列最小最优化算法>
eta = 2.0 * oS.k[i,j] - oS.k[i, i] - oS.k[j, j]
if eta >= 0:
print("eta>=0")
return 0
# 计算出一个新的alphas[j]值
oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej) / eta
# 并使用辅助函数，以及L和H对其进行调整
oS.alphas[j] = clipAlpha(oS.alphas[j], H, L)
# 更新偏差缓存
updateEk(oS, j)
# 检查alpha[j]是否只是轻微的改变，若是是的话，就退出for循环。
if abs(oS.alphas[j] - alphaJold) < 0.00001:
print("j not moving enough")
return 0
# 而后alphas[i]和alphas[j]一样进行改变，虽然改变的大小同样，可是改变的方向正好相反
oS.alphas[i] += oS.labelMat[j] * oS.labelMat[i] * (alphaJold - oS.alphas[j])
# 更新偏差缓存
updateEk(oS, i)
# 在对alpha[i], alpha[j] 进行优化以后，给这两个alpha值设置一个常数b。
# w= Σ[1~n] ai*yi*xi => b = yj Σ[1~n] ai*yi(xi*xj)
# 因此： b1 - b = (y1-y) - Σ[1~n] yi*(a1-a)*(xi*x1)
# 为何减2遍？ 由于是 减去Σ[1~n]，正好2个变量i和j，因此减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[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 smoP(dataMatIn, classLabels, C, toler, maxIter, kTup = ('lin', 0)):
"""
SMO算法外循环
:param dataMatIn: 输入数据集
:param classLabels: 标签
:param C: 松弛变量(常量值)，容许有些数据点能够处于分隔面的错误一侧。
控制最大化间隔和保证大部分的函数间隔小于1.0这两个目标的权重。
能够经过调节该参数达到不一样的结果。
:param toler: 容错率
:param maxIter: 退出前最大的循环次数
:param kTup: 核函数
:return:
"""
# 建立一个 optStruct 对象
oS = optStruct(np.mat(dataMatIn), np.mat(classLabels).transpose(), C, toler, kTup)
iter = 0
entireSet = True
alphaPairsChanged = 0
# 循环遍历：循环maxIter次 而且 （alphaPairsChanged存在能够改变 or 全部行遍历一遍）
# 循环迭代结束 或者 循环遍历全部alpha后，alphaPairs仍是没变化
while iter < maxIter and (alphaPairsChanged > 0 or entireSet):
alphaPairsChanged = 0
# 当entireSet=true or 非边界alpha对没有了；就开始寻找 alpha对，而后决定是否要进行else。
if entireSet:
# 在数据集上遍历全部可能的alpha
for i in range(oS.m):
# 是否存在alpha对，存在就+1
alphaPairsChanged += innerL(i, oS)
print("fullSet, iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
iter += 1
# 对已存在 alpha对，选出非边界的alpha值，进行优化。
else:
# 遍历全部的非边界alpha值，也就是不在边界0或C上的值。
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
# 若是找到alpha对，就优化非边界alpha值，不然，就从新进行寻找，若是寻找一遍 遍历全部的行仍是没找到，就退出循环。
if entireSet:
entireSet = False
elif alphaPairsChanged == 0:
entireSet = True
print("iteration number: %d" % iter)
return oS
def calcWs(alphas, dataArr, classLabels):
"""
基于alpha计算w值
:param alphas: 拉格朗日乘子
:param dataArr: 数据集
:param classLabels:
:return:
wc 回归系数
"""
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
def plotfig_SVM(xArr, yArr, ws, b, alphas):
"""
参考地址：
http://blog.csdn.net/maoersong/article/details/24315633
http://www.cnblogs.com/JustForCS/p/5283489.html
http://blog.csdn.net/kkxgx/article/details/6951959
"""
xMat = np.mat(xArr)
yMat = np.mat(yArr)
# 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最大值，最小值根据原数据集dataArr[:, 0]的大小而定
x = np.arange(-1.0, 10.0, 0.1)
# 根据x.w + b = 0 获得，其式子展开为w0.x1 + w1.x2 + b = 0, 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.plot(xMat[i, 0], xMat[i, 1], 'cx')
else:
ax.plot(xMat[i, 0], xMat[i, 1], 'kp')
# 找到支持向量，并在图中标红
for i in range(70):
if alphas[i] > 0.0:
ax.plot(xMat[i, 0], xMat[i, 1], 'ro')
plt.show()
def predict(data, oS):
r = oS.b
for i in range(oS.m):
r += oS.alphas[i] * oS.labelMat[i] * data * oS.X[i, :].T
return 1 if r > 0 else -1
def score(X_test, y_test, oS):
right_count = 0
for i in range(len(X_test)):
result = predict(X_test[i], oS)
if result == y_test[i]:
right_count += 1
return right_count / len(X_test)
if __name__ == "__main__":
# 获取特征和目标变量
dataArr, labelArr = loadDataSet('testSetRBF2.txt')
X_train, X_test, y_train, y_test = train_test_split(dataArr, labelArr, test_size=0.3, random_state=4)
# print labelArr
# b是常量值， alphas是拉格朗日乘子
# 0.6 0.001 40 0.5666666666666667
oS = smoP(X_train, y_train, 0.6, 0.001, 200, kTup=('rbf', 10))
b = oS.b
alphas = oS.alphas
svInd = np.nonzero(alphas.A > 0)[0]
sVs = np.mat(X_train)[svInd]
labelSV = np.mat(y_train).transpose()[svInd]
print('/n/n/n')
print('b=', b)
print('alphas[alphas>0]=', alphas[alphas > 0])
print('shape(alphas[alphas > 0])=', np.shape(alphas[alphas > 0]))
for i in range(70):
if alphas[i] > 0:
print(dataArr[i], labelArr[i])
# 画图
ws = calcWs(alphas, X_train, y_train)
plotfig_SVM(X_train, y_train, ws, b, alphas)
# print(score(X_test, y_test, oS))
datMat = np.mat(X_test)
lableMat = np.mat(y_test).transpose()
m, n = np.shape(datMat)
right_count = 0
for i in range(m):
kernelEval = kernelTrans(sVs, datMat[i, :], ('rbf', 10))
pre = kernelEval.T * np.multiply(labelSV, alphas[svInd]) + b
sign = 1 if pre > 0 else -1
if sign == labelArr[i]:
right_count += 1
print(right_count / len(X_test))
emmmmmm。。。可能看起来仍是一脸懵逼？？？不要紧，再看看这一份，两份一块儿看能够更好的理解哦
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import random
from sklearn.model_selection import train_test_split
from sklearn.svm import SVC
"""
Author:
ITryagain
Modify:
2019-04-09
"""
def loadDataSet(fileName):
"""
读取数据
Parameters:
fileName - 文件名
Returns:
dataMat - 数据矩阵
labelMat - 数据标签
"""
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 np.array(dataMat), np.array(labelMat)
class SVM:
def __init__(self, max_iter=100, kernel='linear'):
self.max_iter = max_iter
self._kernel = kernel
def init_args(self, features, labels):
self.m, self.n = features.shape
self.X = features
self.Y = labels
self.b = 0.0
# 将Ei保存在一个列表里
self.alpha = np.ones(self.m)
self.E = [self._E(i) for i in range(self.m)]
# 松弛变量
self.C = 1.0
def _KKT(self, i):
y_g = self._g(i) * self.Y[i]
if self.alpha[i] == 0:
return y_g >= 1
elif 0 < self.alpha[i] < self.C:
return y_g == 1
else:
return y_g <= 1
# g(x)预测值，输入xi（X[i]）
def _g(self, i):
r = self.b
for j in range(self.m):
r += self.alpha[j] * self.Y[j] * self.kernel(self.X[i], self.X[j])
return r
# 核函数
def kernel(self, x1, x2):
if self._kernel == 'linear':
return sum([x1[k] * x2[k] for k in range(self.n)])
elif self._kernel == 'poly':
return (sum([x1[k] * x2[k] for k in range(self.n)]) + 1) ** 2
return 0
# E（x）为g(x)对输入x的预测值和y的差
def _E(self, i):
return self._g(i) - self.Y[i]
def _init_alpha(self):
# 外层循环首先遍历全部知足0<a<C的样本点，检验是否知足KKT
index_list = [i for i in range(self.m) if 0 < self.alpha[i] < self.C]
# 不然遍历整个训练集
non_satisfy_list = [i for i in range(self.m) if i not in index_list]
index_list.extend(non_satisfy_list)
for i in index_list:
if self._KKT(i):
continue
E1 = self.E[i]
# 若是E2是+，选择最小的；若是E2是负的，选择最大的
if E1 >= 0:
j = min(range(self.m), key=lambda x: self.E[x])
else:
j = max(range(self.m), key=lambda x: self.E[x])
return i, j
def _compare(self, _alpha, L, H):
if _alpha > H:
return H
elif _alpha < L:
return L
else:
return _alpha
def fit(self, features, labels):
self.init_args(features, labels)
for t in range(self.max_iter):
# train
i1, i2 = self._init_alpha()
# 边界
if self.Y[i1] == self.Y[i2]:
L = max(0, self.alpha[i1] + self.alpha[i2] - self.C)
H = min(self.C, self.alpha[i1] + self.alpha[i2])
else:
L = max(0, self.alpha[i2] - self.alpha[i1])
H = min(self.C, self.C + self.alpha[i2] - self.alpha[i1])
E1 = self.E[i1]
E2 = self.E[i2]
# eta=K11+K22-2K12
eta = self.kernel(self.X[i1], self.X[i1]) + self.kernel(self.X[i2], self.X[i2]) - 2 * self.kernel(
self.X[i1], self.X[i2])
if eta <= 0:
# print('eta <= 0')
continue
alpha2_new_unc = self.alpha[i2] + self.Y[i2] * (E1 - E2) / eta # 此处有修改，根据书上应该是E1 - E2，书上130-131页
alpha2_new = self._compare(alpha2_new_unc, L, H)
alpha1_new = self.alpha[i1] + self.Y[i1] * self.Y[i2] * (self.alpha[i2] - alpha2_new)
b1_new = -E1 - self.Y[i1] * self.kernel(self.X[i1], self.X[i1]) * (alpha1_new - self.alpha[i1]) - self.Y[
i2] * self.kernel(self.X[i2], self.X[i1]) * (alpha2_new - self.alpha[i2]) + self.b
b2_new = -E2 - self.Y[i1] * self.kernel(self.X[i1], self.X[i2]) * (alpha1_new - self.alpha[i1]) - self.Y[
i2] * self.kernel(self.X[i2], self.X[i2]) * (alpha2_new - self.alpha[i2]) + self.b
if 0 < alpha1_new < self.C:
b_new = b1_new
elif 0 < alpha2_new < self.C:
b_new = b2_new
else:
# 选择中点
b_new = (b1_new + b2_new) / 2
# 更新参数
self.alpha[i1] = alpha1_new
self.alpha[i2] = alpha2_new
self.b = b_new
self.E[i1] = self._E(i1)
self.E[i2] = self._E(i2)
return 'train done!'
def predict(self, data):
r = self.b
for i in range(self.m):
r += self.alpha[i] * self.Y[i] * self.kernel(data, self.X[i])
return 1 if r > 0 else -1
def score(self, X_test, y_test):
right_count = 0
for i in range(len(X_test)):
result = self.predict(X_test[i])
if result == y_test[i]:
right_count += 1
return right_count / len(X_test)
def _weight(self):
# linear model
yx = self.Y.reshape(-1, 1) * self.X
self.w = np.dot(yx.T, self.alpha)
return self.w
if __name__ == '__main__':
X, y = loadDataSet('testSetRBF2.txt')
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=4)
svm = SVM(max_iter=200, kernel='poly')
svm.fit(X, y)
print('the score = {}'.format(svm.score(X_test, y_test)))
clf = SVC(gamma='auto')
clf.fit(X_train, y_train)
print(clf.score(X_test, y_test))
# x_ponits = np.arange(4, 8)
# y_ = -(lr_clf.weights[1] * x_ponits + lr_clf.weights[0]) / lr_clf.weights[2]
# plt.plot(x_ponits, y_)
#
# # lr_clf.show_graph()
data_plus = []
data_minus = []
for i in range(len(X)):
if y[i] > 0:
data_plus.append(X[i])
else:
data_minus.append(X[i])
data_plus_np = np.array(data_plus) # 转换为numpy矩阵
data_minus_np = np.array(data_minus) # 转换为numpy矩阵
plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1]) # 正样本散点图
plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1]) # 负样本散点图
plt.show()
输出
the score = 0.7 0.8333333333333334
小节
推了三四天的SVM。。。感受本身快死了，若是不是手头有个比赛须要修改SVM，打死我都不会去推这公式，下次仍是从代码入手好了，数学很差的痛苦orz。。。。
- 1. 机器学习 支持向量机（SVM） 从理论到放弃，从代码到理解
- 2. 支持向量机（SVM）从入门到放弃再到掌握
- 3. 机器学习：支持向量机(svm)
- 4. 机器学习：支持向量机(SVM)
- 5. 机器学习：支持向量机SVM
- 6. 机器学习——支持向量机（SVM）
- 7. 支持向量机SVM-----机器学习
- 8. 【机器学习】支持向量机(SVM)
- 9. 机器学习----SVM支持向量机
- 10. 机器学习-支持向量机SVM
- 更多相关文章...
- • Git 分支管理 - Git 教程
- • CAP理论是什么？ - NoSQL教程
- • 漫谈MySQL的锁机制
- • Docker 清理命令
-
每一个你不满意的现在，都有一个你没有努力的曾经。
- 1. vs2019运行opencv图片显示代码时，窗口乱码
- 2. app自动化 - 元素定位不到？别慌，看完你就能解决
- 3. 在Win8下用cisco ××× Client连接时报Reason 422错误的解决方法
- 4. eclipse快速补全代码
- 5. Eclipse中Java/Html/Css/Jsp/JavaScript等代码的格式化
- 6. idea+spring boot +mabitys(wanglezapin)+mysql (1)
- 7. 勒索病毒发生变种 新文件名将带有“.UIWIX”后缀
- 8. 【原创】Python 源文件编码解读
- 9. iOS9企业部署分发问题深入了解与解决
- 10. 安装pytorch报错CondaHTTPError：******