最近无所事事,学习了SVM,已经学习了一段时间,做了一些仿真实验,看了一下源码,自以为牛逼了,想写点东西。呵呵,其实不是啦,是有了一点心得,所以记录下来,给其他学习SVM的人抛砖引玉而已。算是做点小小贡献。好了,废话就不多说了,正式进入主题。
最近看SVM的源码都是看LIBSVM的源码,里面最重要的就是SVM.cpp文件。那么函数的调用过程以及一些比较重要的函数解释,在这个博客中已经写得比较详细了,我就不进行多说,可以参考(http://blog.csdn.net/flydreamGG/article/details/4470121)。
那么我这篇博文主要补充的是上述博客中一些写得不足的要点,总结一下。
在LibSVM的调用过程中,会调用到一个函数select_working_set(i,j) 。根据上海交通大学模式分析与机器智能实验室(http://www.pami.sjtu.edu.cn/people/gpliu/document/libsvm_src.pdf)的程序解释文档,对该函数进行了一个比较不给力的解释(以下简称上交的结论),文中的解释如下:
看了文档给出的解释,然后再看看SVM.cpp中的解释,你会发现,无论如何都对不上。貌似不是如上述所说的选择方案。后来进一步参考另外一个牛人博客(http://www.cnblogs.com/vivounicorn/archive/2011/08/25/2152824.html),博文给出的解释已经与SVM.cpp中的select_working_set函数差不多了,解释也比较清楚,但是还是有一些差异。最后还是参考了LibSVM作者发表的文章(Working Set Selection Using Second Order Information for TrainingSupport Vector Machines),里面给出了其详细解释,因此问题也最终得以解决。
一般SVM对偶问题为:
SVM收敛的充分必要条件是KKT条件,其表现为:
对(1)式进行求导可以得到:
那么接下来进行一些推导咯:
(1),由(2)和(3)可得:
(2),由(2)和(3)可得:
(3),由(2)和(3)可得:
(4),由(2)和(3)可得:
好了,至此就推导出了上海交大那个不给力的结论。那么我可以发现,(4)和(5)都是b大于某个数,(6)和(7)都是b小于某个数。那么我可以发现,(4)和(5)都是b大于某个数,(6)和(7)都是b小于某个数。因为b是个常量,那么根据上述条件,我们可以得到以下结论,在合理的下,有:
那么我们就是要从中挑选违反上述条件的,来进行重新的迭代和更新,使得所有的都满足上述条件。那么我们可以很容易得到违反条件为:
那么其实上交的结论在LibSVM作者的文章中被进一步简化成了如下形式:
那么,根据(8)的结论中关于i的选择,再看看select_working_set函数中关于i的选择,就可以明白了。嘿嘿。
for(int t=0; t<active_size; t++)
if(y[t]==+1)
{
if(!is_upper_bound(t)) //对应于
if(-G[t]>= Gmax)
{
Gmax= -G[t]; //对应于
Gmax_idx= t;
}
}
else
{
if(!is_lower_bound(t)) //对应于
if(G[t] >=Gmax)
{
Gmax= G[t]; //对应于
Gmax_idx= t;
}
}
补充几句:在SVM.cpp中G存储的就是
好了,i的问题就解决了!
在进入解决j的选择问题前,插入一个定理:
定理一:如果为半正定矩阵,当且仅当待优化乘子为“违反对”时,对于SMO类算法,函数是严格递减的。
LibSVM在做选择的时候,采用的是second order information来选择的。那么我们挑选出了i之后,剩下的任务就是挑选出既是“违反对”同时使目标函数值最小。补充一下:挑选了“违反对”,自然就使得目标函数自然递减了,那么我们挑选目标函数最小,自然使得迭代速度加快啦。这是我们希望看到的结果。
这时候,我们就要用到大名鼎鼎的泰勒展开式啦:
那么我们的优化问题变成了:
由约束条件可知:
补充两句:因为具有取值范围,所以当取到极值的时候,d的取值是有限制的,使得最终的的取值不会超出的取值范围。
那么最终的子优化问题变成了:
那么,经证明上述问题能够取到最小值
证明如下:
结论貌似挺复杂的,其实不然,仔细发现其实(13)就是一元二次函数,就是求极值,变成了一个初中数学题目了。呵呵。根据初中的知识,就能求得该函数的最小值。
源码解释:
for(int j=0; j<active_size; j++)
{
if(y[j]==+1)
{
if(!is_lower_bound(j))
{
//计算
double grad_diff=Gmax+G[j];
if (G[j]>= Gmax2)
Gmax2= G[j];
if (grad_diff> 0) //保证
{
doubleobj_diff;
doublequad_coef=Q_i[i]+QD[j]-2.0*y[i]*Q_i[j];
if(quad_coef > 0)
obj_diff= -(grad_diff*grad_diff)/quad_coef; //代入计算公式最小值
else
obj_diff= -(grad_diff*grad_diff)/TAU;
if(obj_diff <= obj_diff_min)
{
Gmin_idx=j;
obj_diff_min= obj_diff;
}
}
}
}
else
{
if (!is_upper_bound(j))
{
// 计算
doublegrad_diff= Gmax-G[j];
if (-G[j]>= Gmax2)
Gmax2= -G[j];
if (grad_diff> 0) //保证
{
doubleobj_diff;
doublequad_coef=Q_i[i]+QD[j]+2.0*y[i]*Q_i[j];
if (quad_coef> 0)
obj_diff= -(grad_diff*grad_diff)/quad_coef; //代入计算公式最小值
else
obj_diff= -(grad_diff*grad_diff)/TAU;
if(obj_diff <= obj_diff_min)
{
Gmin_idx=j;
obj_diff_min= obj_diff;
}
}
}
}