SVM支持向量机算法及其代码实现.doc

支持向量机算法及其代码实现支持向量机（SVM），起初由vapnik提出时，是作为寻求最优（在一定程度上）二分类器的一种技术。後来它又被拓展到回归和聚类应用。SVM是一种基于核函数的方法，它通过某些核函数把特征向量映射到高维空间，然後建立一个线性判别函数（或者说是一个高维空间中的能够区分训练数据的最优超平面，参考异或那个经典例子）。假如SVM没有明确定义核函数，高维空间中任意两点距离就需要定义。 解是最优的在某种意义上是两类中距离分割面最近的特征向量和分割面的距离最大化。离分割面最近的特征向量被称为”支撑向量”，意即其它向量不影响分割面（决策函数）。 有很多关于SVM的参考文献，这是两篇较好的入门文献。 【Burges98】 C. Burges. A tutorial on support vector machines for pattern recognition, Knowledge Discovery and Data Mining 2(2), 1998. (available online at 1). LIBSVM - A Library for Support Vector Machines. By Chih-Chung Chang and Chih-Jen Lin (2) CvSVM支撑矢量机 class CvSVM : public CvStatModel /继承自基类CvStatModelpublic: / SVM type enum C_SVC=100, NU_SVC=101, ONE_CLASS=102, EPS_SVR=103, NU_SVR=104 ;/SVC是SVM分类器，SVR是SVM回归 / SVM kernel type enum LINEAR=0, POLY=1, RBF=2, SIGMOID=3 ; /提供四种核函数，分别是线性，多项式，径向基，sigmoid型函数。 CvSVM(); virtual CvSVM(); CvSVM( const CvMat* _train_data, const CvMat* _responses, const CvMat* _var_idx=0, const CvMat* _sample_idx=0, CvSVMParams _params=CvSVMParams() ); virtual bool train( const CvMat* _train_data, const CvMat* _responses, const CvMat* _var_idx=0, const CvMat* _sample_idx=0, CvSVMParams _params=CvSVMParams() ); virtual float predict( const CvMat* _sample ) const; virtual int get_support_vector_count() const; virtual const float* get_support_vector(int i) const; virtual void clear(); virtual void save( const char* filename, const char* name=0 ); virtual void load( const char* filename, const char* name=0 ); virtual void write( CvFileStorage* storage, const char* name ); virtual void read( CvFileStorage* storage, CvFileNode* node ); int get_var_count() const return var_idx ? var_idx-cols : var_all; protected: .;CvSVMParamsSVM训练参数struct struct CvSVMParams CvSVMParams(); CvSVMParams( int _svm_type, int _kernel_type, double _degree, double _gamma, double _coef0, double _C, double _nu, double _p, CvMat* _class_weights, CvTermCriteria _term_crit ); int svm_type; int kernel_type; double degree; / for poly double gamma; / for poly/rbf/sigmoid double coef0; / for poly/sigmoid double C; / for CV_SVM_C_SVC, CV_SVM_EPS_SVR and CV_SVM_NU_SVR double nu; / for CV_SVM_NU_SVC, CV_SVM_ONE_CLASS, and CV_SVM_NU_SVR double p; / for CV_SVM_EPS_SVR CvMat* class_weights; / for CV_SVM_C_SVC CvTermCriteria term_crit; / termination criteria;svm_type，SVM的类型： CvSVM:C_SVC - n(n=2)分类器，允许用异常值惩罚因子C进行不完全分类。 CvSVM:NU_SVC - n类似然不完全分类的分类器。参数nu取代了c，其值在区间【0，1】中，nu越大，决策边界越平滑。 CvSVM:ONE_CLASS - 单分类器，所有的训练数据提取自同一个类里，然後SVM建立了一个分界线以分割该类在特征空间中所占区域和其它类在特征空间中所占区域。 CvSVM:EPS_SVR - 回归。 训练集中的特征向量和拟合出来的超平面的距离需要小于p。异常值惩罚因子C被采用。 CvSVM:NU_SVR - 回归；nu 代替了p kernel_type/核类型： CvSVM:LINEAR - 没有任何向映射至高维空间，线性区分（或回归）在原始特征空间中被完成，这是最快的选择。 d(x,y) = x?y = (x,y) CvSVM:POLY - 多项式核: d(x,y) = (gamma*(x?y)+coef0)degree CvSVM:RBF - 径向基，对于大多数情况都是一个较好的选择：d(x,y) = exp(-gamma*|x-y|2) CvSVM:SIGMOID - sigmoid函数被用作核函数: d(x,y) = tanh(gamma*(x?y)+coef0) degree, gamma, coef0：都是核函数的参数，具体的参见上面的核函数的方程。 C, nu, p：在一般的SVM优化求解时的参数。 class_weights: Optional weights, assigned to particular classes. They are multiplied by C and thus affect the mis-classification penalty for different classes. The larger weight, the larger penalty on mis-classification of data from the certain class. term_crit: Termination procedure for iterative SVM training procedure (which solves a partial case of constrained quadratic optimization problem) The structure needs to be initialized and passed to the training method of CvSVM CvSVM:train训练SVM bool CvSVM:train( const CvMat* _train_data, const CvMat* _responses, const CvMat* _var_idx=0, const CvMat* _sample_idx=0, CvSVMParams _params=CvSVMParams() );The method trains SVM model. It follows the conventions of generic train method with the following limitations: only CV_ROW_SAMPLE data layout is supported, the input variables are all ordered, the output variables can be either categorical (_params.svm_type=CvSVM:C_SVC or _params.svm_type=CvSVM:NU_SVC) or ordered (_params.svm_type=CvSVM:EPS_SVR or _params.svm_type=CvSVM:NU_SVR) or not required at all (_params.svm_type=CvSVM:ONE_CLASS), missing measurements are not supported. 所有的参数都被集成在CvSVMParams这个结构中。 CvSVM:get_support_vector*得到支撑矢量和特殊矢量的数 int CvSVM:get_support_vector_count() const;const float* CvSVM:get_support_vector(int i) const;这个方法可以被用来得到支撑矢量的集合。 补充:在WindowsXP+OpenCVRC1平台下整合OpenCV与libSVM 虽然从RC1版开始opencv开始增设ML类，提供对常见的分类器和回归算法的支持。但是尚存在一些问题，比如说例子少（官方许诺说很快会提供一批新例子,见CVS版）。单说SVM这种算法，它自己提供了一套比较完备的函数，但是并不见得优于老牌的libsvm（它也应该参考过libsvm，至于是否效率优于libsvm，我并没有作过测试，官方也没有什么说法，但是libsvm持续开源更新，是公认的现存的开源SVM库中最易上手，性能最好的库）。所以在你的程序里整合opencv和libSVM还是一种比较好的解决方案。在VC中整合有些小地方需要注意，这篇文档主要是提供把图象作为SVM输入时程序遇到的这些小问题的解决方案。希望大家遇到问题时，多去读libSVM的源码，它本身也是开源的，C代码写得也很优秀，数据结构定义得也比较好。 首先是SVM的训练，这部分我并没有整合到VC里，直接使用它提供的python程序，虽然网格搜索这种简易搜索也可以自己写，但是识别时只需要训练生成的SVMmodel文件即可，所以可以和主程序分离开。至于python在windows下的使用，还是要设置一下的，首先要在系统环境变量path里把python的路径设进去，libsvm画训练等高线图还需要gnuplot的支持，打开python源程序(grid.py)，把gnuplot_exe设置成你机器里的安装路径，这样才能正确的运行程序。然後就是改步长和搜索范围，官方建议是先用大步长搜索，搜到最优值後再用小步长在小范围内搜索（我觉得这有可能会陷入局部最优，不过近似出的结果还可以接受）。我用的python版本是2.4，gnuplot4.0。 常用libSVM资料链接 官方站点，有一些tutorial和测试数据 哈工大的机器学习论坛，非常好 上交的一个研究生还写过libsvm2.6版的代码中文注释，源链接找不着了，大家自己搜搜吧，写得很好，上海交通大学模式分析与机器智能实验室。本文来自CSDN博客，转载请标明出处：/xiaoxin_ling/archive/2009/04/16/4083506.aspxSVM 支持向量机 SMO算法 实现 机器学习 如果对SVM原理不是很懂的，可以先看一下入门的视频，对帮助理解很有用的，然后再深入一点可以看看这几篇入门文章，作者写得挺详细，看完以后SVM的基础就了解得差不多了，再然后买本支持向量机导论作者是Nello Cristianini 和 John Shawe-Taylor，电子工业出版社的。然后把书本后面的那个SMO算法实现就基本上弄懂了SVM是怎么一回事，最后再编写一个SVM库出来，比如说像libsvm等工具使用，呵呵，差不多就这样。这些是我学习SVM的整个过程，也算是经验吧。 下面是SVM的简化版SMO算法，我将结合Java代码来解释一下整个SVM的学习训练过程，即所谓的train训练过程。那么什么是SMO算法呢？ SMO算法的目的无非是找出一个函数f(x)，这个函数能让我们把输入的数据x进行分类。既然是分类肯定需要一个评判的标准，比如分出来有两种情况A和B,那么怎么样才能说x是属于A类的，或不是B类的呢？就是需要有个边界，就好像两个国家一样有边界，如果边界越明显，则就越容易区分，因此，我们的目标是最大化边界的宽度，使得非常容易的区分是A类还是B类。 在SVM中，要最大化边界则需要最小化这个数值： w：是参量，值越大边界越明显 C代表惩罚系数，即如果某个x是属于某一类，但是它偏离了该类，跑到边界上后者其他类的地方去了，C越大表明越不想放弃这个点，边界就会缩小代表：松散变量但问题似乎还不好解，又因为SVM是一个凸二次规划问题，凸二次规划问题有最优解，于是问题转换成下列形式（KKT条件）： (1)这里的ai是拉格朗日乘子(问题通过拉格朗日乘法数来求解)对于（a）的情况，表明ai是正常分类，在边界内部（我们知道正确分类的点yi*f(xi)=0）对于（b）的情况，表明了ai是支持向量，在边界上对于（c）的情况，表明了ai是在两条边界之间而最优解需要满足KKT条件，即满足（a）（b）（c）条件都满足以下几种情况出现将会出现不满足：yiui=1但是ai=1但是ai0则是不满足的而原本ai=0yiui=1但是ai=0或者ai=C则表明不满足的，而原本应该是0aiC所以要找出不满足KKT的这些ai，并更新这些ai，但这些ai又受到另外一个约束，即 因此，我们通过另一个方法，即同时更新ai和aj，满足以下等式 就能保证和为0的约束。利用yiai+yjaj=常数，消去ai，可得到一个关于单变量aj的一个凸二次规划问题，不考虑其约束0=aj=C,可以得其解为： (2)这里(3)表示旧值，然后考虑约束0=aj=C可得到a的解析解为： (4)对于那么如何求得ai和aj呢？对于ai，即第一个乘子，可以通过刚刚说的那几种不满足KKT的条件来找，第二个乘子aj可以找满足条件 （5）b的更新：在满足条件：下更新b。（6） 最后更新所有ai，y和b，这样模型就出来了，然后通过函数： （7）输入是x，是一个数组，组中每一个值表示一个特征。输出是A类还是B类。（正类还是负类）以下是主要的代码段： view plaincopy to clipboardprint?01./* 02. * 默认输入参数值 03. * C: regularization parameter 04. * tol: numerical tolerance 05. * max passes 06. */ 07.double C = 1; /对不在界内的惩罚因子 08.double tol = 0.01;/容忍极限值 09.int maxPasses = 5; /表示没有改变拉格朗日乘子的最多迭代次数 10. 11./* 12. * 初始化a, b, passes 13. */ 14. 15.double a = new doublex.length;/拉格朗日乘子 16.this.a = a; 17. 18./将乘子初始化为0 19.for (int i = 0; i x.length; i+) 20. ai = 0; 21. 22.int passes = 0; 23. 24. 25.while (passes maxPasses) 26. /表示改变乘子的次数（基本上是成对改变的） 27. int num_changed_alphas = 0; 28. for (int i = 0; i = 1 and alpha = 0 (正确分类) 36. * yi*f(i) = 1 and 0alpha C (在边界上的支持向量) 37. * yi*f(i) = 0 42. * 如果ri 0并且alpha C 则违反了KKT条件 43. * 因为原本ri 0并且alpha 0则违反了KKT条件 45. * 因为原本ri 0对应的应该是alpha =0 46. */ 47. if (yi * Ei -tol & ai tol & ai 0) 49. 50. /* 51. * ui*yi=1边界上的点 0 ai 0) 60. /参照公式(5) 61. j = findMax(Ei, this.boundAlpha); 62. else 63. /如果边界上没有，就随便选一个j != i的aj 64. j = RandomSelect(i); 65. 66. double Ej = getE(j); 67. 68. /保存当前的ai和aj 69. double oldAi = ai; 70. double oldAj = aj; 71. 72. /* 73. * 计算乘子的范围U, V 74. * 参考公式(4) 75. */ 76. double L, H; 77. if (yi != yj) 78. L = Math.max(0, aj - ai); 79. H = Math.min(C, C - ai + aj); 80. else 81. L = Math.max(0, ai + aj - C); 82. H = Math.min(0, ai + aj); 83. 84. 85. 86. /* 87. * 如果eta等于0或者大于0 则表明a最优值应该在L或者U上 88. */ 89. double eta = 2 * k(i, j) - k(i, i) - k(j, j);/公式(3) 90. 91. if (eta = 0) 92. continue; 93. 94. aj = aj - yj * (Ei - Ej)/ eta;/公式(2) 95. if (0 aj & aj C) 96. this.boundAlpha.add(j); 97. 98. if (aj H) 101. aj = H; 102. 103. if (Math.abs(aj - oldAj) 1e-5) 104. continue; 105. ai = ai + yi * yj * (oldAj - aj); 106. if (0 ai & ai C) 107. this.boundAlpha.add(i); 108. 109. /* 110. * 计算b1， b2 111. * 参照公式(6) 112. */ 113. double b1 = b - Ei - yi * (ai - oldAi) * k(i, i) - yj * (aj - oldAj) * k(i, j); 114. double b2 = b - Ej - yi * (ai - oldAi) * k(i, j) - yj * (aj - oldAj) * k(j, j); 115. 116. if (0 ai & ai C) 117. b = b1; 118. else if (0 aj & aj C) 119. b = b2; 120. else 121. b = (b1 + b2) / 2; 122. 123. num_changed_alphas = num_changed_alphas + 1; 124. 125. 126. if (num_changed_alphas = 0) 127. passes+; 128. else 129. passes = 0; 130. 131. 132.return new SVMModel(a, y, b); /* * 默认输入参数值 * C: regularization parameter * tol: numerical tolerance * max passes */double C = 1; /对不在界内的惩罚因子double tol = 0.01;/容忍极限值int maxPasses = 5; /表示没有改变拉格朗日乘子的最多迭代次数/* * 初始化a, b, passes */double a = new doublex.length;/拉格朗日乘子this.a = a;/将乘子初始化为0for (int i = 0; i x.length; i+) ai = 0;int passes = 0;while (passes maxPasses) /表示改变乘子的次数（基本上是成对改变的）int num_changed_alphas = 0;for (int i = 0; i = 1 and alpha = 0 (正确分类) * yi*f(i) = 1 and 0alpha C (在边界上的支持向量) * yi*f(i) = 0 * 如果ri 0并且alpha C 则违反了KKT条件 * 因为原本ri 0并且alpha 0则违反了KKT条件 * 因为原本ri 0对应的应该是alpha =0 */if (yi * Ei -tol & ai tol & ai 0) /* * ui*yi=1边界上的点 0 ai 0) /参照公式(5)j = findMax(Ei, this.boundAlpha); else /如果边界上没有，就随便选一个j != i的ajj = RandomSelect(i);double Ej = getE(j);/保存当前的ai和ajdouble oldAi = ai;double oldAj = aj;/* * 计算乘子的范围U, V * 参考公式(4) */double L, H;if (yi != yj) L = Math.max(0, aj - ai);H = Math.min(C, C - ai + aj); else L = Math.max(0, ai + aj - C);H = Math.min(0, ai + aj);/* * 如果eta等于0或者大于0 则表明a最优值应