支持向量机整理


参考:

Alexandre KOWALCZYK大神的SVM Tutorial

http://blog.csdn.net/alvine008/article/details/9097111

http://blog.csdn.net/zouxy09/article/details/17292011

http://blog.csdn.net/zy_zhengyang/article/details/45009431


介绍整理了SVM的基本数学推导,SMO算法的基本过程,LibSVM的用法,SMO算法的python实现四大部分(pdf版及测试数据相关电子书百度云 密码:itih)

基本型

支持向量机(Support vector machine)通常用来解决二分类问题,首先考虑一种最简单的情况,即数据线性可分的情况。

  • 划分超平面

考虑下面两点的划分,对于平面上的点集,存在一条直线将其完全划分,但是在实际情况中,样本属性通常是多维的,所以称划分线为划分超平面。

SVM支持向量机推导,工具介绍及python实现-LMLPHP通常将划分超平面表示为:

\[\vec w·\vec x+b=0
\]

注:有的地方将\(\vec w\)表示为\(\vec w^T\) 是考虑到w和x相乘需要列数和行数相等,但是参照线性代数里可省略转置

  • Margin及支持向量

显然存在最优的划分使两类数据集距离最远,此时的泛化性能最好,称最大距离为margin,最优的划分为最优超平面,两边数据集的边界为支持向量。

SVM支持向量机推导,工具介绍及python实现-LMLPHP

首先定义支持向量为:

\[\vec w\vec x+b=s \quad or \quad \vec w\vec x+b=-s
\]

同除以s可得:

\[\vec w\vec x+b=1 \quad or \quad \vec w\vec x+b=-1
\]

则对于右上方的点集y定义为1,左下点集y定义为-1

显然右上方的点集为:\(\vec w\vec x+b \ge1 \quad y=1\),左下方的点集可表示为:\(\vec w\vec x+b1\le-1 \quad y=-1\)

可合并为

\[y(\vec w\vec x+b) \ge1
\]

  • Margin的表示

计算margin如下图,长度用k表示,\(H_1\),\(H_0\)为支持向量

SVM支持向量机推导,工具介绍及python实现-LMLPHP

对于\(H_1\)上的点\(z_0\)满足

\[\vec w\vec z_0+b=1
\]

显然有

\[\vec z_0=\vec x_0 +\vec k
\]

带入化简可得

\[\vec w\vec x_0+b+\vec w \vec k=1
\]

\[\vec w \left (k\frac {\vec w}{||\vec w||} \right )=2
\]

\[k=\frac {2}{||\vec w||}
\]

最大的k值即为所求

  • 有约束的最优化问题

通过上面的计算,可将问题表示为约束条件下的最优化问题

\[max\quad \frac {2}{||\vec w||}
\\
s.t.\quad y_i(\vec w\vec x_i+b) \ge1 \quad i=1...m
\]

目标函数经常表示为 \(min\ \frac{1}{2}||\vec w||^2\)

对偶问题

通常引入拉格朗日函数简化上述问题,称此类问题为对偶问题(The wolfe dual problem)

  • 拉格朗日乘数法

拉格朗日乘数法用来解决等式约束下的最优化问题

将目标函数表示为\(f(\vec w)= \frac {2}{||\vec w||}\) ,约束函数表示为\(g(\vec w,b)=y_i(\vec w\vec x_i+b) \ge1\) ,定义如下拉格朗日函数\(L(\vec w,b,\alpha)\),其中\(\alpha_i\)为拉格朗日算子

\[L(\vec w,b,\alpha)=f(\vec w)-\sum_{i=1}^m \alpha_ig_i(\vec w,b)
\\
L(\vec w,b,\alpha)=\frac{1}{2}||\vec w||^2-\sum_{i=1}^m \alpha_ig_i[y_i(\vec w\vec x_i+b)-1]
\]

求解要求\(\nabla L=0\),可以解得

\[\vec w=\sum_{i=1}^m \alpha_i y_i \vec x_i
\\
\sum_{i=1}^m \alpha_i y_i=0
\]

将\(\vec w\)带入\(L(\vec w,b,\alpha)\)可得

\[f(\vec w,b)=\sum_{i=1}^m\alpha_i-\frac{1}{2} \sum_{i=1}^m\sum_{j=1}^m \alpha_i\alpha_jy_iy_j \vec x_i \vec x_j
\]

所以问题转化为

\[max \quad \sum_{i=1}^m\alpha_i-\frac{1}{2} \sum_{i=1}^m\sum_{j=1}^m \alpha_i\alpha_jy_iy_j \vec x_i \vec x_j
\\
s.t. \quad \alpha_i\ge0 \quad i=1...m
\\
\quad \sum_{i=1}^m \alpha_i y_i=0
\]

  • KKT条件

注意到上面的拉格朗日乘数法用来解决等式优化问题,在对于不等式约束的优化问题时还需附加如下的KKT条件

\[\begin{cases}
\ \alpha_i \ge 0
\\ \ y_if(x_i)-1 \ge 0
\\ \ \alpha_i(y_if(\vec x_i)-1) \ge 0
\end{cases}
\]

软间隔

  • 线性不可分时

实际遇到的问题数据总是有噪音的,也就是说可能存在一两个可忽略的特殊点使margin大大减小,也可能有一两个特殊点导致集合线性不可分(如下图),所以我们需要允许模型犯一定的错误来过滤掉噪音,即允许间隔的调整称为软间隔

SVM支持向量机推导,工具介绍及python实现-LMLPHP

  • 松弛变量

引入松弛变量(slack variable)\(\xi_i\) ,则约束函数可写为

\[y_i(\vec w\vec x_i+b) \ge1-\xi_i
\]

  • 控制出错

同时引入变量C表示对\(\xi_i\)的控制程度,即表示允许模型忽略噪音的程度,显然当C->\(\infty\)退化成基本型的硬间隔,C->0时成为无约束问题

引入C后基本型可表示为

\[max\quad \frac {1}{2}{||\vec w||}^2+C\sum_{i=1}^m \xi_i
\\
s.t.\quad y_i(\vec w\vec x_i+b) \ge1-\xi_i \quad
\\
\xi_i \ge 0 \quad i=1...m
\]

同样由拉格朗日乘数法可得

\[max \quad \sum_{i=1}^m\alpha_i-\frac{1}{2} \sum_{i=1}^m\sum_{j=1}^m \alpha_i\alpha_jy_iy_j\vec x_i \vec x_j
\\
s.t. \quad 0\le \alpha_i\le C \quad i=1...m
\\
\quad \sum_{i=1}^m \alpha_i y_i=0
\]

核方法

  • 升维

考虑更通常的情况,当点集线性不可分时,更有趣的方法是升维,如下图在二维平面内不能找到一个划分超平面

SVM支持向量机推导,工具介绍及python实现-LMLPHP

然而当在更高维的情况比如三维,我们就可以找到划分超平面来解决问题

SVM支持向量机推导,工具介绍及python实现-LMLPHP

这种方法通常成为核方法(kernel)

  • 构造核函数

比如我们可以简单的\((\vec x,\vec y)\) ->\((\vec x,\vec x\vec y,\vec y)\) 将点集升维,但是当样本空间很大时代价太高,考虑上面的目标函数

SVM支持向量机推导,工具介绍及python实现-LMLPHP

我们可以利用目标函数中的点积来构建核函数避免计算量的增加

常用核函数

  1. 线性核:\(k(\vec x_i,\vec y_i)=x_i·y_i\)
  2. 多项式核:\(k(\vec x,\vec y)=(\vec x_i·\vec y_i)^d\) ,\(d\ge 1\)为多项式的次数
  3. 高斯核(RBF,常用):\(k(\vec x,\vec y)=exp(-\frac {||x_i-x_j||^2}{2 \sigma ^2})\),\(\sigma\)为高斯核的带宽

SMO算法

  • SMO概念

对于上述问题的求解,SMO算法每次迭代优化两个\(\alpha\) (若只选择1个,由约束条件可以直接解出\(\alpha\),就不是最优化问题了),把原始求解N个参数二次规划问题分解成很多个子二次规划问题分别求解,每个子问题只需要求解2个参数,方法类似于坐标上升,节省时间成本和降低了内存需求。每次启发式选择两个变量进行优化,不断循环,直到达到函数最优值。

  1. 视为一个二元函数

SMO算法选择同时优化两个参数,固定其他N-2个参数,假设选择的变量为α1,α2,固定其他参数α3,α4,...,αN,由于参数α3,α4,...,αN的固定,可以简化目标函数为只关于α1,α2的二元函数,Const表示常数项(不包含变量α1,α2的项)。

\[minW(\alpha_i,\alpha_j)=\frac {1}{2}K_{11}\alpha_1^2+\frac {1}{2}K_{22}\alpha_2^2+y_1y_2K_{12}\alpha_1 \alpha_2-(\alpha_1+\alpha_2)+y_1v_1\alpha_1+y_2v_2\alpha_2+Const
\\
其中,v_i=\sum _{j=3}^m \alpha_jy_jK(x_i,x_j) \quad i=1,2
\]

  1. 视为一元函数求极值点

由等式约束:\(\alpha_i y1+\alpha_2 y2=-\sum _{i=3}^{m} \alpha_i y_i=\xi\),可见\(\xi\)为定值,将\(\alpha_1 y1+\alpha_2 y2=\xi\)两边同时乘以\(y_1\),且\(y_1^2=1\),得

\[\alpha_1=(\xi-y_2\alpha_2)y_1
\]

代回上述二元函数可得关于\(\alpha_2\)的一元函数,求导令其等于0可得

\[\frac {\partial W(\alpha_2)}{\partial \alpha_2}=(K_{11}+K_{22}-2K_{12})\alpha_2-K_{11}\xi y_2+K_{12}\xi y_2+y_1y_2-1-v_1y_1+v_2y_2=0
\]

将上述式子带入\(\frac {\partial W(\alpha_2)}{\partial \alpha_2}=0\)可得

\[\alpha_2^{new}=\alpha_2^{old}+\frac{y_2(E_1-E_2)}{\eta},其中\eta=K_{11}+K{22}-2K_{12}
\]

  1. 对上述解进行修剪

考虑到约束条件

\[0\le \alpha_i \le C
\\
\alpha_1y_1+\alpha_2y_2=\xi
\]

SVM支持向量机推导,工具介绍及python实现-LMLPHP

定义\(\alpha_2\)区间为[L,H],则有

当\(y_1\not=y_2\)时,\(L=max(0,\alpha_2^{old}-\alpha_1^{old})\);\(H=min(C,C+\alpha_2^{old}-\alpha_1^{old})\);

当\(y_1=y_2\)时,\(L=max(0,\alpha_2^{old}+\alpha_1^{old}-C)\);\(H=min(C,\alpha_2^{old}+\alpha_1^{old})\);

  1. 求解\(\alpha_1^{new}\)

由\(\alpha_1^{old} y1+\alpha_2^{old}y2=\alpha_1^{new} y1+\alpha_2^{new}y2\) ,可得

\[\alpha_1^{new}=\alpha_1^{old} +y_1y_2(\alpha_2^{old}-\alpha_2^{new})
\]

  1. 考虑临界情况

大部分情况下由\(\eta=K_{11}+K_{22}-2K_{12}>0\),但当如下两种情况时\(\alpha_2\)取临界值L或H

  • η<0,当核函数K不满足Mercer定理时,矩阵K非正定;
  • η=0,样本x1与x2输入特征相同;

也可以理解成对一元函数求二阶导就是\(\eta=K_{11}+K_{22}-2K_{12}>0\)

  • 当\(\eta\)<0时,目标函数为凸函数,没有极小值,极值在定义域边界处取得
  • 当\(\eta\)=0时,目标函数为单调函数,同样在边界处取极值
  1. 启发式选择变量

上述分析是在从N个变量中已经选出两个变量进行优化的方法,下面分析如何高效地选择两个变量进行优化,使得目标函数下降的最快。

  • \(\alpha_1\)

第一个变量的选择称为外循环

先遍历非边界样本集(0<\(\alpha_1\)<C)中违反KKT的\(\alpha_1\)作为第一个变量,同样依据相关规则选择第二个变量,对此两个变量进行优化。

当遍历完非边界样本集后,再次回到边界样本集中寻找,即在边界样本集与非边界样本集上来回切换,寻找违反KKT条件的\(\alpha_1\)作为第一个变量。直到遍历整个样本集后,没有违反KKT条件αi,然后退出。

边界上的样本对应的\(\alpha_1\)=0或者\(\alpha_1\)=C,在优化过程中很难变化,然而非边界样本0<\(\alpha_1\)<C会随着对其他变量的优化会有大的变化。

  • \(\alpha_2\)

第二个变量的选择过程为内循环

第二个变量的选择希望能使\(\alpha_2\)有较大的变化,通常为每个样本的\(E_i\)保存在一个列表中,选择最大的\(|E_1-E_2|\)来近似最大化步长,由于α2是依赖于\(|E_1-E_2|\),当\(E_1\)为正时,那么选择最小的\(E_i\)作为\(E_2\),如果\(E_1\)为负,选择最大\(E_i\)作为\(E_2\)

  1. 更新b

每完成对两个变量的优化后,要对b的值进行更新,因为b的值关系到f(x)的计算,即关系到下次优化时\(E_I\)的计算。

SVM支持向量机推导,工具介绍及python实现-LMLPHP

LibSVM工具包求解

LibSVM是台湾林智仁(Chih-Jen Lin)教授2001年开发的一套支持向量机的库,这套库运算速度快,可以很方便的对数据做分类或回归。由于libSVM程序小,运用灵活,输入参数少,并且是开源的,易于扩展,因此成为目前国内应用最多的SVM的库。(以下参考博客)

  • 基本介绍

可以从https://www.csie.ntu.edu.tw/~cjlin/libsvm/下载,解压到任意目录,有如下几个目录

Java——主要是应用于java平台;

Python——是用来参数优选的工具,稍后介绍;

svm-toy——一个可视化的工具,用来展示训练数据和分类界面,里面是源码,其编译后的程序在windows文件夹下;

tools——主要包含四个python文件,用来数据集抽样(subset),参数优选(grid),集成测试(easy),数据检查(checkdata);

windows——包含libSVM四个exe程序包,一般直接用这四个就够了

其他.h和.cpp文件都是程序的源码,可以编译出相应的.exe文件。

  • 自带数据测试
  1. 命令行下进入libsvm下的windows目录
  2. 输入命令 .\svm-train ..\heart_scale test_model

SVM支持向量机推导,工具介绍及python实现-LMLPHP

heart_scale是目录下的已经存在的样本文件,以后改成自己的文件名就可以了, test_model——是创建模型

其中,#iter为迭代次数,

nu是你选择的核函数类型的参数

obj为SVM文件转换为的二次规划求解得到的最小值

rho为判决函数的偏置项b

nSV为标准支持向量个数(0<a[i]<c)

nBSV为边界上的支持向量个数(a[i]=c),

Total nSV为支持向量总个数(对于二分类来说,因为只有一个分类模型Total nSV = nSV,但是对于多类,这个是各个分类模型的nSV之和)。

在目录下,还可以看到产生了一个train.model文件,可以用记事本打开,记录了训练后的结果。

svm_type c_svc//所选择的svm类型,默认为c_svc

kernel_type rbf//训练采用的核函数类型,此处为RBF核

gamma 0.0769231//RBF核的参数γ

nr_class 2//类别数,此处为两分类问题

total_sv 132//支持向量总个数

rho 0.424462//判决函数的偏置项b

label 1 -1//原始文件中的类别标识

nr_sv 64 68//每个类的支持向量机的个数

SV//以下为各个类的权系数及相应的支持向量

  • libsvm使用规范
  1. libsvm的数据格式

Label 1:value 2:value(value为0时序号可直接省略)

  1. svm-train的用法

svmtrain我们在前面已经接触过,他主要实现对训练数据集的训练,并可以获得SVM模型。

用法: svmtrain [options] training_set_file [model_file]

其中,options为操作参数,可用的选项即表示的涵义如下所示:

-s设置svm类型:

0 – C-SVC

1 – v-SVC

2 – one-class-SVM

3 –ε-SVR

4 – n - SVR

-t设置核函数类型,默认值为2

0 --线性核:u'**v*

1 --多项式核:(gu'*v+coef0)degree*

2 -- RBF核:exp(-γ*||u-v||2)

3 -- sigmoid核:tanh(γu'*v+coef*0)

-d degree:设置多项式核中degree的值,默认为3

-gγ:设置核函数中γ的值,默认为1/k,k为特征(或者说是属性)数;

-r coef 0:设置核函数中的coef 0,默认值为0;

-c cost:设置C-SVC、ε-SVR、n - SVR中从惩罚系数C,默认值为1;

-n v:设置v-SVC、one-class-SVM与n - SVR中参数n,默认值0.5;

-pε:设置v-SVR的损失函数中的e,默认值为0.1;

-m cachesize:设置cache内存大小,以MB为单位,默认值为40;

-eε:设置终止准则中的可容忍偏差,默认值为0.001;

-h shrinking:是否使用启发式,可选值为0或1,默认值为1;

-b概率估计:是否计算SVC或SVR的概率估计,可选值0或1,默认0;

-wi weight:对各类样本的惩罚系数C加权,默认值为1;

-v n:n折交叉验证模式;

model_file:可选项,为要保存的结果文件,称为模型文件,以便在预测时使用。

默认情况下,只需要给函数提供一个样本文件名就可以了,但为了能保存结果,还是要提供一个结果文件名,比如:test.model,则命令为:svmtrain test.txt test.model

  • svm-predict的用法

svmpredict是根据训练获得的模型,对数据集合进行预测。

用法:svmpredict [options] test_file model_file output_file

其中,options为操作参数,可用的选项即表示的涵义如下所示:

-b probability_estimates——是否需要进行概率估计预测,可选值为0或者1,默认值为0。

model_file ——是由svmtrain产生的模型文件;

test_file——是要进行预测的数据文件,格式也要符合libsvm格式,即使不知道label的值,也要任意填一个,svmpredict会在output_file中给出正确的label结果,如果知道label的值,就会输出正确率;

output_file ——是svmpredict的输出文件,表示预测的结果值。

Python实现SMO算法

#smo
from numpy import *
import time # calulate kernel value
def calcKernelValue(matrix_x, sample_x, kernelOption):
kernelType = kernelOption[0]
numSamples = matrix_x.shape[0]
kernelValue = mat(zeros((numSamples, 1))) if kernelType == 'linear':
kernelValue = matrix_x * sample_x.T
elif kernelType == 'rbf':
sigma = kernelOption[1]
if sigma == 0:
sigma = 1.0
for i in xrange(numSamples):
diff = matrix_x[i, :] - sample_x
kernelValue[i] = exp(diff * diff.T / (-2.0 * sigma**2))
else:
raise NameError('Not support kernel type! You can use linear or rbf!')
return kernelValue # calculate kernel matrix given train set and kernel type
def calcKernelMatrix(train_x, kernelOption):
numSamples = train_x.shape[0]
kernelMatrix = mat(zeros((numSamples, numSamples)))
for i in xrange(numSamples):
kernelMatrix[:, i] = calcKernelValue(train_x, train_x[i, :], kernelOption)
return kernelMatrix # define a struct just for storing variables and data
class SVMStruct:
def __init__(self, dataSet, labels, C, toler, kernelOption):
self.train_x = dataSet # each row stands for a sample
self.train_y = labels # corresponding label
self.C = C # slack variable
self.toler = toler # termination condition for iteration
self.numSamples = dataSet.shape[0] # number of samples
self.alphas = mat(zeros((self.numSamples, 1))) # Lagrange factors for all samples
self.b = 0
self.errorCache = mat(zeros((self.numSamples, 2)))
self.kernelOpt = kernelOption
self.kernelMat = calcKernelMatrix(self.train_x, self.kernelOpt) # calculate the error for alpha k
def calcError(svm, alpha_k):
output_k = float(multiply(svm.alphas, svm.train_y).T * svm.kernelMat[:, alpha_k] + svm.b)
error_k = output_k - float(svm.train_y[alpha_k])
return error_k # update the error cache for alpha k after optimize alpha k
def updateError(svm, alpha_k):
error = calcError(svm, alpha_k)
svm.errorCache[alpha_k] = [1, error] # select alpha j which has the biggest step
def selectAlpha_j(svm, alpha_i, error_i):
svm.errorCache[alpha_i] = [1, error_i] # mark as valid(has been optimized)
candidateAlphaList = nonzero(svm.errorCache[:, 0].A)[0] # mat.A return array
maxStep = 0; alpha_j = 0; error_j = 0 # find the alpha with max iterative step
if len(candidateAlphaList) > 1:
for alpha_k in candidateAlphaList:
if alpha_k == alpha_i:
continue
error_k = calcError(svm, alpha_k)
if abs(error_k - error_i) > maxStep:
maxStep = abs(error_k - error_i)
alpha_j = alpha_k
error_j = error_k
# if came in this loop first time, we select alpha j randomly
else:
alpha_j = alpha_i
while alpha_j == alpha_i:
alpha_j = int(random.uniform(0, svm.numSamples))
error_j = calcError(svm, alpha_j) return alpha_j, error_j # the inner loop for optimizing alpha i and alpha j
def innerLoop(svm, alpha_i):
error_i = calcError(svm, alpha_i) ### check and pick up the alpha who violates the KKT condition
## satisfy KKT condition
# 1) yi*f(i) >= 1 and alpha == 0 (outside the boundary)
# 2) yi*f(i) == 1 and 0<alpha< C (on the boundary)
# 3) yi*f(i) <= 1 and alpha == C (between the boundary)
## violate KKT condition
# because y[i]*E_i = y[i]*f(i) - y[i]^2 = y[i]*f(i) - 1, so
# 1) if y[i]*E_i < 0, so yi*f(i) < 1, if alpha < C, violate!(alpha = C will be correct)
# 2) if y[i]*E_i > 0, so yi*f(i) > 1, if alpha > 0, violate!(alpha = 0 will be correct)
# 3) if y[i]*E_i = 0, so yi*f(i) = 1, it is on the boundary, needless optimized
if (svm.train_y[alpha_i] * error_i < -svm.toler) and (svm.alphas[alpha_i] < svm.C) or\
(svm.train_y[alpha_i] * error_i > svm.toler) and (svm.alphas[alpha_i] > 0): # step 1: select alpha j
alpha_j, error_j = selectAlpha_j(svm, alpha_i, error_i)
alpha_i_old = svm.alphas[alpha_i].copy()
alpha_j_old = svm.alphas[alpha_j].copy() # step 2: calculate the boundary L and H for alpha j
if svm.train_y[alpha_i] != svm.train_y[alpha_j]:
L = max(0, svm.alphas[alpha_j] - svm.alphas[alpha_i])
H = min(svm.C, svm.C + svm.alphas[alpha_j] - svm.alphas[alpha_i])
else:
L = max(0, svm.alphas[alpha_j] + svm.alphas[alpha_i] - svm.C)
H = min(svm.C, svm.alphas[alpha_j] + svm.alphas[alpha_i])
if L == H:
return 0 # step 3: calculate eta (the similarity of sample i and j)
eta = 2.0 * svm.kernelMat[alpha_i, alpha_j] - svm.kernelMat[alpha_i, alpha_i] \
- svm.kernelMat[alpha_j, alpha_j]
if eta >= 0:
return 0 # step 4: update alpha j
svm.alphas[alpha_j] -= svm.train_y[alpha_j] * (error_i - error_j) / eta # step 5: clip alpha j
if svm.alphas[alpha_j] > H:
svm.alphas[alpha_j] = H
if svm.alphas[alpha_j] < L:
svm.alphas[alpha_j] = L # step 6: if alpha j not moving enough, just return
if abs(alpha_j_old - svm.alphas[alpha_j]) < 0.00001:
updateError(svm, alpha_j)
return 0 # step 7: update alpha i after optimizing aipha j
svm.alphas[alpha_i] += svm.train_y[alpha_i] * svm.train_y[alpha_j] \
* (alpha_j_old - svm.alphas[alpha_j]) # step 8: update threshold b
b1 = svm.b - error_i - svm.train_y[alpha_i] * (svm.alphas[alpha_i] - alpha_i_old) \
* svm.kernelMat[alpha_i, alpha_i] \
- svm.train_y[alpha_j] * (svm.alphas[alpha_j] - alpha_j_old) \
* svm.kernelMat[alpha_i, alpha_j]
b2 = svm.b - error_j - svm.train_y[alpha_i] * (svm.alphas[alpha_i] - alpha_i_old) \
* svm.kernelMat[alpha_i, alpha_j] \
- svm.train_y[alpha_j] * (svm.alphas[alpha_j] - alpha_j_old) \
* svm.kernelMat[alpha_j, alpha_j]
if (0 < svm.alphas[alpha_i]) and (svm.alphas[alpha_i] < svm.C):
svm.b = b1
elif (0 < svm.alphas[alpha_j]) and (svm.alphas[alpha_j] < svm.C):
svm.b = b2
else:
svm.b = (b1 + b2) / 2.0 # step 9: update error cache for alpha i, j after optimize alpha i, j and b
updateError(svm, alpha_j)
updateError(svm, alpha_i) return 1
else:
return 0 # the main training procedure
def trainSVM(train_x, train_y, C, toler, maxIter, kernelOption = ('rbf', 1.0)):
# calculate training time
startTime = time.time() # init data struct for svm
svm = SVMStruct(mat(train_x), mat(train_y), C, toler, kernelOption) # start training
entireSet = True
alphaPairsChanged = 0
iterCount = 0
# Iteration termination condition:
# Condition 1: reach max iteration
# Condition 2: no alpha changed after going through all samples,
# in other words, all alpha (samples) fit KKT condition
while (iterCount < maxIter) and ((alphaPairsChanged > 0) or entireSet):
alphaPairsChanged = 0 # update alphas over all training examples
if entireSet:
for i in xrange(svm.numSamples):
alphaPairsChanged += innerLoop(svm, i)
print '---iter:%d entire set, alpha pairs changed:%d' % (iterCount, alphaPairsChanged)
iterCount += 1
# update alphas over examples where alpha is not 0 & not C (not on boundary)
else:
nonBoundAlphasList = nonzero((svm.alphas.A > 0) * (svm.alphas.A < svm.C))[0]
for i in nonBoundAlphasList:
alphaPairsChanged += innerLoop(svm, i)
print '---iter:%d non boundary, alpha pairs changed:%d' % (iterCount, alphaPairsChanged)
iterCount += 1 # alternate loop over all examples and non-boundary examples
if entireSet:
entireSet = False
elif alphaPairsChanged == 0:
entireSet = True print 'Congratulations, training complete! Took %fs!' % (time.time() - startTime)
return svm # testing your trained svm model given test set
def testSVM(svm, test_x, test_y):
test_x = mat(test_x)
test_y = mat(test_y)
numTestSamples = test_x.shape[0]
supportVectorsIndex = nonzero(svm.alphas.A > 0)[0]
supportVectors = svm.train_x[supportVectorsIndex]
supportVectorLabels = svm.train_y[supportVectorsIndex]
supportVectorAlphas = svm.alphas[supportVectorsIndex]
matchCount = 0
for i in xrange(numTestSamples):
kernelValue = calcKernelValue(supportVectors, test_x[i, :], svm.kernelOpt)
predict = kernelValue.T * multiply(supportVectorLabels, supportVectorAlphas) + svm.b
if sign(predict) == sign(test_y[i]):
matchCount += 1
accuracy = float(matchCount) / numTestSamples
return accuracy

数据读写及输出的控制程序如下

from numpy import *
import SVM
import csv ################## test svm #####################
## step 1: load data
print "step 1: load data..."
dataSet = []
labels = [] in_file=csv.reader(open('spammer_order.csv','r')) maxm=[0]*35
minm=[0]*35
index=1
for in_line in in_file:
if index<3:
if index==1:
for i in range(4,32):
maxm[i]=float(in_line[i])
if index==2:
for i in range(4,32):
minm[i]=float(in_line[i])
else:
out_line=[]
if in_line[32]=='yes':
labels.append('1')
else:
labels.append('-1') for i in range(4,32):
if i!=6:
if i==4:
if in_line[4]=='female':
out_line.append('1')
else:
out_line.append('-1')
else:
max_value=maxm[i]
min_value=minm[i]
value=float(in_line[i])
n=2*(value-min_value)/(max_value-min_value)-1
out_line.append(str(n))
dataSet.append(map(float,out_line))
index=index+1
#3 860 1718
p1=3+(860-3)/3*2
p=860
p2=860+(1718-860)/3*2
labels=map(float,labels)
train_data=mat(dataSet[0:p1]+dataSet[p:p2])
train_labels=mat(labels[0:p1]+labels[p:p2]).T
test_data=mat(dataSet[p1:p]+dataSet[p2:])
test_labels=mat(labels[p1:p]+labels[p:p2]).T train_x = train_data
train_y = train_labels
test_x = test_data
test_y = test_labels ## step 2: training...
print "step 2: training..."
C = 1
toler = 0.001
maxIter = 428
svmClassifier = SVM.trainSVM(train_x, train_y, C, toler, maxIter, kernelOption = ('rbf', 0.440787)) ## step 3: testing
print "step 3: testing..."
accuracy = SVM.testSVM(svmClassifier, test_x, test_y) ## step 4: show the result
print "step 4: show the result..."
print 'The classify accuracy is: %.3f%%' % (accuracy * 100)
05-11 16:01