完整版和简化版的smo不通之处在于alpha的选择方式上，alpha的更改和代数运算和简化版本，完整版的alpha的选择方式采用了启发式进行选择，前面文章已经详解了，这里再简单的叙述一下：
Platt SMO算法是通过一个外循环来选择第一个alpha值的，并且其选择过程会在两种方式之间进行交替即：一种方式是在所有数据集上进行单遍扫描，另一种方式则是在非边界alpha中实现单遍扫描，而所谓非边界alpha指的就是那些不等于边界0或者C的alpha值。对整个数据集扫描相当容易，而实现非边界alpha值的扫描时，首先需要建立这些alpha值的列表，然后在对这个表进行遍历，同时该步骤跳过那些已知不会改变的alpha值。
在选择第一个alpha值后，算法会通过一个内循环来选择第二个alpha值，在优化过程中，会通过最大步长的方式来获的第二个alpha值，在简易的smo中，我会在选择j后计算误差率E，但这里会建立一个全局的缓存用于保存误差值，并从中选择使得步长或者说Ei-Ej最大的alpha值，为什么是这样呢？看下式：
从上式可以清楚的看到的步长取决于和，而又无法人为改变，所以能改变的只有，使其差值越大，说明前进步长越大，下面给出代码，代码均已详细注释：
# 完整版smo算法实现# 下面的函数就是参数方面的存储定义# 简化版的smo都是在函数体内定义，这里单独拿出来，不同的是处理E值，增加了ecache项class optStruct:def __init__(self, dataMatIn, classLabels, C, toler):# dataMatIn、classLabels为输入数据和标签，C为松弛因子，toler为容忍度，简单来说就是alpha多大才认为不在变化了# kTup为核函数的信息self.X = dataMatInself.labelMat = classLabelsself.C = Cself.tol = tolerself.m = np.shape(dataMatIn)[0]self.alphas = np.mat(np.zeros((self.m, 1))) # 给alpha占位self.b = 0self.eCache = np.mat(np.zeros((self.m, 2)))# 首先这是一个矩阵，两列，第一列是表面E是否有效的标志位，第二列给出的是实际的E值# 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)# 其实就是f(x) = ∑alpha_i*y_i*<x_i,x> + b,前面理论也详细讲了这一点，就是预测x的类别，只是单独拿出来了# 然后计算 差值 Ekdef calcEk(oS, k):fXk = float(np.multiply(oS.alphas,oS.labelMat).T*(oS.X*oS.X[k,:].T)) + oS.bEk = fXk - float(oS.labelMat[k])return Ek# 用来选择第二个alpha参数，或者说内循环的alpha值，其目标是选择合适的alpha值使其步长最大def selectJ(i, oS, Ei): # this is the second choice -heurstic, and calcs EjmaxK = -1 # 定义使的计算E的差值最大maxDeltaE = 0 # E值的最大变化量Ej = 0 # 使的差值最大的那个EjoS.eCache[i] = [1, Ei] # 这是把计算出来的差值记录并标为有效validEcacheList = np.nonzero(oS.eCache[:, 0].A)[0] # 刚开始有点不理解.A是什么意思，查了一下，解释一下# 情况下面的示例,其中.A的意思是把矩阵转换成数组，因为oS.eCache[:, 0].A，取的是一列的数据E值有效位，转换后还是一列# 而nonzero(oS.eCache[:, 0].A)[0]返回的是非零的目录列表值并取出赋给validEcacheList，共下面处理if (len(validEcacheList)) > 1:for k in validEcacheList: # 索引validEcacheList的标号，求最大差值if k == i: continue # 如果相等，说明选取两个alpha一样，直接跳出本次循环，进行下一次循环Ek = calcEk(oS, k) # 计算EkdeltaE = abs(Ei - Ek) # 求两个Ek的差值，选择最大的步长if (deltaE > maxDeltaE): # 比较，选择最大的差值maxK = k;maxDeltaE = deltaE;Ej = Ekreturn maxK, Ejelse:# 这个是为第一次准备的，不满足上面的条件即(len(validEcacheList)) > 1，因为第一次的E都为0，# 所以长度都为0，不满足条件，来到这里，进行随机抽取alpha 的 ijj = selectJrand(i, oS.m)Ej = calcEk(oS, j)return j, Ej'''In [12]: a = [[0,0,3],[0,5,0],[7,0,9]]In [13]: b = np.mat(a)In [14]: c = b.AIn [15]: print('type(a): ',type(a),'type(b): ',type(b),'type(c): ',type(c))type(a): <class 'list'> type(b): <class 'numpy.matrixlib.defmatrix.matrix'> type(c): <class 'numpy.ndarray'>In [16]: bOut[16]:matrix([[0, 0, 3],[0, 5, 0],[7, 0, 9]])In [17]: cOut[17]:array([[0, 0, 3],[0, 5, 0],[7, 0, 9]])In [18]: v =np.nonzero(c)[0]In [19]: vOut[19]: array([0, 1, 2, 2], dtype=int64)'''# 每一次循环计算的E都会缓存，即用于更新def updateEk(oS, k): # after any alpha has changed update the new value in the cacheEk = calcEk(oS, k)oS.eCache[k] = [1, Ek]# 完整版的SMO优化例程，和前面的简化版类似def innerL(i, oS):Ei = calcEk(oS, i) # 计算第一个alpha的差值# 寻找非边界的alpha，为寻找第二个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)):j,Ej = selectJ(i, oS, Ei) # 挑选第二个alpha，并得到Ej# 保存旧的alpha，为后面比较使用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, 大家还记得上篇的参数更新里的: K_11 + K_22 - 2K_12吗，其实就是这个eta，有点不同，# 是因为假设条件不一样，但是原理相同eta = 2.0 * oS.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].Tif eta >= 0:print("eta>=0")return 0# 更新alphaoS.alphas[j] -= oS.labelMat[j]*(Ei - Ej)/etaoS.alphas[j] = clipAlpha(oS.alphas[j],H,L)# 把选出的第二个alpha对应的E添加到EcacheupdateEk(oS, j) #added this for the Ecache# 比较第二alpha的值变化是否在容忍度以内，如果在这属于不在变化的alpha了，则返回，反之继续执行if (abs(oS.alphas[j] - alphaJold) < 0.00001):print("j not moving enough")return 0# 更新另一个alpha，即第一个alphaoS.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# 计算b值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,:].Tb2 = 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# 根据条件更新bif (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2else: oS.b = (b1 + b2)/2.0return 1else:return 0# 外循环即第一个alpha的选择def smoPK(dataMatIn, classLabels, C, toler, maxIter): #full Platt SMO# 把数据传入数据结构中oS = optStruct(np.mat(dataMatIn), np.mat(classLabels).transpose(), C, toler)iter = 0entireSet = TruealphaPairsChanged = 0# 执行while的条件是1.迭代次数达到最大，2.alpha值存在变化,或者整个开始启动即entireSet = Truewhile (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):alphaPairsChanged = 0if entireSet: #go over allfor i in range(oS.m): # 遍历所有的alpha值alphaPairsChanged += innerL(i,oS) # 从上面我们知道进入innerL是第一个alpha选取就为i# 第二个跟剧最大的E差值进行选择也是在innerl完成选择# 如果寻找满足条件则返回1，反之返回0，所以alphaPairsChanged是记录变化的alpha的变化次数print("fullSet, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged))iter += 1 # 整个遍历完以后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 += 1if entireSet:entireSet = False # 开始时先遍历整个alpha，然后迭代完成一次以后，开始转向非边界进行遍历elif (alphaPairsChanged == 0):entireSet = True # 如果遍历完所有的非边界都没有可更新的alpha，则继续转为全局遍历print("iteration number: %d" % iter)return oS.b,oS.alphas# 计算w，大家还记的求w的式子吗？# W = ∑alpha_i*y_i*x_idef 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)# 虽然该过程遍历所有的样本点，但是通过前面的原理详解可知，大部分的alpha值为0，只有支持向量的值不为零# 所以这个计算还是很简单的return w
下面给出测试代码：
b, alphas = SMO_simplify.smoPK(dataarr, labelarr, 0.6, 0.001, 40)print('b: ',b)print('alphas: ',alphas[alphas>0])ws = SMO_simplify.calcWs(alphas, dataarr, labelarr)print('ws: ',ws)# 对数据进行预测datmat = np.mat(dataarr)print('datmat[0]*np.mat(ws) + b = ', datmat[0]*np.mat(ws) + b)print('labelarr[0]',labelarr[0])err=0# 计算全体数据的错误率for i in range(len(dataarr)):y = datmat[i]*np.mat(ws) + bif y > 0:y = 1else:y = -1if y != labelarr[i]:err += 1;print('err: ', err)print('错误率为: ', err/len(dataarr))
测试 结果：
non-bound, iter: 1 i:46, pairs changed 0non-bound, iter: 1 i:55, pairs changed 0non-bound, iter: 1 i:94, pairs changed 0iteration number: 2j not moving enoughfullSet, iter: 2 i:0, pairs changed 0fullSet, iter: 2 i:1, pairs changed 0fullSet, iter: 2 i:2, pairs changed 0j not moving enoughb: [[-2.89901748]]alphas: [[0.06961952 0.0169055 0.0169055 0.0272699 0.04522972 0.02726990.0243898 0.06140181 0.06140181]]ws: [[ 0.65307162][-0.17196128]]datmat[0]*np.mat(ws) + b = [[-0.92555695]]labelarr[0] -1.0err: 0错误率为: 0.0
简单测试100%正确，这是线性分类，如果是非线性呢？下一节将详细讨论核函数实战代码
