机器学习模型自我代码复现：SVM(SMO算法)

根据模型的数学原理进行简单的代码自我复现以及使用测试，仅作自我学习用。模型原理此处不作过多赘述，仅罗列自己将要使用到的部分公式。

如文中或代码有错误或是不足之处，还望能不吝指正。

SVM是一种分类算法，其思路是：找出一个（超）平面，以此分割2个类别。满足这一要求的（超）平面有许多。我们还期望，2个类别中离超平面距离最近的那几个点到超平面的距离尽可能的大，使得分割效果更“明显”，效果更好，离超平面距离最近的那几个点被称为“支持向量”。

将上述语言转化为一个公式，就是希望找到

$min_{\omega,b}\frac{1}{2} ||\omega||^2$

$s.t. y_i(w^Tx_i+b>=1,i=1,2,...,m)$

其中， $f(x_i) = w^Tx_i+b$ 为对应超平面； $min_{\omega,b}\frac{1}{2} ||\omega||^2$ 为间隔最大化的等价条件，此处省略推导过程。对其使用拉格朗日乘子法，得

$L(\omega,b,\alpha)=\frac{1}{2}||\omega||^2 + \Sigma_{i=1}^m\alpha_i(1-y1(\omega^Tx_i+b))$ 其中，α_i为拉格朗日乘子。

同时满足KKT条件：

$\begin{equation} \left\{ \begin{array}{lr} \alpha_i>=0, & \\ y_if(x_i)-1>=0\\ \alpha_i(y_if(x_i)-1)=0 & \end{array} \right. \end{equation}$

对于拉格朗日乘子法的解而言，如果有某个点的拉格朗日乘子 $\alpha_i$ 不为0，就代表其是“支持向量”，实际只有它们是对最终结果产生影响的。

根据SMO算法的思想，我们每次首先选择一个违反KKT条件最严重的 $\alpha_i$ ，之后寻找一个与样本 $x_i$ 距离最远的样本 $x_j$ 对应的 $\alpha_j$ 。因为这2个变量有很大的差别，所以带给函数值的变化更大。

以下是 $\alpha_i$ $\alpha_j$ 的求解过程，笔者正在自行理解中，不再赘述推理过程。

$\eta = K_{ii} + K_{jj} - 2K_{ij} \\ E_i = (\Sigma_k\alpha_ky_kK_{ik}+b)-y_i \\ alpha_j^{*} = alpha_2^{old} + \frac{y_2(E_i-E_j)}{\eta}$

$\alpha_j^{new}= \begin{equation} \left\{ \begin{array}{lr} H , \alpha_j^*>H, & \\ \alpha_j^*,L<=\alpha_j^*<=H \\ L,\alpha_j^*<L & \end{array} \right. \end{equation}$

其中，当y_i=y_j时， $L = max(0,\alpha_j^{old}-\alpha_i^{old}) H = min(C,C+\alpha_j^{old}-\alpha_i^{old})$

当y_i!=y_j时 $L = max(0,\alpha_j^{old}-\alpha_i^{old}-C) H = min(C,\alpha_j^{old}+\alpha_i^{old})$

$\alpha_i^{new} = \alpha_i^{old}+y_iy_j(\alpha_j^{old}-\alpha_j^{new})$

同时更新偏置项b

$b1 = -E_i - y_i*(\alpha_i^{new}-\alpha_i^{old})*K_{ii}-y_jK_{ij}(\alpha_j^{new}-\alpha_j^{old})+b$

$b2 = -E_j - y_i*(\alpha_i^{new}-\alpha_i^{old})*K_{ij}-y_jK_{jj}(\alpha_j^{new}-\alpha_j^{old})+b$

当 $0<\alpha_i^{new}<C$ 时，取b=b1

当 $0<\alpha_j^{new}<C$ 时，取b=b2

同时满足以上条件时,b=b1=b2

否则取b = (b1+b2)/2

最终，

$W = \sum \alpha_iy_ix_i$

此处的代码参考了GitHub上的GitHub - xinlianghu/svm: 用Python实现SVM多分类器

应该是以《统计学习方法》为基础的代码。

import numpy as np
import pandas as pd
import randomclass SVM:def __init__(self,data,labels,C,toler,gaussian_kernel=True) -> None:"""初始化data:自变量labels:因变量C:软间隔的惩罚因子toler:偏离KKT条件的容错率gaussian_kernel:是否使用高斯核Es:误差缓存m:数据量b:偏置"""self.kernel = gaussian_kernelself.data = dataself.labels = labelsself.C = Cself.toler = tolerself.gaussian_kernel = gaussian_kernelself.m = data.shape[0]self.alphas = np.zeros(self.m)self.Es = np.zeros((self.m,2)) #Ei是否有效 Ei当前值self.b = 0self.iter = 0def K(self,i,j,gamma = 1/2):"""核函数计算，gaussian_kernel==False时默认为线性核函数"""if self.gaussian_kernel:return np.exp(-gamma*np.linalg.norm(self.data[i,:]-self.data[j,:])**2)else:return sum(self.data[i,:]*self.data[j,:])def calc_Ek(self,k):"""根据拉格朗日因子alpha与核函数计算第k个误差项"""gxi = np.dot(self.alphas*self.labels,np.array([self.K(k,j) for j in range(self.m)]))+self.breturn gxi-float(self.labels[k])def train(self,max_Iter):"""训练此函数中首先选择alpha_i，之后再内循环中选择alpha_jmax_Iter:最大迭代次数"""last_changed = 0entrySet = True #是否遍历整个数据集last_changed = 0 #有几对alpha被改变了while self.iter<max_Iter and ((last_changed>0) or entrySet):last_changed = 0.0if entrySet:for i in range(self.m):last_changed+=self.innerL(i)self.iter += 1entrySet = Falseelse:nonBounds = np.nonzero((self.alphas>0)*(self.alphas<self.C))[0]for i in nonBounds:last_changed+=self.innerL(i)self.iter+=1if last_changed==0:entrySet = True#计算权重self.w = np.zeros((1,self.data.shape[1]))for i in range(self.m):if self.alphas[i]>0:self.w+=self.labels[i]*self.alphas[i]*self.data[i,:]def innerL(self,k):"""内循环根据给定条件选定alpha_j,与之前的alpha_i一起更新"""Ek = self.calc_Ek(k)last_changed = 0if (self.labels[k]*Ek < -self.toler and self.alphas[k]<self.C) or (self.labels[k]*Ek>self.toler and self.alphas[k]>0):self.Es[k] = [1,Ek]j,Ej = self.selectJ(k)res = self.update(k,j,Ek,Ej)last_changed+=resreturn last_changeddef update(self,i,j,Ei,Ej):"""更新alpha_i,alpha_j,b"""#首先更新alpha_jif(self.labels[i]==self.labels[j]):L = max(0,self.alphas[j]+self.alphas[i]-self.C)H = min(self.C,self.alphas[j]+self.alphas[i])else:L = max(0,self.alphas[j]-self.alphas[i])H = min(self.C,self.C+self.alphas[j]-self.alphas[i])if L==H:return 0eta = self.K(i,i)+self.K(j,j)-2*self.K(i,j)if eta == 0: #这里不考虑eta（分母）==0的情况return 0origin_alpha_j = self.alphas[j].copy()origin_alpha_i = self.alphas[i].copy()self.alphas[j] = self.alphas[j]+self.labels[j]*(self.labels[j]*(Ei-Ej)/eta)if self.alphas[j]>H:self.alphas[j]=Helif self.alphas[j]<L:self.alphas[j]=Lself.Es[j] = [1,self.calc_Ek(j)]if abs(self.alphas[j]-origin_alpha_j)<0.00001:return 0#更新alpha_iself.alphas[i] = self.alphas[i]+self.labels[i]*self.labels[j]*(origin_alpha_j-self.alphas[j])self.Es[i] = [1,self.calc_Ek(i)]#更新bbi = -self.Es[i][1]-self.labels[i]*(self.alphas[i]-origin_alpha_i)*self.K(i,i)-self.labels[j]*(self.alphas[j]-origin_alpha_j)*self.K(i,j)+self.bbj = -self.Es[j][1]-self.labels[i]*(self.alphas[i]-origin_alpha_i)*self.K(i,j)-self.labels[j]*(self.alphas[j]-origin_alpha_j)*self.K(j,j)+self.bprint(self.Es[i][1],self.Es[j][1],self.alphas[i]-origin_alpha_i,self.alphas[j]-origin_alpha_j)if self.C>self.alphas[i] and self.alphas[i]>0:self.b = bielif self.C>self.alphas[j] and self.alphas[j]>0:self.b = bjelse:self.b = (bi+bj)/2.0return 1def selectJ(self,i):"""根据i选择alpha_j"""Ei = self.Es[i][1]maxdelta = 0j=0Ej = 0validEsList = np.nonzero(self.Es[:,0])[0]#初次循环时，随机选择一个j，之后循环选择误差相距的alpha_jif len(validEsList)>1:for k in validEsList:if k == i:continueEk = self.calc_Ek(k)delta = abs(Ei-Ek)if delta>maxdelta:maxdelta = deltaj=kEj = Ekelse:j = iwhile j==i:j = int(np.random.uniform(0,self.m))Ej = self.calc_Ek(j)return j,Ej

使用自造数据集进行检验

x1 = np.random.randn(100)*2+10
x2 = np.random.randn(100)*2+20
x3 = np.random.randn(100)*2+10
x4 = np.random.randn(100)*2+20data = np.r_[np.c_[x1,x3],np.c_[x2,x4]]
labels = np.array([1]*100+[-1]*100)svm = SVM(data,labels,C=0.225,toler=0.001,gaussian_kernel=False)
svm.train(max_Iter=40)x_pre = np.linspace(5,25,200)
y_pre = -(svm.w[0][0] * x_pre + svm.b)/(svm.w[0][1])from matplotlib import pyplot as plt
import matplotlib as mplmpl.rcParams["font.family"]="SimHei"
mpl.rcParams["axes.unicode_minus"]=Falseplt.scatter(x1,x3)
plt.scatter(x2,x4)
plt.scatter(x_pre,y_pre)
plt.show()

从图像上看，的确找到了能够分割2个类别的线。

但是，依然有2个问题存在：一是b的值似乎对惩罚数C十分敏感，从理论上来说不应该这样的。

二是我不是很理解Es矩阵中的“有效”概念。我的理解是，不能随意取Es，否则各个 $\alpha$ 都有机会大于0。关于这两点我还得再研究一下。

机器学习模型自我代码复现：SVM(SMO算法)相关推荐

机器学习模型自我代码复现：GBDT
根据模型的数学原理进行简单的代码自我复现以及使用测试,仅作自我学习用.模型原理此处不作过多赘述,仅罗列自己将要使用到的部分公式. 如文中或代码有错误或是不足之处,还望能不吝指正. 集成学习,通过构建并 ...
机器学习模型自我代码复现：使用numpy复现CNN
根据模型的数学原理进行简单的代码自我复现以及使用测试,仅作自我学习用.模型原理此处不作过多赘述. 如文中或代码有错误或是不足之处,还望能不吝指正. 本文侧重于使用numpy重新写出一个CNN模型,故而 ...
支持向量机(SVM) SMO算法详解
1.寻找最大间隔训练样本集:D = { (x1, y1) , (x2, y2) , ... ,(xm, ym) } , yi ϵ { -1, +1} 划分超平面的线性方程:wTx + b = 0( ...
SVM SMO算法代码详细剖析
前言一:本文要结合SVM理论部分来看即笔者另一篇: SVM原理从头到尾详细推导二:有了理论部分下面就是直接代码啦,本文用四部分进行介绍:最简版的SMO,改进版platt SMO,核函数,sklea ...
[比赛记录] 主流机器学习模型模板代码+经验分享[xgb, lgb, Keras, LR]
向AI转型的程序员都关注了这个号??? 大数据挖掘DT数据分析公众号: datadw 最近打各种比赛,在这里分享一些General Model,稍微改改就能用的 XGBoost调参大全: http ...
主流机器学习模型模板代码+经验分享[xgb, lgb, Keras, LR]
刷比赛利器,感谢分享的人. 摘要最近打各种比赛,在这里分享一些General Model,稍微改改就能用的环境: python 3.5.2 XGBoost调参大全: http://blog.csd ...
论文解读+代码复现【AIDD】贝叶斯、决策树、随机森林+2种机器学习模型在癌症治疗药物发现中的应用
AIDD(AI Drug Discovery & Design):是近年来非常火热的技术应用,且已经介入到新药设计到研发的大部分环节当中,为新药发现与开发带来了极大的助力.倾向于机器对数据库信 ...
SVM问题的求解方法SMO算法
西瓜书学习笔记:支持向量机(6.1-6.2)笔记 1.SMO算法思路讲解列表就是西瓜书上面的公式6.11 2.SMO算法简单实现的点就是支持向量机上面的点. 30:43分钟开始讲code exam ...
手把手教你使用Flask轻松部署机器学习模型（附代码链接） | CSDN博文精选
作者 | Abhinav Sagar 翻译 | 申利彬校对 | 吴金笛来源 | 数据派THU(ID:DatapiTHU) 本文旨在让您把训练好的机器学习模型通过Flask API 投入到生产环境 ...

机器学习模型自我代码复现：SVM(SMO算法)

机器学习模型自我代码复现：SVM(SMO算法)相关推荐

最新文章

热门文章