CS_ALearn/ailearning
1
0
Fork 0
mirror of https://github.com/apachecn/ailearning.git synced 2026-11-06 07:17:43 +08:00
Code Issues Projects Releases Wiki Activity

Merge pull request #83 from jiangzhonglian/master

Browse Source 
添加6.SVM视频中2个问题的注释
This commit is contained in:
片刻
2017-04-20 00:33:45 +08:00
committed by GitHub
parent c25cad8f09 11a3822879
commit 7463424061
7 changed files with 306 additions and 154 deletions
Download Patch File Download Diff File Expand all files Collapse all files

42
docs/15.大数据与MapReduce.md
View File

@@ -0,0 +1,42 @@
# 大数据与MapReduce
> 本章内容
* MapReduce
* Python中Hadoop流的使用
* 使用mrjob库将MapReduce自动化
* 利用Pegasos算法并行训练支持向量机
## MapReduce：分布式计算的框架
```
优点：可在短时间内完成大量工作。
缺点：算法必须经过重写，需要对系统工程有一定的理解。
适用数据类型：数值型和标称型数据。
```
* MapReduce集群的示意图
![MapReduce集群的示意图](/images/15.BigData_MapReduce/MR_1_cluster.jpg)
> 关于MapRduce的学习要点
* 主节点控制MapReduce的作业流程
* MapReduce的作业可以分成map任务和reduce任务
* map任务之间不做数据交流，reduce任务也一样
* 在map和reduce阶段中间，有一个sort和combine阶段
* 数据被重复存放在不同的机器上，以防止某个机器实效
* mapper和reducer传输的数据形式为key/value对
## Hadoop流
## 在Amazon网络服务商运行Hadoop程序
## MapReduce上的机器学习
## 在Python中使用mrjob来自动化MapReduce
## 示例：分布式SVM的Pegasos算法
## 你真的需要MapReduce吗？
## 本章小节

51
docs/6.支持向量机.md
View File

@@ -11,7 +11,7 @@
* 机(Machine)就是表示一种算法，而不是表示机器。
* SVM有很多种实现，最流行的一种实现是： `序列最小优化(Sequential Minimal Optimization, SMO)算法`。
* 下面还会介绍一种称为`核函数(kernel)`的方式将SVM扩展到更多数据集上。
* 实战项目：回顾第1章中手写识别的案例，并考察其能否通过SVM来提高识别的效果。
* 实战项目：回顾第2章中手写识别的案例，并考察其能否通过SVM来提高识别的效果。
* 注意：`SVM几何含义比较直观，但其算法实现较复杂，牵扯大量数学公式的推导。`
```
@@ -53,55 +53,56 @@ This is the simplest kind of SVM (Called an LSVM) Support Vectors are those data
```
* 选择D会比B、C分隔的效果要好很多，原因是上述的5个结论。
* 所有的点看作地雷吧，那么我们(超平面)得找到最近所有的地雷，并保证我们离它最远。
![线性可分](/images/6.SVM/SVM_3_linearly-separable.jpg)
### 怎么寻找最大间隔
> 点到超平面的距离
* 分隔超平面`函数间距`: $$y(x)=w^Tx+b$$
* 分类的结果： $$f(x)=sign(w^Tx+b)$$ (sign表示>0为1，<0为-1，=0为0)
* 点到超平面的`几何间距`: $$d(x)=(w^Tx+b)/||w||$$ [ ||w||表示w矩阵的二范式=> sqrt(w*w^T), 点到超平面的距离也是类似的，需要复习一下向量的知识 ]
* 分隔超平面`函数间距`: \\(y(x)=w^Tx+b\\)
* 分类的结果： \\(f(x)=sign(w^Tx+b)\\) (sign表示>0为1，<0为-1，=0为0)
* 点到超平面的`几何间距`: \\(d(x)=(w^Tx+b)/||w||\\) （||w||表示w矩阵的二范式=> \\(\sqrt{w*w^T}\\), 点到超平面的距离也是类似的，需要复习一下向量的知识）
![点到直线的几何距离](/images/6.SVM/SVM_4_point2line-distance.png)
> 拉格朗日乘子法的使用
* 类别标签用-1、1，是为了后期方便`lable*(w^Tx+b)`的标识和距离计算；如果`lable*(w^Tx+b)>0`表示预测正确，否则预测错误。
* 类别标签用-1、1，是为了后期方便 \\(lable*(w^Tx+b)\\) 的标识和距离计算；如果 \\(lable*(w^Tx+b)>0\\) 表示预测正确，否则预测错误。
* 现在目标很明确，就是要找到`w`和`b`，因此我们必须要找到最小间隔的数据点，也就是前面所说的`支持向量`。
* 也就说，让最小的距离取最大.(最小的距离：就是最小间隔的数据点；最大：就是最大间距，为了找出最优超平面--最终就是支持向量)
* 怎么理解呢？ 例如： 所有的点看作地雷吧，那么我们(超平面)得找到最近所有的地雷，并保证我们离它最远。
* 目标函数：`arg max{min(lable*(w^Tx+b)/||w||)}`
* 1.如果`lable*(w^Tx+b)>0`表示预测正确，`||w||`可以理解为归一化，我们始终可以找到一个阈值让`lable*(w^Tx+b)>=1`
* 2.所以令`lable*(w^Tx+b)=1`，我们本质上是求 `arg max[关于w, b] (1/||w||)`；也就说，我们约束(前提)条件是`lable*(w^Tx+b)>=1`
* 新的目标函数求解： `arg max[关于w, b] (1/||w||)`
* => 就是求: `arg min[关于w, b] (||w||)` (求矩阵会比较麻烦，如果x只是1/2*X^2的偏导数，那么。。同样是求最小值)
* => 就是求: `arg min[关于w, b] (1/2*||w||^2)` (二次函数求导，求极值，平方也方便计算)
* 目标函数：\\(arg: max\ \{min\ [lable*(w^Tx+b)/||w||]\}\\)
* 1.如果 \\(lable*(w^Tx+b)>0\\) 表示预测正确，也称`函数间隔`，\\(||w||\\) 可以理解为归一化，也称`几何间隔`，我们始终可以找到一个阈值让 \\(lable*(w^Tx+b)>=1\\)
* 2.所以令 \\(lable*(w^Tx+b)=1\\)，我们本质上是求 \\(arg: max\{关于w, b\}\ (1/||w||)\\)；也就说，我们约束(前提)条件是: \\(lable*(w^Tx+b)=1\\)
* 新的目标函数求解： \\(arg: max\{关于w, b\}\ (1/||w||)\\)
* => 就是求: \\(arg: min\{关于w, b\}\ (||w||)\\) (求矩阵会比较麻烦，如果x只是1/2*X^2的偏导数，那么。。同样是求最小值)
* => 就是求: \\(arg: min\{关于w, b\}\ (\frac{1}{2}*||w||^2)\\) (二次函数求导，求极值，平方也方便计算)
* 本质上就是求线性不等式的二次优化问题(求分隔超平面，等价于求解相应的凸二次规划问题。)
* 通过拉格朗日乘子法，求二次优化问题
* 假设需要求极值的目标函数 (objective function) 为 f(x,y)，限制条件为 φ(x,y)=M
* 设g(x,y)=M-φ(x,y) # 临时φ(x,y)表示下文中`label*(w^Tx+b)`
* 假设需要求极值的目标函数 (objective function) 为 f(x,y)，限制条件为 φ(x,y)=M # M=1
* 设g(x,y)=M-φ(x,y) # 临时φ(x,y)表示下文中 \\(label*(w^Tx+b)\\)
* 定义一个新函数: F(x,y,λ)=f(x,y)+λg(x,y)
* a为λ，代表要引入的拉格朗日乘子(Lagrange multiplier)
* 那么： `L(w,b,a)=1/2||w||^2 + Σ a*[1-label*(w^Tx+b)]` => `L(w,b,a)=1/2||w||^2 - Σ a*[label*(w^Tx+b)-1]`
* 因为：`label*(w^Tx+b)>=1, a>0` , 所以`a*[label*(w^Tx+b)-1]>=0`, Σ 0~n的加和也都大于0
* `max[关于a] (L(w,b,a)) = 1/2*||w||^2`
* 如果求： `min(1/2*||w||^2)`, 也就是要求： `arg min[关于w, b]{ max[关于a] (L(w,b,a)) }`
* 那么： \\(L(w,b,\alpha)=\frac{1}{2} * ||w||^2 + \sum_{i=1}^{n} \alpha * [1 - label * (w^Tx+b)]\\)
* 因为：\\(label*(w^Tx+b)>=1, \alpha>0\\) , 所以 \\(\alpha*[label*(w^Tx+b)-1]>=0\\), \\(\sum_{i=1}^{n} \alpha * [label * (w^Tx+b) - 1]>=0\\)
* \\(max\{关于\alpha\}\ L(w,b,\alpha) = \frac{1}{2} *||w||^2\\)
* 如果求： \\(min\{关于w, b\}\ \frac{1}{2} *||w||^2\\) , 也就是要求： \\(min\{关于w, b\}\ [max\{关于\alpha\}\ L(w,b,\alpha)]\\)
* 现在转化到对偶问题的求解
* `min[关于w, b]{ max[关于w, b] (L(w,b,a)) } >= max[关于w, b] { min[关于w, b] (L(w,b,a)) }`
* \\(min\{关于w, b\}\ [max\{关于\alpha\}\ L(w,b,\alpha)]\\) >= \\(max\{关于\alpha\}\ [min\{关于w, b\}\ L(w,b,\alpha)]\\)
* 现在分2步
* 先求：`min[关于w, b] (L(w,b,a))=1/2||w||^2 + Σ a*[1-label*(w^Tx+b)]`
* 先求： \\(min\{关于w, b\}\ L(w,b,\alpha)=\frac{1}{2} * ||w||^2 + \sum_{i=1}^{n} \alpha * [1 - label * (w^Tx+b)]\\)
* 就是求`L(w,b,a)`关于[w, b]的偏导数, 得到`w和b的值`，并化简为：`L和a的方程`。
* 参考： 1.【计算拉格朗日函数的对偶函数-下图是参考小象学院】 2.如果公式推导还是不懂，也可以参考《统计学习方法》李航-P103<学习的对偶算法>
* 参考： 如果公式推导还是不懂，也可以参考《统计学习方法》李航-P103<学习的对偶算法>
* ![计算拉格朗日函数的对偶函数](/images/6.SVM/SVM_5_Lagrangemultiplier.png)
* 终于得到课本上的公式：`max(Σ(i=>n) a[i] - 1/2*Σ(i,j=>n) a[i]*a[j]]*label(i)*label(j)*<xi*xj内积>)`
* 约束条件： `a>=0 and Σ a[i]*label(i)=0`
* 终于得到课本上的公式： \\(max\ \{\alpha\}\ [\sum_{i=1}^{m} \alpha - \frac{1}{2} \sum_{i, j=1}^{m} label_i·label_j·\alpha_i·\alpha_j·<x_i, x_j>]\\)
* 约束条件： \\(a>=0,\ and\ \sum_{i=1}^{m} a_i·label_i=0\\)
> 松弛变量(slack variable)
* 我们知道几乎所有的数据都不那么干净, 通过引入松弛变量来允许数据点可以处于分隔面错误的一侧。
* 约束条件： `C>=a>=0 and Σ a[i]*label(i)=0`
* 约束条件： \\(C>=a>=0,\ and\ \sum_{i=1}^{m} a_i·label_i=0\\)
* 这里常量C用于控制“最大化间隔”和“保证大部分点的函数间隔小于1.0” 这两个目标的权重。
* 常量C是一个常数，我们通过调节钙参数得到不同的结果。一旦求出了所有的alpha，那么分隔超平面就可以通过这些alpha来表示。
* 常量C是一个常数，我们通过调节该参数得到不同的结果。一旦求出了所有的alpha，那么分隔超平面就可以通过这些alpha来表示。
* 这一结论十分直接，SVM中的主要工作就是要求解 alpha.
> SVM应用的一般框架
@@ -130,7 +131,7 @@ SVM的一般流程
* 这里指的合适必须要符合一定的条件
* 1.这两个alpha必须要在间隔边界之外
* 2.这两个alpha还没有进行过区间化处理或者不在边界上。
* 之所以要同时改变2个alpha；原因，我们有一个约束条件：`Σ a[i]*label(i)=0`；如果只是修改一个alpha，很可能导致约束条件失效。
* 之所以要同时改变2个alpha；原因，我们有一个约束条件： \\(\sum_{i=1}^{m} a_i·label_i=0\\)；如果只是修改一个alpha，很可能导致约束条件失效。
```
SMO伪代码大致如下：

BIN
images/15.BigData_MapReduce/MR_1_cluster.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 68 KiB

BIN
images/6.SVM/SVM_5_Lagrangemultiplier.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 80 KiB

After

Width:  |  Height:  |  Size: 102 KiB

12
src/python/6.SVM/svm-complete.py
View File

@@ -101,7 +101,7 @@ def calcEk(oS, k):
k 具体的某一行
Returns:
Ek 预测结果与真实结果比对，计算误差Ek
"""
fXk = float(multiply(oS.alphas, oS.labelMat).T * oS.K[:, k] + oS.b)
Ek = fXk - float(oS.labelMat[k])
@@ -142,7 +142,7 @@ def selectJ(i, oS, Ei): # this is the second choice -heurstic, and calcs Ej
maxK = -1
maxDeltaE = 0
Ej = 0
# # 首先将输入值Ei在缓存中设置成为有效的。这里的有效意味着它已经计算好了。
# 首先将输入值Ei在缓存中设置成为有效的。这里的有效意味着它已经计算好了。
oS.eCache[i] = [1, Ei]
# print 'oS.eCache[%s]=%s' % (i, oS.eCache[i])
@@ -247,8 +247,7 @@ def innerL(i, oS):
return 0
# eta是alphas[j]的最优修改量，如果eta==0，需要退出for循环的当前迭代过程
# 如果ETA为0，那么计算新的alphas[j]就比较麻烦了, 为什么呢？ 因为2个值一样。
# 2ab <= a^2 + b^2
# 参考《统计学习方法》李航-P125~P128<序列最小最优化算法>
eta = 2.0 * oS.K[i, j] - oS.K[i, i] - oS.K[j, j] # changed for kernel
if eta >= 0:
print("eta>=0")
@@ -274,6 +273,7 @@ def innerL(i, oS):
# 在对alpha[i], alpha[j] 进行优化之后，给这两个alpha值设置一个常数b。
# w= Σ[1~n] ai*yi*xi => b = yi- Σ[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]):
@@ -515,8 +515,8 @@ if __name__ == "__main__":
# plotfig_SVM(dataArr, labelArr, ws, b, alphas)
# # 有核函数的测试
# testRbf(1)
testRbf(0.8)
# 项目实战
# 示例：手写识别问题回顾
testDigits(('rbf', 20))
# testDigits(('rbf', 20))

348
src/python/6.SVM/svm-complete_Non-Kernel.py
View File

@@ -11,6 +11,18 @@ from numpy import *
import matplotlib.pyplot as plt
class optStruct:
def __init__(self, dataMatIn, classLabels, C, toler): # Initialize the structure with the parameters
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))) # first column is valid flag
def loadDataSet(fileName):
"""loadDataSet（对文件进行逐行解析，从而得到第行的类标签和整个数据矩阵）
@@ -61,6 +73,222 @@ def clipAlpha(aj, H, L):
return aj
def calcEk(oS, k):
"""calcEk（求 Ek误差：预测值-真实值的差）
该过程在完整版的SMO算法中陪出现次数较多，因此将其单独作为一个方法
Args:
oS optStruct对象
k 具体的某一行
Returns:
Ek 预测结果与真实结果比对，计算误差Ek
"""
fXk = float(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): # this is the second choice -heurstic, and calcs Ej
"""selectJ（返回最优的j和Ej）
内循环的启发式方法。
选择第二个(内循环)alpha的alpha值
这里的目标是选择合适的第二个alpha值以保证每次优化中采用最大步长。
该函数的误差与第一个alpha值Ei和下标i有关。
Args:
i 具体的第i一行
oS optStruct对象
Ei 预测结果与真实结果比对，计算误差Ei
Returns:
j 随机选出的第j一行
Ej 预测结果与真实结果比对，计算误差Ej
"""
maxK = -1
maxDeltaE = 0
Ej = 0
# 首先将输入值Ei在缓存中设置成为有效的。这里的有效意味着它已经计算好了。
oS.eCache[i] = [1, Ei]
# print 'oS.eCache[%s]=%s' % (i, oS.eCache[i])
# print 'oS.eCache[:, 0].A=%s' % oS.eCache[:, 0].A.T
# """
# # 返回非0的：行列值
# nonzero(oS.eCache[:, 0].A)= (
# 行： array([ 0, 2, 4, 5, 8, 10, 17, 18, 20, 21, 23, 25, 26, 29, 30, 39, 46,52, 54, 55, 62, 69, 70, 76, 79, 82, 94, 97]),
# 列： array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,0, 0, 0, 0, 0])
# )
# """
# print 'nonzero(oS.eCache[:, 0].A)=', nonzero(oS.eCache[:, 0].A)
# # 取行的list
# print 'nonzero(oS.eCache[:, 0].A)[0]=', nonzero(oS.eCache[:, 0].A)[0]
# 非零E值的行的list列表，所对应的alpha值
validEcacheList = nonzero(oS.eCache[:, 0].A)[0]
if (len(validEcacheList)) > 1:
for k in validEcacheList: # 在所有的值上进行循环，并选择其中使得改变最大的那个值
if k == i:
continue # don't calc for i, waste of time
# 求 Ek误差：预测值-真实值的差
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)
# 求 Ek误差：预测值-真实值的差
Ej = calcEk(oS, j)
return j, Ej
def updateEk(oS, k): # after any alpha has changed update the new value in the cache
"""updateEk（计算误差值并存入缓存中。）
在对alpha值进行优化之后会用到这个值。
Args:
oS optStruct对象
k 某一列的行号
"""
# 求 误差：预测值-真实值的差
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]
def innerL(i, oS):
"""innerL
内循环代码
Args:
i 具体的某一行
oS optStruct对象
Returns:
0 找不到最优的值
1 找到了最优的值，并且oS.Cache到缓存中
"""
# 求 Ek误差：预测值-真实值的差
Ei = calcEk(oS, i)
# 约束条件 (KKT条件是解决最优化问题的时用到的一种方法。我们这里提到的最优化问题通常是指对于给定的某一函数，求其在指定作用域上的全局最小值。)
# 0<=alphas[i]<=C，但由于0和C是边界值，我们无法进行优化，因为需要增加一个alphas和降低一个alphas。
# 表示发生错误的概率：labelMat[i]*Ei 如果超出了 toler， 才需要优化。至于正负号，我们考虑绝对值就对了。
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.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
# 计算出一个新的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.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(dataMatIn, classLabels, C, toler, maxIter):
"""
完整SMO算法外循环，与smoSimple有些类似，但这里的循环退出条件更多一些
Args:
dataMatIn 数据集
classLabels 类别标签
C 松弛变量(常量值)，允许有些数据点可以处于分隔面的错误一侧。
控制最大化间隔和保证大部分的函数间隔小于1.0这两个目标的权重。
可以通过调节该参数达到不同的结果。
toler 容错率
maxIter 退出前最大的循环次数
Returns:
b 模型的常量值
alphas 拉格朗日乘子
"""
# 创建一个 optStruct 对象
oS = optStruct(mat(dataMatIn), mat(classLabels).transpose(), C, toler)
iter = 0
entireSet = True
alphaPairsChanged = 0
# 循环遍历：循环maxIter次 并且 （alphaPairsChanged存在可以改变 or 所有行遍历一遍）
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 = 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 # toggle entire set loop
elif (alphaPairsChanged == 0):
entireSet = True
print("iteration number: %d" % iter)
return oS.b, oS.alphas
def calcWs(alphas, dataArr, classLabels):
"""
基于alpha计算w值
@@ -81,126 +309,6 @@ def calcWs(alphas, dataArr, classLabels):
return w
'''
#######********************************
Non-Kernel VErsions below
#######********************************
'''
class optStruct:
def __init__(self, dataMatIn, classLabels, C, toler): # Initialize the structure with the parameters
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))) # first column is valid flag
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
def selectJ(i, oS, Ei): # this is the second choice -heurstic, and calcs Ej
maxK = -1
maxDeltaE = 0
Ej = 0
oS.eCache[i] = [1, Ei] # set valid #choose the alpha that gives the maximum delta E
validEcacheList = nonzero(oS.eCache[:, 0].A)[0]
if (len(validEcacheList)) > 1:
for k in validEcacheList: # loop through valid Ecache values and find the one that maximizes delta E
if k == i: continue # don't calc for i, waste of time
Ek = calcEk(oS, k)
deltaE = abs(Ei - Ek)
if (deltaE > maxDeltaE):
maxK = k
maxDeltaE = deltaE
Ej = Ek
return maxK, Ej
else: # in this case (first time around) we don't have any valid eCache values
j = selectJrand(i, oS.m)
Ej = calcEk(oS, j)
return j, Ej
def updateEk(oS, k): # after any alpha has changed update the new value in the cache
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]
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) # this has been changed from selectJrand
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) # added this for the Ecache
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]) # update i by the same amount as j
updateEk(oS, i) # added this for the Ecache #the update is in the oppostie direction
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(dataMatIn, classLabels, C, toler, maxIter): # full Platt SMO
oS = optStruct(mat(dataMatIn), mat(classLabels).transpose(), C, toler)
iter = 0
entireSet = True
alphaPairsChanged = 0
while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):
alphaPairsChanged = 0
if entireSet: # go over all
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: # go over non-bound (railed) alphas
nonBoundIs = 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 # toggle entire set loop
elif (alphaPairsChanged == 0):
entireSet = True
print("iteration number: %d" % iter)
return oS.b, oS.alphas
def plotfig_SVM(xArr, yArr, ws, b, alphas):
"""
参考地址：

7
src/python/6.SVM/svm-simple.py
View File

@@ -127,8 +127,7 @@ def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
continue
# eta是alphas[j]的最优修改量，如果eta==0，需要退出for循环的当前迭代过程
# 如果ETA为0，那么计算新的alphas[j]就比较麻烦了, 为什么呢？ 因为2个值一样。
# 2ab <= a^2 + b^2
# 参考《统计学习方法》李航-P125~P128<序列最小最优化算法>
eta = 2.0 * dataMatrix[i, :]*dataMatrix[j, :].T - dataMatrix[i, :]*dataMatrix[i,:].T - dataMatrix[j, :]*dataMatrix[j, :].T
if eta >= 0:
print("eta>=0")
@@ -145,8 +144,9 @@ def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
# 然后alphas[i]和alphas[j]同样进行改变，虽然改变的大小一样，但是改变的方向正好相反
alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])
# 在对alpha[i], alpha[j] 进行优化之后，给这两个alpha值设置一个常数b。
# w= Σ[1~n] ai*yi*xi => b = yi- Σ[1~n] ai*yi(xi*xj)
# 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 = 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]):
@@ -158,6 +158,7 @@ def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
alphaPairsChanged += 1
print("iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
# 在for循环外，检查alpha值是否做了更新，如果在更新则将iter设为0后继续运行程序
# 知道更新完毕后，iter次循环无变化，才推出循环。
if (alphaPairsChanged == 0):
iter += 1
else: