精华内容
下载资源
问答
  • 说来惭愧,断更快半个月了,本打算是一周一篇的。...本篇文章参考了诸多大牛的文章写成的,对于什么是SVM做出了生动的阐述,同时也进行了线性SVM的理论推导,以及最后的编程实践,公式较多,还需静下心来一点一点推导。

    转载请注明作者和出处: https://zhuanlan.zhihu.com/ml-jack
    机器学习知乎专栏:https://zhuanlan.zhihu.com/ml-jack
    CSDN博客专栏:http://blog.csdn.net/column/details/16415.html
    Github代码获取:https://github.com/Jack-Cherish/Machine-Learning/
    Python版本: Python3.x
    运行平台: Windows
    IDE: Sublime text3

    一 前言

    说来惭愧,断更快半个月了,本打算是一周一篇的。感觉SVM瞬间难了不少,推导耗费了很多时间,同时身边的事情也不少,忙了许久。本篇文章参考了诸多大牛的文章写成的,对于什么是SVM做出了生动的阐述,同时也进行了线性SVM的理论推导,以及最后的编程实践,公式较多,还需静下心来一点一点推导。

    本文出现的所有代码,均可在我的github上下载,欢迎Follow、Star:https://github.com/Jack-Cherish/Machine-Learning

    二 什么是SVM?

    SVM的英文全称是Support Vector Machines,我们叫它支持向量机。支持向量机是我们用于分类的一种算法。让我们以一个小故事的形式,开启我们的SVM之旅吧。

    在很久以前的情人节,一位大侠要去救他的爱人,但天空中的魔鬼和他玩了一个游戏。

    魔鬼在桌子上似乎有规律放了两种颜色的球,说:“你用一根棍分开它们?要求:尽量在放更多球之后,仍然适用。”

    于是大侠这样放,干的不错?

    然后魔鬼,又在桌上放了更多的球,似乎有一个球站错了阵营。显然,大侠需要对棍做出调整。

    SVM就是试图把棍放在最佳位置,好让在棍的两边有尽可能大的间隙。这个间隙就是球到棍的距离。

    现在好了,即使魔鬼放了更多的球,棍仍然是一个好的分界线。

    魔鬼看到大侠已经学会了一个trick(方法、招式),于是魔鬼给了大侠一个新的挑战。

    现在,大侠没有棍可以很好帮他分开两种球了,现在怎么办呢?当然像所有武侠片中一样大侠桌子一拍,球飞到空中。然后,凭借大侠的轻功,大侠抓起一张纸,插到了两种球的中间。

    现在,从空中的魔鬼的角度看这些球,这些球看起来像是被一条曲线分开了。

    再之后,无聊的大人们,把这些球叫做data,把棍子叫做classifier, 找到最大间隙的trick叫做optimization,拍桌子叫做kernelling, 那张纸叫做hyperplane

    更为直观地感受一下吧(需要翻墙):https://www.youtube.com/watch?v=3liCbRZPrZA

    概述一下:

    当一个分类问题,数据是线性可分的,也就是用一根棍就可以将两种小球分开的时候,我们只要将棍的位置放在让小球距离棍的距离最大化的位置即可,寻找这个最大间隔的过程,就叫做最优化。但是,现实往往是很残酷的,一般的数据是线性不可分的,也就是找不到一个棍将两种小球很好的分类。这个时候,我们就需要像大侠一样,将小球拍起,用一张纸代替小棍将小球进行分类。想要让数据飞起,我们需要的东西就是核函数(kernel),用于切分小球的纸,就是超平面。

    也许这个时候,你还是似懂非懂,没关系。根据刚才的描述,可以看出,问题是从线性可分延伸到线性不可分的。那么,我们就按照这个思路,进行原理性的剖析。

    三 线性SVM

    先看下线性可分的二分类问题。

    上图中的(a)是已有的数据,红色和蓝色分别代表两个不同的类别。数据显然是线性可分的,但是将两类数据点分开的直线显然不止一条。上图的(b)和©分别给出了B、C两种不同的分类方案,其中黑色实线为分界线,术语称为“决策面”。每个决策面对应了一个线性分类器。虽然从分类结果上看,分类器A和分类器B的效果是相同的。但是他们的性能是有差距的,看下图:

    在"决策面"不变的情况下,我又添加了一个红点。可以看到,分类器B依然能很好的分类结果,而分类器C则出现了分类错误。显然分类器B的"决策面"放置的位置优于分类器C的"决策面"放置的位置,SVM算法也是这么认为的,它的依据就是分类器B的分类间隔比分类器C的分类间隔大。这里涉及到第一个SVM独有的概念"分类间隔"。在保证决策面方向不变且不会出现错分样本的情况下移动决策面,会在原来的决策面两侧找到两个极限位置(越过该位置就会产生错分现象),如虚线所示。虚线的位置由决策面的方向和距离原决策面最近的几个样本的位置决定。而这两条平行虚线正中间的分界线就是在保持当前决策面方向不变的前提下的最优决策面。两条虚线之间的垂直距离就是这个最优决策面对应的分类间隔。显然每一个可能把数据集正确分开的方向都有一个最优决策面(有些方向无论如何移动决策面的位置也不可能将两类样本完全分开),而不同方向的最优决策面的分类间隔通常是不同的,那个具有“最大间隔”的决策面就是SVM要寻找的最优解。而这个真正的最优解对应的两侧虚线所穿过的样本点,就是SVM中的支持样本点,称为"支持向量"。

    1 数学建模

    求解这个"决策面"的过程,就是最优化。一个最优化问题通常有两个基本的因素:1)目标函数,也就是你希望什么东西的什么指标达到最好;2)优化对象,你期望通过改变哪些因素来使你的目标函数达到最优。在线性SVM算法中,目标函数显然就是那个"分类间隔",而优化对象则是决策面。所以要对SVM问题进行数学建模,首先要对上述两个对象(“分类间隔"和"决策面”)进行数学描述。按照一般的思维习惯,我们先描述决策面。

    数学建模的时候,先在二维空间建模,然后再推广到多维。

    (1)"决策面"方程

    我们都知道二维空间下一条直线的方式如下所示:

    现在我们做个小小的改变,让原来的x轴变成x1,y轴变成x2。

    移项得:

    将公式向量化得:

    进一步向量化,用w列向量和x列向量和标量γ进一步向量化:

    其中,向量w和x分别为:

    这里w1=a,w2=-1。我们都知道,最初的那个直线方程a和b的几何意义,a表示直线的斜率,b表示截距,a决定了直线与x轴正方向的夹角,b决定了直线与y轴交点位置。那么向量化后的直线的w和r的几何意义是什么呢?

    现在假设:

    可得:

    在坐标轴上画出直线和向量w:

    蓝色的线代表向量w,红色的先代表直线y。我们可以看到向量w和直线的关系为垂直关系。这说明了向量w也控制这直线的方向,只不过是与这个直线的方向是垂直的。标量γ的作用也没有变,依然决定了直线的截距。此时,我们称w为直线的法向量。

    二维空间的直线方程已经推导完成,将其推广到n为空间,就变成了超平面方程。(一个超平面,在二维空间的例子就是一个直线)但是它的公式没变,依然是:

    不同之处在于:

    我们已经顺利推导出了"决策面"方程,它就是我们的超平面方程,之后,我们统称其为超平面方程。

    ###(2)"分类间隔"方程

    现在,我们依然对于一个二维平面的简单例子进行推导。

    我们已经知道间隔的大小实际上就是支持向量对应的样本点到决策面的距离的二倍。那么图中的距离d我们怎么求?我们高中都学过,点到直线的距离距离公式如下:

    公式中的直线方程为Ax0+By0+C=0,点P的坐标为(x0,y0)。

    现在,将直线方程扩展到多维,求得我们现在的超平面方程,对公式进行如下变形:

    这个d就是"分类间隔"。其中||w||表示w的二范数,求所有元素的平方和,然后再开方。比如对于二维平面:

    那么,

    我们目的是为了找出一个分类效果好的超平面作为分类器。分类器的好坏的评定依据是分类间隔W=2d的大小,即分类间隔W越大,我们认为这个超平面的分类效果越好。此时,求解超平面的问题就变成了求解分类间隔W最大化的为题。W的最大化也就是d最大化的。

    ###(3)约束条件

    看起来,我们已经顺利获得了目标函数的数学形式。但是为了求解w的最大值。我们不得不面对如下问题:

    • 我们如何判断超平面是否将样本点正确分类?
    • 我们知道相求距离d的最大值,我们首先需要找到支持向量上的点,怎么在众多的点中选出支持向量上的点呢?

    上述我们需要面对的问题就是约束条件,也就是说我们优化的变量d的取值范围受到了限制和约束。事实上约束条件一直是最优化问题里最让人头疼的东西。但既然我们已经知道了这些约束条件确实存在,就不得不用数学语言对他们进行描述。但SVM算法通过一些巧妙的小技巧,将这些约束条件融合到一个不等式里面。

    这个二维平面上有两种点,我们分别对它们进行标记:

    • 红颜色的圆点标记为1,我们人为规定其为正样本;
    • 蓝颜色的五角星标记为-1,我们人为规定其为负样本。

    对每个样本点xi加上一个类别标签yi:

    如果我们的超平面方程能够完全正确地对上图的样本点进行分类,就会满足下面的方程:

    如果我们要求再高一点,假设决策面正好处于间隔区域的中轴线上,并且相应的支持向量对应的样本点到决策面的距离为d,那么公式进一步写成:

    上述公式的解释就是,对于所有分类标签为1的样本点,它们到直线的距离都大于等于d(支持向量上的样本点到超平面的距离)。对于所有分类标签为-1的样本点,它们到直线的距离都小于等于d。公式两边都除以d,就可以得到:

    其中,

    因为||w||和d都是标量。所上述公式的两个矢量,依然描述一条直线的法向量和截距。

    上述两个公式,都是描述一条直线,数学模型代表的意义是一样的。现在,让我们对wd和γd重新起个名字,就叫它们w和γ。因此,我们就可以说:“对于存在分类间隔的两类样本点,我们一定可以找到一些超平面面,使其对于所有的样本点均满足下面的条件:”

    上述方程即给出了SVM最优化问题的约束条件。这时候,可能有人会问了,为什么标记为1和-1呢?因为这样标记方便我们将上述方程变成如下形式:

    正是因为标签为1和-1,才方便我们将约束条件变成一个约束方程,从而方便我们的计算。

    (4)线性SVM优化问题基本描述

    现在整合一下思路,我们已经得到我们的目标函数:

    我们的优化目标是是d最大化。我们已经说过,我们是用支持向量上的样本点求解d的最大化的问题的。那么支持向量上的样本点有什么特点呢?

    你赞同这个观点吗?所有支持向量上的样本点,都满足如上公式。如果不赞同,请重看"分类间隔"方程推导过程。

    现在我们就可以将我们的目标函数进一步化简:

    因为,我们只关心支持向量上的点。随后我们求解d的最大化问题变成了||w||的最小化问题。进而||w||的最小化问题等效于

    为什么要做这样的等效呢?这是为了在进行最优化的过程中对目标函数求导时比较方便,但这绝对不影响最优化问题最后的求解。我们将最终的目标函数和约束条件放在一起进行描述:

    这里n是样本点的总个数,缩写s.t.表示"Subject to",是"服从某某条件"的意思。上述公式描述的是一个典型的不等式约束条件下的二次型函数优化问题,同时也是支持向量机的基本数学模型。

    ###(5)求解准备

    我们已经得到支持向量机的基本数学模型,接下来的问题就是如何根据数学模型,求得我们想要的最优解。在学习求解方法之前,我们得知道一点,想用我下面讲述的求解方法有一个前提,就是我们的目标函数必须是凸函数。理解凸函数,我们还要先明确另一个概念,凸集。在凸几何中,凸集(convex set)是在)凸组合下闭合的放射空间的子集。看一幅图可能更容易理解:

    左右量图都是一个集合。**如果集合中任意2个元素连线上的点也在集合中,那么这个集合就是凸集。**显然,上图中的左图是一个凸集,上图中的右图是一个非凸集。

    凸函数的定义也是如此,其几何意义表示为函数任意两点连线上的值大于对应自变量处的函数值。若这里凸集C即某个区间L,那么,设函数f为定义在区间L上的函数,若对L上的任意两点x1,x2和任意的实数λ,λ属于(0,1),总有:

    则函数f称为L上的凸函数,当且仅当其上镜图(在函数图像上方的点集)为一个凸集。再看一幅图,也许更容易理解:

    像上图这样的函数,它整体就是一个非凸函数,我们无法获得全局最优解的,只能获得局部最优解。比如红框内的部分,如果单独拿出来,它就是一个凸函数。对于我们的目标函数:

    很显然,它是一个凸函数。所以,可以使用我接下来讲述的方法求取最优解。

    通常我们需要求解的最优化问题有如下几类:

    • 无约束优化问题,可以写为:

    - 有等式约束的优化问题,可以写为:

    - 有不等式约束的优化问题,可以写为:

    对于第(a)类的优化问题,尝尝使用的方法就是费马大定理(Fermat),即使用求取函数f(x)的导数,然后令其为零,可以求得候选最优值,再在这些候选值中验证;如果是凸函数,可以保证是最优解。这也就是我们高中经常使用的求函数的极值的方法。

    对于第(b)类的优化问题,常常使用的方法就是拉格朗日乘子法(Lagrange Multiplier) ,即把等式约束h_i(x)用一个系数与f(x)写为一个式子,称为拉格朗日函数,而系数称为拉格朗日乘子。通过拉格朗日函数对各个变量求导,令其为零,可以求得候选值集合,然后验证求得最优值。

    对于第©类的优化问题,常常使用的方法就是KKT条件。同样地,我们把所有的等式、不等式约束与f(x)写为一个式子,也叫拉格朗日函数,系数也称拉格朗日乘子,通过一些条件,可以求出最优值的必要条件,这个条件称为KKT条件。

    必要条件和充要条件如果不理解,可以看下面这句话:

    • A的必要条件就是A可以推出的结论
    • A的充分条件就是可以推出A的前提

    了解到这些,现在让我们再看一下我们的最优化问题:

    现在,我们的这个对优化问题属于哪一类?很显然,它属于第©类问题。因为,在学习求解最优化问题之前,我们还要学习两个东西:拉格朗日函数和KKT条件。

    (6)拉格朗日函数

    首先,我们先要从宏观的视野上了解一下拉格朗日对偶问题出现的原因和背景。

    我们知道我们要求解的是最小化问题,所以一个直观的想法是如果我能够构造一个函数,使得该函数在可行解区域内与原目标函数完全一致,而在可行解区域外的数值非常大,甚至是无穷大,那么这个没有约束条件的新目标函数的优化问题就与原来有约束条件的原始目标函数的优化问题是等价的问题。这就是使用拉格朗日方程的目的,它将约束条件放到目标函数中,从而将有约束优化问题转换为无约束优化问题。

    随后,人们又发现,使用拉格朗日获得的函数,使用求导的方法求解依然困难。进而,需要对问题再进行一次转换,即使用一个数学技巧:拉格朗日对偶。

    所以,显而易见的是,我们在拉格朗日优化我们的问题这个道路上,需要进行下面二个步骤:

    • 将有约束的原始目标函数转换为无约束的新构造的拉格朗日目标函数
    • 使用拉格朗日对偶性,将不易求解的优化问题转化为易求解的优化

    下面,进行第一步:将有约束的原始目标函数转换为无约束的新构造的拉格朗日目标函数

    公式变形如下:

    其中αi是拉格朗日乘子,αi大于等于0,是我们构造新目标函数时引入的系数变量(我们自己设置)。现在我们令:

    当样本点不满足约束条件时,即在可行解区域外

    此时,我们将αi设置为正无穷,此时θ(w)显然也是正无穷。

    当样本点满足约束条件时,即在可行解区域内:

    此时,显然θ(w)为原目标函数本身。我们将上述两种情况结合一下,就得到了新的目标函数:

    此时,再看我们的初衷,就是为了建立一个在可行解区域内与原目标函数相同,在可行解区域外函数值趋近于无穷大的新函数,现在我们做到了。

    现在,我们的问题变成了求新目标函数的最小值,即:

    这里用p*表示这个问题的最优值,且和最初的问题是等价的。

    接下来,我们进行第二步:将不易求解的优化问题转化为易求解的优化

    我们看一下我们的新目标函数,先求最大值,再求最小值。这样的话,我们首先就要面对带有需要求解的参数w和b的方程,而αi又是不等式约束,这个求解过程不好做。所以,我们需要使用拉格朗日函数对偶性,将最小和最大的位置交换一下,这样就变成了:

    交换以后的新问题是原始问题的对偶问题,这个新问题的最优值用d来表示。而且d<=p*。我们关心的是d=p的时候,这才是我们要的解。需要什么条件才能让d=p呢?

    • 首先必须满足这个优化问题是凸优化问题。
    • 其次,需要满足KKT条件。

    凸优化问题的定义是:**求取最小值的目标函数为凸函数的一类优化问题。**目标函数是凸函数我们已经知道,这个优化问题又是求最小值。所以我们的最优化问题就是凸优化问题。

    接下里,就是探讨是否满足KKT条件了。

    ###(7)KKT条件

    我们已经使用拉格朗日函数对我们的目标函数进行了处理,生成了一个新的目标函数。通过一些条件,可以求出最优值的必要条件,这个条件就是接下来要说的KKT条件。一个最优化模型能够表示成下列标准形式:

    KKT条件的全称是Karush-Kuhn-Tucker条件,KKT条件是说最优值条件必须满足以下条件:

    • 条件一:经过拉格朗日函数处理之后的新目标函数L(w,b,α)对α求导为零:
    • 条件二:h(x) = 0;
    • 条件三:α*g(x) = 0;

    对于我们的优化问题:

    显然,条件二已经满足了。另外两个条件为啥也满足呢?

    这里原谅我省略一系列证明步骤,感兴趣的可以移步这里:http://blog.csdn.net/xianlingmao/article/details/7919597

    这里已经给出了很好的解释。现在,凸优化问题和KKT都满足了,问题转换成了对偶问题。而求解这个对偶学习问题,可以分为三个步骤:首先要让L(w,b,α)关于w和b最小化,然后求对α的极大,最后利用SMO算法求解对偶问题中的拉格朗日乘子。

    ###(8)对偶问题求解

    第一步:

    根据上述推导已知:

    首先固定α,要让L(w,b,α)关于w和b最小化,我们分别对w和b偏导数,令其等于0,即:

    将上述结果带回L(w,b,α)得到:

    从上面的最后一个式子,我们可以看出,此时的L(w,b,α)函数只含有一个变量,即αi。

    第二步:

    现在内侧的最小值求解完成,我们求解外侧的最大值,从上面的式子得到

    现在我们的优化问题变成了如上的形式。对于这个问题,我们有更高效的优化算法,即序列最小优化(SMO)算法。我们通过这个优化算法能得到α,再根据α,我们就可以求解出w和b,进而求得我们最初的目的:找到超平面,即"决策平面"。

    总结一句话:我们为啥使出吃奶的劲儿进行推导?因为我们要将最初的原始问题,转换到可以使用SMO算法求解的问题,这是一种最流行的求解方法。为啥用这种求解方法?因为它牛逼啊!

    2 SMO算法

    现在,我们已经得到了可以用SMO算法求解的目标函数,但是对于怎么编程实现SMO算法还是感觉无从下手。那么现在就聊聊如何使用SMO算法进行求解。

    ###(1)Platt的SMO算法

    1996年,John Platt发布了一个称为SMO的强大算法,用于训练SVM。SM表示序列最小化(Sequential Minimal Optimizaion)。Platt的SMO算法是将大优化问题分解为多个小优化问题来求解的。这些小优化问题往往很容易求解,并且对它们进行顺序求解的结果与将它们作为整体来求解的结果完全一致的。在结果完全相同的同时,SMO算法的求解时间短很多。

    SMO算法的目标是求出一系列alpha和b,一旦求出了这些alpha,就很容易计算出权重向量w并得到分隔超平面。

    SMO算法的工作原理是:每次循环中选择两个alpha进行优化处理。一旦找到了一对合适的alpha,那么就增大其中一个同时减小另一个。这里所谓的"合适"就是指两个alpha必须符合以下两个条件,条件之一就是两个alpha必须要在间隔边界之外,而且第二个条件则是这两个alpha还没有进进行过区间化处理或者不在边界上。

    ###(2)SMO算法的解法

    先来定义特征到结果的输出函数为:

    接着,我们回忆一下原始优化问题,如下:

    求导得:

    将上述公式带入输出函数中:

    与此同时,拉格朗日对偶后得到最终的目标化函数:

    将目标函数变形,在前面增加一个符号,将最大值问题转换成最小值问题:

    实际上,对于上述目标函数,是存在一个假设的,即数据100%线性可分。但是,目前为止,我们知道几乎所有数据都不那么"干净"。这时我们就可以通过引入所谓的松弛变量(slack variable),来允许有些数据点可以处于超平面的错误的一侧。这样我们的优化目标就能保持仍然不变,但是此时我们的约束条件有所改变:

    根据KKT条件可以得出其中αi取值的意义为:

    • 对于第1种情况,表明αi是正常分类,在边界内部;
    • 对于第2种情况,表明αi是支持向量,在边界上;
    • 对于第3种情况,表明αi是在两条边界之间。

    而最优解需要满足KKT条件,即上述3个条件都得满足,以下几种情况出现将会不满足:

    也就是说,如果存在不能满足KKT条件的αi,那么需要更新这些αi,这是第一个约束条件。此外,更新的同时还要受到第二个约束条件的限制,即:

    因为这个条件,我们同时更新两个α值,因为只有成对更新,才能保证更新之后的值仍然满足和为0的约束,假设我们选择的两个乘子为α1和α2:

    其中, ksi为常数。因为两个因子不好同时求解,所以可以先求第二个乘子α2的解(α2 new),得到α2的解(α2 new)之后,再用α2的解(α2 new)表示α1的解(α1 new )。为了求解α2 new ,得先确定α2 new的取值范围。假设它的上下边界分别为H和L,那么有:

    接下来,综合下面两个条件:

    当y1不等于y2时,即一个为正1,一个为负1的时候,可以得到:

    所以有:

    此时,取值范围如下图所示:

    当y1等于y2时,即两个都为正1或者都为负1,可以得到:

    所以有:

    此时,取值范围如下图所示:

    如此,根据y1和y2异号或同号,可以得出α2 new的上下界分别为:

    这个界限就是编程的时候需要用到的。已经确定了边界,接下来,就是推导迭代式,用于更新 α值。

    我们已经知道,更新α的边界,接下来就是讨论如何更新α值。我们依然假设选择的两个乘子为α1和α2。固定这两个乘子,进行推导。于是目标函数变成了:

    [点击放大图片](https://img-blog.csdn.net/20170923181141690?watermark/2/text/aHR0cDovL2Jsb2cuY3Nkbi5uZXQvYzQwNjQ5NTc2Mg==/font/5a6L5L2T/fontsize/400/fill/I0JBQkFCMA==/dissolve/70/gravity/SouthEast)

    为了描述方便,我们定义如下符号:

    最终目标函数变为:

    我们不关心constant的部分,因为对于α1和α2来说,它们都是常数项,在求导的时候,直接变为0。对于这个目标函数,如果对其求导,还有个未知数α1,所以要推导出α1和α2的关系,然后用α2代替α1,这样目标函数就剩一个未知数了,我们就可以求导了,推导出迭代公式。所以现在继续推导α1和α2的关系。注意第一个约束条件:

    我们在求α1和α2的时候,可以将α3,α4,…,αn和y3,y4,…,yn看作常数项。因此有:

    我们不必关心常数B的大小,现在将上述等式两边同时乘以y1,得到(y1y1=1):

    其中γ为常数By1,我们不关心这个值,s=y1y2。接下来,我们将得到的α1带入W(α2)公式得:

    这样目标函数中就只剩下α2了,我们对其求偏导(注意:s=y1y2,所以s的平方为1,y1的平方和y2的平方均为1):

    继续化简,将s=y1y2带入方程。

    我们令:

    Ei为误差项,η为学习速率。

    再根据我们已知的公式:

    将α2 new继续化简得:

    这样,我们就得到了最终需要的迭代公式。这个是没有经过剪辑是的解,需要考虑约束:

    根据之前推导的α取值范围,我们得到最终的解析解为:

    又因为:

    消去γ得:

    这样,我们就知道了怎样计算α1和α2了,也就是如何对选择的α进行更新。

    当我们更新了α1和α2之后,需要重新计算阈值b,因为b关系到了我们f(x)的计算,也就关系到了误差Ei的计算。

    我们要根据α的取值范围,去更正b的值,使间隔最大化。当α1 new在0和C之间的时候,根据KKT条件可知,这个点是支持向量上的点。因此,满足下列公式:

    公式两边同时乘以y1得(y1y1=1):

    因为我们是根据α1和α2的值去更新b,所以单独提出i=1和i=2的时候,整理可得:

    其中前两项为:

    将上述两个公式,整理得:

    同理可得b2 new为:

    当b1和b2都有效的时候,它们是相等的,即:

    当两个乘子都在边界上,则b阈值和KKT条件一致。当不满足的时候,SMO算法选择他们的中点作为新的阈值:

    最后,更新所有的α和b,这样模型就出来了,从而即可求出我们的分类函数。

    现在,让我们梳理下SMO算法的步骤:

    • 步骤1:计算误差:

    - 步骤2:计算上下界L和H:

    - 步骤3:计算η:

    - 步骤4:更新αj:

    - 步骤5:根据取值范围修剪αj:

    - 步骤6:更新αi:

    - 步骤7:更新b1和b2:

    - 步骤8:根据b1和b2更新b:

    四 编程求解线性SVM

    已经梳理完了SMO算法实现步骤,接下来按照这个思路编写代码,进行实战练习。

    (1)可视化数据集

    我们先使用简单的数据集进行测试,数据集下载地址:https://github.com/Jack-Cherish/Machine-Learning/blob/master/SVM/testSet.txt

    编写程序可视化数据集,看下它是长什么样的:

    # -*- coding:UTF-8 -*-
    import matplotlib.pyplot as plt
    import numpy as np
    
    """
    函数说明:读取数据
    
    Parameters:
        fileName - 文件名
    Returns:
        dataMat - 数据矩阵
        labelMat - 数据标签
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-21
    """
    def loadDataSet(fileName):
        dataMat = []; labelMat = []
        fr = open(fileName)
        for line in fr.readlines():                                     #逐行读取,滤除空格等
            lineArr = line.strip().split('\t')
            dataMat.append([float(lineArr[0]), float(lineArr[1])])      #添加数据
            labelMat.append(float(lineArr[2]))                          #添加标签
        return dataMat,labelMat
    
    """
    函数说明:数据可视化
    
    Parameters:
        dataMat - 数据矩阵
        labelMat - 数据标签
    Returns:
        无
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-21
    """
    def showDataSet(dataMat, labelMat):
        data_plus = []                                  #正样本
        data_minus = []                                 #负样本
        for i in range(len(dataMat)):
            if labelMat[i] > 0:
                data_plus.append(dataMat[i])
            else:
                data_minus.append(dataMat[i])
        data_plus_np = np.array(data_plus)              #转换为numpy矩阵
        data_minus_np = np.array(data_minus)            #转换为numpy矩阵
        plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1])   #正样本散点图
        plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1]) #负样本散点图
        plt.show()
    
    if __name__ == '__main__':
        dataMat, labelMat = loadDataSet('testSet.txt')
        showDataSet(dataMat, labelMat)
    

    运行程序,查看结果:

    这就是我们使用的二维数据集,显然线性可分。现在我们使用简化版的SMO算法进行求解。

    (2)简化版SMO算法

    按照上述已经推导的步骤编写代码:

    # -*- coding:UTF-8 -*-
    from time import sleep
    import matplotlib.pyplot as plt
    import numpy as np
    import random
    import types
    
    """
    函数说明:读取数据
    
    Parameters:
        fileName - 文件名
    Returns:
        dataMat - 数据矩阵
        labelMat - 数据标签
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-21
    """
    def loadDataSet(fileName):
        dataMat = []; labelMat = []
        fr = open(fileName)
        for line in fr.readlines():                                     #逐行读取,滤除空格等
            lineArr = line.strip().split('\t')
            dataMat.append([float(lineArr[0]), float(lineArr[1])])      #添加数据
            labelMat.append(float(lineArr[2]))                          #添加标签
        return dataMat,labelMat
    
    
    """
    函数说明:随机选择alpha
    
    Parameters:
        i - alpha
        m - alpha参数个数
    Returns:
        j -
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-21
    """
    def selectJrand(i, m):
        j = i                                 #选择一个不等于i的j
        while (j == i):
            j = int(random.uniform(0, m))
        return j
    
    """
    函数说明:修剪alpha
    
    Parameters:
        aj - alpha值
        H - alpha上限
        L - alpha下限
    Returns:
        aj - alpah值
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-21
    """
    def clipAlpha(aj,H,L):
        if aj > H:
            aj = H
        if L > aj:
            aj = L
        return aj
    
    """
    函数说明:简化版SMO算法
    
    Parameters:
        dataMatIn - 数据矩阵
        classLabels - 数据标签
        C - 松弛变量
        toler - 容错率
        maxIter - 最大迭代次数
    Returns:
        无
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-23
    """
    def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
        #转换为numpy的mat存储
        dataMatrix = np.mat(dataMatIn); labelMat = np.mat(classLabels).transpose()
        #初始化b参数,统计dataMatrix的维度
        b = 0; m,n = np.shape(dataMatrix)
        #初始化alpha参数,设为0
        alphas = np.mat(np.zeros((m,1)))
        #初始化迭代次数
        iter_num = 0
        #最多迭代matIter次
        while (iter_num < maxIter):
            alphaPairsChanged = 0
            for i in range(m):
                #步骤1:计算误差Ei
                fXi = float(np.multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + b
                Ei = fXi - float(labelMat[i])
                #优化alpha,更设定一定的容错率。
                if ((labelMat[i]*Ei < -toler) and (alphas[i] < C)) or ((labelMat[i]*Ei > toler) and (alphas[i] > 0)):
                    #随机选择另一个与alpha_i成对优化的alpha_j
                    j = selectJrand(i,m)
                    #步骤1:计算误差Ej
                    fXj = float(np.multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[j,:].T)) + b
                    Ej = fXj - float(labelMat[j])
                    #保存更新前的aplpha值,使用深拷贝
                    alphaIold = alphas[i].copy(); alphaJold = alphas[j].copy();
                    #步骤2:计算上下界L和H
                    if (labelMat[i] != labelMat[j]):
                        L = max(0, alphas[j] - alphas[i])
                        H = min(C, C + alphas[j] - alphas[i])
                    else:
                        L = max(0, alphas[j] + alphas[i] - C)
                        H = min(C, alphas[j] + alphas[i])
                    if L==H: print("L==H"); continue
                    #步骤3:计算eta
                    eta = 2.0 * dataMatrix[i,:]*dataMatrix[j,:].T - dataMatrix[i,:]*dataMatrix[i,:].T - dataMatrix[j,:]*dataMatrix[j,:].T
                    if eta >= 0: print("eta>=0"); continue
                    #步骤4:更新alpha_j
                    alphas[j] -= labelMat[j]*(Ei - Ej)/eta
                    #步骤5:修剪alpha_j
                    alphas[j] = clipAlpha(alphas[j],H,L)
                    if (abs(alphas[j] - alphaJold) < 0.00001): print("alpha_j变化太小"); continue
                    #步骤6:更新alpha_i
                    alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])
                    #步骤7:更新b_1和b_2
                    b1 = b - Ei- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[i,:].T - labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[i,:]*dataMatrix[j,:].T
                    b2 = b - Ej- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[j,:].T - labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[j,:]*dataMatrix[j,:].T
                    #步骤8:根据b_1和b_2更新b
                    if (0 < alphas[i]) and (C > alphas[i]): b = b1
                    elif (0 < alphas[j]) and (C > alphas[j]): b = b2
                    else: b = (b1 + b2)/2.0
                    #统计优化次数
                    alphaPairsChanged += 1
                    #打印统计信息
                    print("第%d次迭代 样本:%d, alpha优化次数:%d" % (iter_num,i,alphaPairsChanged))
            #更新迭代次数
            if (alphaPairsChanged == 0): iter_num += 1
            else: iter_num = 0
            print("迭代次数: %d" % iter_num)
        return b,alphas
    
    """
    函数说明:分类结果可视化
    
    Parameters:
        dataMat - 数据矩阵
        w - 直线法向量
        b - 直线解决
    Returns:
        无
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-23
    """
    def showClassifer(dataMat, w, b):
        #绘制样本点
        data_plus = []                                  #正样本
        data_minus = []                                 #负样本
        for i in range(len(dataMat)):
            if labelMat[i] > 0:
                data_plus.append(dataMat[i])
            else:
                data_minus.append(dataMat[i])
        data_plus_np = np.array(data_plus)              #转换为numpy矩阵
        data_minus_np = np.array(data_minus)            #转换为numpy矩阵
        plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1], s=30, alpha=0.7)   #正样本散点图
        plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1], s=30, alpha=0.7) #负样本散点图
        #绘制直线
        x1 = max(dataMat)[0]
        x2 = min(dataMat)[0]
        a1, a2 = w
        b = float(b)
        a1 = float(a1[0])
        a2 = float(a2[0])
        y1, y2 = (-b- a1*x1)/a2, (-b - a1*x2)/a2
        plt.plot([x1, x2], [y1, y2])
        #找出支持向量点
        for i, alpha in enumerate(alphas):
            if abs(alpha) > 0:
                x, y = dataMat[i]
                plt.scatter([x], [y], s=150, c='none', alpha=0.7, linewidth=1.5, edgecolor='red')
        plt.show()
    
    
    """
    函数说明:计算w
    
    Parameters:
        dataMat - 数据矩阵
        labelMat - 数据标签
        alphas - alphas值
    Returns:
        无
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-09-23
    """
    def get_w(dataMat, labelMat, alphas):
        alphas, dataMat, labelMat = np.array(alphas), np.array(dataMat), np.array(labelMat)
        w = np.dot((np.tile(labelMat.reshape(1, -1).T, (1, 2)) * dataMat).T, alphas)
        return w.tolist()
    
    
    if __name__ == '__main__':
        dataMat, labelMat = loadDataSet('testSet.txt')
        b,alphas = smoSimple(dataMat, labelMat, 0.6, 0.001, 40)
        w = get_w(dataMat, labelMat, alphas)
        showClassifer(dataMat, w, b)
    

    程序运行结果:

    其中,中间的蓝线为求出来的分类器,用红圈圈出的点为支持向量点。

    五 总结

    • 本文主要进行了线性SVM的推导,并通过编程实现一个简化版SMO算法;
    • 本文的简化版SMO算法在选取α的时候,没有选择启发式的选择方法,并且没有两个乘子的计算没有进行优化,所以算法比较耗时,下一篇文章会讲解相应的优化方法;
    • 本文讨论的是线性SVM,没有使用核函数,下一篇文章将会讲解如何应用核函数,将SVM应用于非线性数据集;
    • 如有问题,请留言。如有错误,还望指正,谢谢!

    PS: 如果觉得本篇本章对您有所帮助,欢迎关注、评论、赞!

    参考资料:

    • [1] 五岁小孩也能看懂的SVM:https://www.zhihu.com/question/21094489/answer/8627319
    • [2] 五岁小孩也能看懂的SVM :https://www.reddit.com/r/MachineLearning/comments/15zrpp/please_explain_support_vector_machines_svm_like_i/
    • [3] pluskid大牛博客:http://blog.pluskid.org/?page_id=683
    • [4] 陈东岳老师文章:https://zhuanlan.zhihu.com/p/24638007
    • [5] 深入理解拉格朗日乘子法和KKT条件:http://blog.csdn.net/xianlingmao/article/details/7919597
    • [6] 充分条件和必要条件:https://www.zhihu.com/question/30469121
    • [7] 凸函数:https://zh.wikipedia.org/wiki/%E5%87%B8%E5%87%BD%E6%95%B0
    • [8]《机器学习实战》第6章内容。
    • [9] SVM之SMO算法:http://www.cnblogs.com/zangrunqiang/p/5515872.html

    1

    展开全文
  • 《机器学习实战》6.2支持向量机SVM基础之处理复杂的非线性SVM 搜索微信公众号:‘AI-ming3526’或者’计算机视觉这件小事’ 获取更多人工智能、机器学习干货 csdn:https://blog.csdn.net/baidu_31657889/ github:...

    《机器学习实战》6.2支持向量机SVM基础之处理复杂的非线性SVM

    搜索微信公众号:‘AI-ming3526’或者’计算机视觉这件小事’ 获取更多人工智能、机器学习干货
    csdn:https://blog.csdn.net/baidu_31657889/
    github:https://github.com/aimi-cn/AILearners

    本文出现的所有代码,均可在github上下载,不妨来个Star把谢谢~:Github代码地址

    一、非线性SVM

    1、核函数(kernel) 使用

    我们已经了解到,SVM如何处理线性可分的情况,而对于非线性的情况,SVM的处理方式就是选择一个核函数。简而言之:在线性不可分的情况下,SVM通过某种事先选择的非线性映射(核函数)将输入变量映到一个高维特征空间,将其变成在高维空间线性可分,在这个高维空间中构造最优分类超平面。也就是上一节最开始的大侠,凭借大侠的轻功,大侠抓起一张纸,插到了两种球的中间。让这些球感觉是被曲线分开了一样。

    根据上篇文章,线性可分的情况下,可知最终的超平面方程为:
    在这里插入图片描述

    将上述公式用内积来表示:
    在这里插入图片描述

    对于线性不可分,我们使用一个非线性映射,将数据映射到特征空间,在特征空间中使用线性学习器,分类函数变形如下:
    在这里插入图片描述

    其中ϕ是从输入空间(X)到某个特征空间(F)的映射,这意味着建立非线性学习器分为两步:

    • 首先使用一个非线性映射将数据变换到一个特征空间F;
    • 然后在特征空间使用线性分类学习器学习。

    如果有一种方法可以在特征空间中直接计算内积 <ϕ(xi),ϕ(x)>,就像在原始输入点的函数中一样,就有可能将两个步骤融合到一起建立一个分线性的学习器,这样直接计算的方法称为核函数方法

    这里直接给出一个定义:核是一个函数k,对所有x,z∈X,满足k(x,z)=<ϕ(xi),ϕ(x)>,这里ϕ(·)是从原始输入空间X到内积空间F的映射。

    简而言之:如果不是用核技术,就会先计算线性映ϕ(x1)和ϕ(x2),然后计算这它们的内积,使用了核技术之后,先把ϕ(x1)和ϕ(x2)的一般表达式<ϕ(x1),ϕ(x2)>=k(<ϕ(x1),ϕ(x2) >)计算出来,这里的<·,·>表示内积,k(·,·)就是对应的核函数,这个表达式往往非常简单,所以计算非常方便。

    这种将内积替换成核函数的方式被称为核技巧(kernel trick)

    2、非线性数据处理

    已经知道了核技巧是什么,但是为什么要这样做呢?我们先举一个简单的例子,进行说明。假设二维平面x-y上存在若干点,其中点集A服从 {x,y|x2+y2=1},点集B服从{x,y|x2+y2=9},那么这些点在二维平面上的分布是这样的:
    在这里插入图片描述

    蓝色的是点集A,红色的是点集B,他们在xy平面上并不能线性可分,即用一条直线分割(虽然肉眼是可以识别的)。采用映射(x,y)->(x,y,x2+y2)后,在三维空间的点的分布为:
    在这里插入图片描述

    可见红色和蓝色的点被映射到了不同的平面,在更高维空间中是线性可分的(用一个平面去分割)。

    上述例子中的样本点的分布遵循圆的分布。继续推广到椭圆的一般样本形式:
    在这里插入图片描述

    上图的两类数据分布为两个椭圆的形状,这样的数据本身就是不可分的。不难发现,这两个半径不同的椭圆是加上了少量的噪音生成得到的。所以,一个理想的分界应该也是一个椭圆,而不是一个直线。如果用X1和X2来表示这个二维平面的两个坐标的话,我们知道这个分界椭圆可以写为:
    在这里插入图片描述

    这个方程就是高中学过的椭圆一般方程。注意上面的形式,如果我们构造另外一个五维的空间,其中五个坐标的值分别为:
    在这里插入图片描述

    那么,显然我们可以将这个分界的椭圆方程写成如下形式:
    在这里插入图片描述

    这个关于新的坐标Z1,Z2,Z3,Z4,Z5的方程,就是一个超平面方程,它的维度是5。也就是说,如果我们做一个映射 ϕ : 二维 → 五维,将 X1,X2按照上面的规则映射为 Z1,Z2,··· ,Z5,那么在新的空间中原来的数据将变成线性可分的,从而使用之前我们推导的线性分类算法就可以进行处理了。

    我们举个简单的计算例子,现在假设已知的映射函数为:
    在这里插入图片描述

    这个是一个从2维映射到5维的例子。如果没有使用核函数,根据上一小节的介绍,我们需要先结算映射后的结果,然后再进行内积运算。那么对于两个向量a1=(x1,x2)和a2=(y1,y2)有:
    在这里插入图片描述

    另外,如果我们不进行映射计算,直接运算下面的公式:
    在这里插入图片描述

    你会发现,这两个公式的计算结果是相同的。区别在于什么呢?

    • 一个是根据映射函数,映射到高维空间中,然后再根据内积的公式进行计算,计算量大;
    • 另一个则直接在原来的低维空间中进行计算,而不需要显式地写出映射后的结果,计算量小。

    其实,在这个例子中,核函数就是:
    在这里插入图片描述

    我们通过k(x1,x2)的低维运算得到了先映射再内积的高维运算的结果,这就是核函数的神奇之处,它有效减少了我们的计算量。在这个例子中,我们对一个2维空间做映射,选择的新的空间是原始空间的所以一阶和二阶的组合,得到了5维的新空间;如果原始空间是3维的,那么我们会得到19维的新空间,这个数目是呈爆炸性增长的。如果我们使用ϕ(·)做映射计算,难度非常大,而且如果遇到无穷维的情况,就根本无从计算了。所以使用核函数进行计算是非常有必要的。

    3、核技巧的实现

    通过核技巧的转变,我们的分类函数变为:
    在这里插入图片描述

    我们的对偶问题变成了:
    在这里插入图片描述

    这样,我们就避开了高纬度空间中的计算。当然,我们刚刚的例子是非常简单的,我们可以手动构造出来对应映射的核函数出来,如果对于任意一个映射,要构造出对应的核函数就很困难了。因此,通常,人们会从一些常用的核函数中进行选择,根据问题和数据的不同,选择不同的参数,得到不同的核函数。接下来,要介绍的就是一个非常流行的核函数,那就是径向基核函数

    径向基核函数是SVM中常用的一个核函数。径向基核函数采用向量作为自变量的函数,能够基于向量举例运算输出一个标量。径向基核函数的高斯版本的公式如下:
    在这里插入图片描述

    其中,σ是用户自定义的用于确定到达率(reach)或者说函数值跌落到0的速度参数。上述高斯核函数将数据从原始空间映射到无穷维空间。关于无穷维空间,我们不必太担心。高斯核函数只是一个常用的核函数,使用者并不需要确切地理解数据到底是如何表现的,而且使用高斯核函数还会得到一个理想的结果。如果σ选得很大的话,高次特征上的权重实际上衰减得非常快,所以实际上(数值上近似一下)相当于一个低维的子空间;反过来,如果σ选得很小,则可以将任意的数据映射为线性可分——当然,这并不一定是好事,因为随之而来的可能是非常严重的过拟合问题。不过,总的来说,通过调控参数σ,高斯核实际上具有相当高的灵活性,也是使用最广泛的核函数之一。

    二、编程实现非线性SVM

    接下来,我们将使用testSetRBF.txt和testSetRBF2.txt,前者作为训练集,后者作为测试集。数据集下载地址:数据集下载

    1、可视化数据集

    Github代码地址

    先编写程序简单的看下数据集:

    #!/usr/bin/env python
    # -*- encoding: utf-8 -*-
    '''
    @File    :   svm2.py
    @Time    :   2019/06/17 17:13:43
    @Author  :   xiao ming 
    @Version :   1.0
    @Contact :   xiaoming3526@gmail.com
    @Desc    :   SVM支持向量机处理非线性数据
    @github  :   https://github.com/aimi-cn/AILearners
    '''
    
    import matplotlib.pyplot as plt
    import numpy as np
    import random
    
    '''
    @description: 读取数据
    @param:fileName - 文件名
    @return: dataMat - 数据矩阵
            labelMat - 数据标签
    '''
    def loadDataSet(filename):
        dataMat = []; labelMat = []
        fr = open(filename)
        #逐行读取,滤除空格等
        for line in fr.readlines():
            lineArr = line.strip().split('\t')
            dataMat.append([float(lineArr[0]), float(lineArr[1])])
            labelMat.append(float(lineArr[2]))
        return dataMat,labelMat
    
    '''
    @description: 数据可视化
    @param:dataMat - 数据矩阵
            labelMat - 数据标签
    @return: None
    '''
    def showDataSet(dataMat,labelMat):
        #正样本
        data_plus = []     
        #负样本                             
        data_minus = []                                 
        for i in range(len(dataMat)):
            if labelMat[i] > 0:
                data_plus.append(dataMat[i])
            else:
                data_minus.append(dataMat[i])
        #转换为numpy矩阵
        data_plus_np = np.array(data_plus)              
        data_minus_np = np.array(data_minus)           
        #正样本散点图
        plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1])   
        #负样本散点图
        plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1]) 
        plt.show()
    
    if __name__ == "__main__":
        dataMat,labelMat = loadDataSet('C:/Users/Administrator/Desktop/blog/github/AILearners/data/ml/jqxxsz/6.SVM/testSetRBF.txt')
        showDataSet(dataMat,labelMat)   
    

    运行结果:
    在这里插入图片描述

    可见,数据明显是线性不可分的。下面我们根据公式,编写核函数,并增加初始化参数kTup用于存储核函数有关的信息,同时我们只要将之前的内积运算变成核函数的运算即可。最后编写testRbf()函数,用于测试。在svm2.py文件中,编写代码如下:代码地址

    if __name__ == "__main__":
        start = time.clock()
        testRbf()
        end = time.clock()
        t=end-start
        print("Runtime is:",t)
    

    运行结果如下图所示:
    在这里插入图片描述

    可以看到,训练集错误率为9%,测试集错误率都是18%,训练耗时0.97s 。可以尝试更换不同的K1参数以观察测试错误率、训练错误率、支持向量个数随k1的变化情况。你会发现K1过大,会出现过拟合的情况,即训练集错误率低,但是测试集错误率高。

    三、Sklearn构建SVM分类器

    你的老板要求:你写的那个手写识别程序非常好,但是它占用内存太大。顾客无法通过无线的方式下载我们的应用。
    所以:我们可以考虑使用支持向量机,保留支持向量就行(knn需要保留所有的向量),就可以获得非常好的效果。

    最早的那篇KNN识别手写数字的文章介绍了KNN的算法和数据集:数据集介绍

    使用的数据集还是kNN用到的数据集(testDigits和trainingDigits):下载地址

    首先,我们先使用自己用python写的代码进行训练。创建文件svm_demo02.py文件,编写代码如下:代码地址

    SMO算法实现部分跟上文是一样的,我们新创建了img2vector()、loadImages()、testDigits()函数,它们分别用于二进制图形转换、图片加载、训练SVM分类器。我们自己的SVM分类器是个二类分类器,所以在设置标签的时候,将9作为负类,其余的0-8作为正类,进行训练。这是一种’ovr’思想,即one vs rest,就是对一个类别和剩余所有的类别进行分类。如果想实现10个数字的识别,一个简单的方法是,训练出10个分类器。这里简单起见,只训练了一个用于分类9和其余所有数字的分类器,运行结果如下:
    在这里插入图片描述

    可以看到,虽然我们进行了所谓的"优化",但是训练仍然很耗时,迭代10次,花费了537.8s。因为我们没有多进程、没有设置自动的终止条件,总之一句话,需要优化的地方太多了。尽管如此,我们训练后得到的结果还是不错的,可以看到训练集错误率为0,测试集错误率也仅为0.01%。

    然后我们来看这个超级方便的方法–sklearn.svm.SVC

    官方文档:地址

    部分参数说明:SVC这个函数

    • C:惩罚项,float类型,可选参数,默认为1.0,C越大,即对分错样本的惩罚程度越大,因此在训练样本中准确率越高,但是泛化能力降低,也就是对测试数据的分类准确率降低。相反,减小C的话,容许训练样本中有一些误分类错误样本,泛化能力强。对于训练样本带有噪声的情况,一般采用后者,把训练样本集中错误分类的样本作为噪声。
    • kernel:核函数类型,str类型,默认为’rbf’。可选参数为:
      • ‘linear’:线性核函数
      • ‘poly’:多项式核函数
      • ‘rbf’:径像核函数/高斯核
      • ‘sigmod’:sigmod核函数
      • ‘precomputed’:核矩阵
      • precomputed表示自己提前计算好核函数矩阵,这时候算法内部就不再用核函数去计算核矩阵,而是直接用你给的核矩阵,核矩阵需要为n*n的。

    SVC很是强大,我们不用理解算法实现的具体细节,不用理解算法的优化方法。同时,它也满足我们的多分类需求。创建文件svm-svc.py文件,编写代码如下:代码地址

    #!/usr/bin/env python
    # -*- encoding: utf-8 -*-
    '''
    @File    :   svm-svc.py
    @Time    :   2019/06/17 21:20:49
    @Author  :   xiao ming 
    @Version :   1.0
    @Contact :   xiaoming3526@gmail.com
    @Desc    :   sklearn.svm.SVC实现手写体识别
    @github  :   https://github.com/aimi-cn/AILearners
    '''
    
    # here put the import lib
    import numpy as np
    import operator
    from os import listdir
    from sklearn.svm import SVC
    import time
    
    def img2vector(filename):
        """
        将32x32的二进制图像转换为1x1024向量。
        Parameters:
            filename - 文件名
        Returns:
            returnVect - 返回的二进制图像的1x1024向量
        """
        #创建1x1024零向量
        returnVect = np.zeros((1, 1024))
        #打开文件
        fr = open(filename)
        #按行读取
        for i in range(32):
            #读一行数据
            lineStr = fr.readline()
            #每一行的前32个元素依次添加到returnVect中
            for j in range(32):
                returnVect[0, 32*i+j] = int(lineStr[j])
        #返回转换后的1x1024向量
        return returnVect
     
    def handwritingClassTest():
        """
        手写数字分类测试
        Parameters:
            无
        Returns:
            无
        """
        #测试集的Labels
        hwLabels = []
        #返回trainingDigits目录下的文件名
        trainingFileList = listdir('C:/Users/Administrator/Desktop/blog/github/AILearners/data/ml/jqxxsz/2.KNN/trainingDigits')
        #返回文件夹下文件的个数
        m = len(trainingFileList)
        #初始化训练的Mat矩阵,测试集
        trainingMat = np.zeros((m, 1024))
        #从文件名中解析出训练集的类别
        for i in range(m):
            #获得文件的名字
            fileNameStr = trainingFileList[i]
            #获得分类的数字
            classNumber = int(fileNameStr.split('_')[0])
            #将获得的类别添加到hwLabels中
            hwLabels.append(classNumber)
            #将每一个文件的1x1024数据存储到trainingMat矩阵中
            trainingMat[i,:] = img2vector('C:/Users/Administrator/Desktop/blog/github/AILearners/data/ml/jqxxsz/2.KNN/trainingDigits/%s' % (fileNameStr))
        clf = SVC(C=200,kernel='rbf')
        clf.fit(trainingMat,hwLabels)
        #返回testDigits目录下的文件列表
        testFileList = listdir('C:/Users/Administrator/Desktop/blog/github/AILearners/data/ml/jqxxsz/2.KNN/testDigits')
        #错误检测计数
        errorCount = 0.0
        #测试数据的数量
        mTest = len(testFileList)
        #从文件中解析出测试集的类别并进行分类测试
        for i in range(mTest):
            #获得文件的名字
            fileNameStr = testFileList[i]
            #获得分类的数字
            classNumber = int(fileNameStr.split('_')[0])
            #获得测试集的1x1024向量,用于训练
            vectorUnderTest = img2vector('C:/Users/Administrator/Desktop/blog/github/AILearners/data/ml/jqxxsz/2.KNN/testDigits/%s' % (fileNameStr))
            #获得预测结果
            # classifierResult = classify0(vectorUnderTest, trainingMat, hwLabels, 3)
            classifierResult = clf.predict(vectorUnderTest)
            print("分类返回结果为%d\t真实结果为%d" % (classifierResult, classNumber)).decode('utf-8').encode('gb2312')
            if(classifierResult != classNumber):
                errorCount += 1.0
        print("总共错了%d个数据\n错误率为%f%%" % (errorCount, errorCount/mTest * 100)).decode('utf-8').encode('gb2312')
     
    if __name__ == '__main__':
        start = time.clock()
        handwritingClassTest()
        end = time.clock()
        t=end-start
        print("Runtime is:",t)
    

    代码和kNN的实现是差不多的,就是换了个分类器而已。运行结果如下:

    在这里插入图片描述
    可以看到,训练和测试的时间总共加起来才5.9s。而且,测试集的错误率仅为1.37%。试着改变SVC的参数,慢慢体会一下吧~

    总结一下:

    SVM的优缺点

    优点

    -可用于线性/非线性分类,也可以用于回归,泛化错误率低,也就是说具有良好的学习能力,且学到的结果具有很好的推广性。

    • 可以解决小样本情况下的机器学习问题,可以解决高维问题,可以避免神经网络结构选择和局部极小点问题。
    • SVM是最好的现成的分类器,现成是指不加修改可直接使用。并且能够得到较低的错误率,SVM可以对训练集之外的数据点做很好的分类决策。

    缺点

    • 对参数调节和和函数的选择敏感。

    AIMI-CN AI学习交流群【1015286623】 获取更多AI资料
    扫码加群:

    在这里插入图片描述

    分享技术,乐享生活:我们的公众号计算机视觉这件小事每周推送“AI”系列资讯类文章,欢迎您的关注!

    在这里插入图片描述

    展开全文
  • 上篇文章讲解的是线性SVM的推导过程以及简化版SMO算法的代码实现。本篇文章将讲解SMO算法的优化方法以及非线性SVM

    转载请注明作者和出处: http://blog.csdn.net/c406495762
    **机器学习知乎专栏:**https://zhuanlan.zhihu.com/ml-jack
    **CSDN博客专栏:**http://blog.csdn.net/column/details/16415.html
    **Github代码获取:**https://github.com/Jack-Cherish/Machine-Learning/
    Python版本: Python3.x
    运行平台: Windows
    IDE: Sublime text3

    一 前言

    上篇文章讲解的是线性SVM的推导过程以及简化版SMO算法的代码实现。本篇文章将讲解SMO算法的优化方法以及非线性SVM。

    本文出现的所有代码,均可在我的github上下载,欢迎Follow、Star:https://github.com/Jack-Cherish/Machine-Learning

    更多精彩内容,尽在微信公众号,欢迎您的关注:

    在这里插入图片描述


    二 SMO算法优化

    在几百个点组成的小规模数据集上,简化版SMO算法的运行是没有什么问题的,但是在更大的数据集上的运行速度就会变慢。简化版SMO算法的第二个α的选择是随机的,针对这一问题,我们可以使用启发式选择第二个α值,来达到优化效果。

    1 启发选择方式

    下面这两个公式想必已经不再陌生:

    在实现SMO算法的时候,先计算η,再更新a_j。为了加快第二个α_j乘子的迭代速度,需要让直线的斜率增大,对于α_j的更新公式,其中η值没有什么文章可做,于是只能令:

    因此,我们可以明确自己的优化方法了:

    • 最外层循环,首先在样本中选择违反KKT条件的一个乘子作为最外层循环,然后用"启发式选择"选择另外一个乘子并进行这两个乘子的优化
    • 在非边界乘子中寻找使得|E_i - E_j|最大的样本
    • 如果没有找到,则从整个样本中随机选择一个样本

    接下来,让我们看看完整版SMO算法如何实现。

    2 完整版SMO算法

    完整版Platt SMO算法是通过一个外循环来选择违反KKT条件的一个乘子,并且其选择过程会在这两种方式之间进行交替:

    • 在所有数据集上进行单遍扫描
    • 在非边界α中实现单遍扫描

    非边界α指的就是那些不等于边界0或C的α值,并且跳过那些已知的不会改变的α值。所以我们要先建立这些α的列表,用于才能出α的更新状态。

    在选择第一个α值后,算法会通过"启发选择方式"选择第二个α值。

    3 编写代码

    我们首先构建一个仅包含init方法的optStruct类,将其作为一个数据结构来使用,方便我们对于重要数据的维护。代码思路和之前的简化版SMO算法是相似的,不同之处在于增加了优化方法,如果上篇文章已经看懂,我想这个代码会很好理解。创建一个svm-smo.py文件,编写代码如下:

    # -*-coding:utf-8 -*-
    import matplotlib.pyplot as plt
    import numpy as np
    import random
    
    """
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-10-03
    """
    
    class optStruct:
        """
        数据结构,维护所有需要操作的值
        Parameters:
            dataMatIn - 数据矩阵
            classLabels - 数据标签
            C - 松弛变量
            toler - 容错率
        """
        def __init__(self, dataMatIn, classLabels, C, toler):
            self.X = dataMatIn                                #数据矩阵
            self.labelMat = classLabels                        #数据标签
            self.C = C                                         #松弛变量
            self.tol = toler                                 #容错率
            self.m = np.shape(dataMatIn)[0]                 #数据矩阵行数
            self.alphas = np.mat(np.zeros((self.m,1)))         #根据矩阵行数初始化alpha参数为0   
            self.b = 0                                         #初始化b参数为0
            self.eCache = np.mat(np.zeros((self.m,2)))         #根据矩阵行数初始化虎误差缓存,第一列为是否有效的标志位,第二列为实际的误差E的值。
    
    def loadDataSet(fileName):
        """
        读取数据
        Parameters:
            fileName - 文件名
        Returns:
            dataMat - 数据矩阵
            labelMat - 数据标签
        """
        dataMat = []; labelMat = []
        fr = open(fileName)
        for line in fr.readlines():                                     #逐行读取,滤除空格等
            lineArr = line.strip().split('\t')
            dataMat.append([float(lineArr[0]), float(lineArr[1])])      #添加数据
            labelMat.append(float(lineArr[2]))                          #添加标签
        return dataMat,labelMat
    
    def calcEk(oS, k):
        """
        计算误差
        Parameters:
            oS - 数据结构
            k - 标号为k的数据
        Returns:
            Ek - 标号为k的数据误差
        """
        fXk = float(np.multiply(oS.alphas,oS.labelMat).T*(oS.X*oS.X[k,:].T) + oS.b)
        Ek = fXk - float(oS.labelMat[k])
        return Ek
    
    def selectJrand(i, m):
        """
        函数说明:随机选择alpha_j的索引值
    
        Parameters:
            i - alpha_i的索引值
            m - alpha参数个数
        Returns:
            j - alpha_j的索引值
        """
        j = i                                 #选择一个不等于i的j
        while (j == i):
            j = int(random.uniform(0, m))
        return j
    
    def selectJ(i, oS, Ei):
        """
        内循环启发方式2
        Parameters:
            i - 标号为i的数据的索引值
            oS - 数据结构
            Ei - 标号为i的数据误差
        Returns:
            j, maxK - 标号为j或maxK的数据的索引值
            Ej - 标号为j的数据误差
        """
        maxK = -1; maxDeltaE = 0; Ej = 0                         #初始化
        oS.eCache[i] = [1,Ei]                                      #根据Ei更新误差缓存
        validEcacheList = np.nonzero(oS.eCache[:,0].A)[0]        #返回误差不为0的数据的索引值
        if (len(validEcacheList)) > 1:                            #有不为0的误差
            for k in validEcacheList:                           #遍历,找到最大的Ek
                if k == i: continue                             #不计算i,浪费时间
                Ek = calcEk(oS, k)                                #计算Ek
                deltaE = abs(Ei - Ek)                            #计算|Ei-Ek|
                if (deltaE > maxDeltaE):                        #找到maxDeltaE
                    maxK = k; maxDeltaE = deltaE; Ej = Ek
            return maxK, Ej                                        #返回maxK,Ej
        else:                                                   #没有不为0的误差
            j = selectJrand(i, oS.m)                            #随机选择alpha_j的索引值
            Ej = calcEk(oS, j)                                    #计算Ej
        return j, Ej                                             #j,Ej
    
    def updateEk(oS, k):
        """
        计算Ek,并更新误差缓存
        Parameters:
            oS - 数据结构
            k - 标号为k的数据的索引值
        Returns:
            无
        """
        Ek = calcEk(oS, k)                                        #计算Ek
        oS.eCache[k] = [1,Ek]                                    #更新误差缓存
    
    
    def clipAlpha(aj,H,L):
        """
        修剪alpha_j
        Parameters:
            aj - alpha_j的值
            H - alpha上限
            L - alpha下限
        Returns:
            aj - 修剪后的alpah_j的值
        """
        if aj > H:
            aj = H
        if L > aj:
            aj = L
        return aj
    
    def innerL(i, oS):
        """
        优化的SMO算法
        Parameters:
            i - 标号为i的数据的索引值
            oS - 数据结构
        Returns:
            1 - 有任意一对alpha值发生变化
            0 - 没有任意一对alpha值发生变化或变化太小
        """
        #步骤1:计算误差Ei
        Ei = calcEk(oS, i)
        #优化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)):
            #使用内循环启发方式2选择alpha_j,并计算Ej
            j,Ej = selectJ(i, oS, Ei)
            #保存更新前的aplpha值,使用深拷贝
            alphaIold = oS.alphas[i].copy(); alphaJold = oS.alphas[j].copy();
            #步骤2:计算上下界L和H
            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
            #步骤3:计算eta
            eta = 2.0 * oS.X[i,:] * oS.X[j,:].T - oS.X[i,:] * oS.X[i,:].T - oS.X[j,:] * oS.X[j,:].T
            if eta >= 0:
                print("eta>=0")
                return 0
            #步骤4:更新alpha_j
            oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej)/eta
            #步骤5:修剪alpha_j
            oS.alphas[j] = clipAlpha(oS.alphas[j],H,L)
            #更新Ej至误差缓存
            updateEk(oS, j)
            if (abs(oS.alphas[j] - alphaJold) < 0.00001):
                print("alpha_j变化太小")
                return 0
            #步骤6:更新alpha_i
            oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])
            #更新Ei至误差缓存
            updateEk(oS, i)
            #步骤7:更新b_1和b_2
            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,:].T
            b2 = 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
            #步骤8:根据b_1和b_2更新b
            if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1
            elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2
            else: oS.b = (b1 + b2)/2.0
            return 1
        else:
            return 0
    
    def smoP(dataMatIn, classLabels, C, toler, maxIter):
        """
        完整的线性SMO算法
        Parameters:
            dataMatIn - 数据矩阵
            classLabels - 数据标签
            C - 松弛变量
            toler - 容错率
            maxIter - 最大迭代次数
        Returns:
            oS.b - SMO算法计算的b
            oS.alphas - SMO算法计算的alphas
        """
        oS = optStruct(np.mat(dataMatIn), np.mat(classLabels).transpose(), C, toler)                    #初始化数据结构
        iter = 0                                                                                         #初始化当前迭代次数
        entireSet = True; alphaPairsChanged = 0
        while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):                            #遍历整个数据集都alpha也没有更新或者超过最大迭代次数,则退出循环
            alphaPairsChanged = 0
            if entireSet:                                                                                #遍历整个数据集                           
                for i in range(oS.m):       
                    alphaPairsChanged += innerL(i,oS)                                                    #使用优化的SMO算法
                    print("全样本遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter,i,alphaPairsChanged))
                iter += 1
            else:                                                                                         #遍历非边界值
                nonBoundIs = np.nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]                        #遍历不在边界0和C的alpha
                for i in nonBoundIs:
                    alphaPairsChanged += innerL(i,oS)
                    print("非边界遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter,i,alphaPairsChanged))
                iter += 1
            if entireSet:                                                                                #遍历一次后改为非边界遍历
                entireSet = False
            elif (alphaPairsChanged == 0):                                                                #如果alpha没有更新,计算全样本遍历
                entireSet = True 
            print("迭代次数: %d" % iter)
        return oS.b,oS.alphas                                                                             #返回SMO算法计算的b和alphas
    
    
    def showClassifer(dataMat, classLabels, w, b):
        """
        分类结果可视化
        Parameters:
            dataMat - 数据矩阵
            w - 直线法向量
            b - 直线解决
        Returns:
            无
        """
        #绘制样本点
        data_plus = []                                  #正样本
        data_minus = []                                 #负样本
        for i in range(len(dataMat)):
            if classLabels[i] > 0:
                data_plus.append(dataMat[i])
            else:
                data_minus.append(dataMat[i])
        data_plus_np = np.array(data_plus)              #转换为numpy矩阵
        data_minus_np = np.array(data_minus)            #转换为numpy矩阵
        plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1], s=30, alpha=0.7)   #正样本散点图
        plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1], s=30, alpha=0.7) #负样本散点图
        #绘制直线
        x1 = max(dataMat)[0]
        x2 = min(dataMat)[0]
        a1, a2 = w
        b = float(b)
        a1 = float(a1[0])
        a2 = float(a2[0])
        y1, y2 = (-b- a1*x1)/a2, (-b - a1*x2)/a2
        plt.plot([x1, x2], [y1, y2])
        #找出支持向量点
        for i, alpha in enumerate(alphas):
            if abs(alpha) > 0:
                x, y = dataMat[i]
                plt.scatter([x], [y], s=150, c='none', alpha=0.7, linewidth=1.5, edgecolor='red')
        plt.show()
    
    
    def calcWs(alphas,dataArr,classLabels):
        """
        计算w
        Parameters:
            dataArr - 数据矩阵
            classLabels - 数据标签
            alphas - alphas值
        Returns:
            w - 计算得到的w
        """
        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)
        return w
    
    if __name__ == '__main__':
        dataArr, classLabels = loadDataSet('testSet.txt')
        b, alphas = smoP(dataArr, classLabels, 0.6, 0.001, 40)
        w = calcWs(alphas,dataArr, classLabels)
        showClassifer(dataArr, classLabels, w, b)
    

    完整版SMO算法(左图)与简化版SMO算法(右图)运行结果对比如下图所示:

    图中画红圈的样本点为支持向量上的点,是满足算法的一种解。完整版SMO算法覆盖整个数据集进行计算,而简化版SMO算法是随机选择的。可以看出,完整版SMO算法选出的支持向量样点更多,更接近理想的分隔超平面。

    对比两种算法的运算时间,我的测试结果是完整版SMO算法的速度比简化版SMO算法的速度快6倍左右

    其实,优化方法不仅仅是简单的启发式选择,还有其他优化方法,SMO算法速度还可以进一步提高。但是鉴于文章进度,这里不再进行展开。感兴趣的朋友,可以移步这里进行理论学习:http://www.cnblogs.com/zangrunqiang/p/5515872.html


    三 非线性SVM

    1 核技巧

    我们已经了解到,SVM如何处理线性可分的情况,而对于非线性的情况,SVM的处理方式就是选择一个核函数。简而言之:在线性不可分的情况下,SVM通过某种事先选择的非线性映射(核函数)将输入变量映到一个高维特征空间,将其变成在高维空间线性可分,在这个高维空间中构造最优分类超平面。

    根据上篇文章,线性可分的情况下,可知最终的超平面方程为:

    将上述公式用内积来表示:

    对于线性不可分,我们使用一个非线性映射,将数据映射到特征空间,在特征空间中使用线性学习器,分类函数变形如下:

    其中ϕ从输入空间(X)到某个特征空间(F)的映射,这意味着建立非线性学习器分为两步:

    • 首先使用一个非线性映射将数据变换到一个特征空间F;
    • 然后在特征空间使用线性学习器分类。

    如果有一种方法可以在特征空间中直接计算内积<ϕ(x_i),ϕ(x)>,就像在原始输入点的函数中一样,就有可能将两个步骤融合到一起建立一个分线性的学习器,这样直接计算的方法称为核函数方法。

    这里直接给出一个定义:核是一个函数k,对所有x,z∈X,满足k(x,z)=<ϕ(x_i),ϕ(x)>,这里ϕ(·)是从原始输入空间X到内积空间F的映射。

    简而言之:如果不是用核技术,就会先计算线性映ϕ(x_1)和ϕ(x_2),然后计算这它们的内积,使用了核技术之后,先把ϕ(x_1)和ϕ(x_2)的一般表达式<ϕ(x_1),ϕ(x_2)>=k(<ϕ(x_1),ϕ(x_2) >)计算出来,这里的<·,·>表示内积,k(·,·)就是对应的核函数,这个表达式往往非常简单,所以计算非常方便。

    这种将内积替换成核函数的方式被称为核技巧(kernel trick)。

    2 非线性数据处理

    已经知道了核技巧是什么,但是为什么要这样做呢?我们先举一个简单的例子,进行说明。假设二维平面x-y上存在若干点,其中点集A服从{x,y|x2+y2=1},点集B服从{x,y|x2+y2=9},那么这些点在二维平面上的分布是这样的:

    蓝色的是点集A,红色的是点集B,他们在xy平面上并不能线性可分,即用一条直线分割( 虽然肉眼是可以识别的) 。采用映射(x,y)->(x,y,x2+y2)后,在三维空间的点的分布为:

    可见红色和蓝色的点被映射到了不同的平面,在更高维空间中是线性可分的(用一个平面去分割)。

    上述例子中的样本点的分布遵循圆的分布。继续推广到椭圆的一般样本形式:

    上图的两类数据分布为两个椭圆的形状,这样的数据本身就是不可分的。不难发现,这两个半径不同的椭圆是加上了少量的噪音生成得到的。所以,一个理想的分界应该也是一个椭圆,而不是一个直线。如果用X1和X2来表示这个二维平面的两个坐标的话,我们知道这个分界椭圆可以写为:

    这个方程就是高中学过的椭圆一般方程。注意上面的形式,如果我们构造另外一个五维的空间,其中五个坐标的值分别为:

    那么,显然我们可以将这个分界的椭圆方程写成如下形式:

    这个关于新的坐标Z1,Z2,Z3,Z4,Z5的方程,就是一个超平面方程,它的维度是5。也就是说,如果我们做一个映射 ϕ : 二维 → 五维,将 X1,X2按照上面的规则映射为 Z1,Z2,··· ,Z5,那么在新的空间中原来的数据将变成线性可分的,从而使用之前我们推导的线性分类算法就可以进行处理了。

    我们举个简单的计算例子,现在假设已知的映射函数为:

    这个是一个从2维映射到5维的例子。如果没有使用核函数,根据上一小节的介绍,我们需要先结算映射后的结果,然后再进行内积运算。那么对于两个向量a1=(x1,x2)和a2=(y1,y2)有:

    另外,如果我们不进行映射计算,直接运算下面的公式:

    你会发现,这两个公式的计算结果是相同的。区别在于什么呢?

    • 一个是根据映射函数,映射到高维空间中,然后再根据内积的公式进行计算,计算量大;
    • 另一个则直接在原来的低维空间中进行计算,而不需要显式地写出映射后的结果,计算量小。

    其实,在这个例子中,核函数就是:

    我们通过k(x1,x2)的低维运算得到了先映射再内积的高维运算的结果,这就是核函数的神奇之处,它有效减少了我们的计算量。在这个例子中,我们对一个2维空间做映射,选择的新的空间是原始空间的所以一阶和二阶的组合,得到了5维的新空间;如果原始空间是3维的,那么我们会得到19维的新空间,这个数目是呈爆炸性增长的。如果我们使用ϕ(·)做映射计算,难度非常大,而且如果遇到无穷维的情况,就根本无从计算了。所以使用核函数进行计算是非常有必要的。

    3 核技巧的实现

    通过核技巧的转变,我们的分类函数变为:

    我们的对偶问题变成了:

    这样,我们就避开了高纬度空间中的计算。当然,我们刚刚的例子是非常简单的,我们可以手动构造出来对应映射的核函数出来,如果对于任意一个映射,要构造出对应的核函数就很困难了。因此,通常,人们会从一些常用的核函数中进行选择,根据问题和数据的不同,选择不同的参数,得到不同的核函数。接下来,要介绍的就是一个非常流行的核函数,那就是径向基核函数。

    径向基核函数是SVM中常用的一个核函数。径向基核函数采用向量作为自变量的函数,能够基于向量举例运算输出一个标量。径向基核函数的高斯版本的公式如下:

    其中,σ是用户自定义的用于确定到达率(reach)或者说函数值跌落到0的速度参数。上述高斯核函数将数据从原始空间映射到无穷维空间。关于无穷维空间,我们不必太担心。高斯核函数只是一个常用的核函数,使用者并不需要确切地理解数据到底是如何表现的,而且使用高斯核函数还会得到一个理想的结果。如果σ选得很大的话,高次特征上的权重实际上衰减得非常快,所以实际上(数值上近似一下)相当于一个低维的子空间;反过来,如果σ选得很小,则可以将任意的数据映射为线性可分——当然,这并不一定是好事,因为随之而来的可能是非常严重的过拟合问题。不过,总的来说,通过调控参数σ,高斯核实际上具有相当高的灵活性,也是使用最广泛的核函数之一。


    四 编程实现非线性SVM

    接下来,我们将使用testSetRBF.txt和testSetRBF2.txt,前者作为训练集,后者作为测试集。数据集下载地址:https://github.com/Jack-Cherish/Machine-Learning/tree/master/SVM

    1 可视化数据集

    我们先编写程序简单看下数据集:

    # -*-coding:utf-8 -*-
    import matplotlib.pyplot as plt
    import numpy as np
    
    def showDataSet(dataMat, labelMat):
        """
        数据可视化
        Parameters:
            dataMat - 数据矩阵
            labelMat - 数据标签
        Returns:
            无
        """
        data_plus = []                                  #正样本
        data_minus = []                                 #负样本
        for i in range(len(dataMat)):
            if labelMat[i] > 0:
                data_plus.append(dataMat[i])
            else:
                data_minus.append(dataMat[i])
        data_plus_np = np.array(data_plus)              #转换为numpy矩阵
        data_minus_np = np.array(data_minus)            #转换为numpy矩阵
        plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1])   #正样本散点图
        plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1]) #负样本散点图
        plt.show()
    
    if __name__ == '__main__':
        dataArr,labelArr = loadDataSet('testSetRBF.txt')                        #加载训练集
        showDataSet(dataArr, labelArr)
    

    程序运行结果:

    可见,数据明显是线性不可分的。下面我们根据公式,编写核函数,并增加初始化参数kTup用于存储核函数有关的信息,同时我们只要将之前的内积运算变成核函数的运算即可。最后编写testRbf()函数,用于测试。创建svmMLiA.py文件,编写代码如下:

    # -*-coding:utf-8 -*-
    import matplotlib.pyplot as plt
    import numpy as np
    import random
    
    """
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-10-03
    """
    
    class optStruct:
        """
        数据结构,维护所有需要操作的值
        Parameters:
            dataMatIn - 数据矩阵
            classLabels - 数据标签
            C - 松弛变量
            toler - 容错率
            kTup - 包含核函数信息的元组,第一个参数存放核函数类别,第二个参数存放必要的核函数需要用到的参数
        """
        def __init__(self, dataMatIn, classLabels, C, toler, kTup):
            self.X = dataMatIn                                #数据矩阵
            self.labelMat = classLabels                        #数据标签
            self.C = C                                         #松弛变量
            self.tol = toler                                 #容错率
            self.m = np.shape(dataMatIn)[0]                 #数据矩阵行数
            self.alphas = np.mat(np.zeros((self.m,1)))         #根据矩阵行数初始化alpha参数为0   
            self.b = 0                                         #初始化b参数为0
            self.eCache = np.mat(np.zeros((self.m,2)))         #根据矩阵行数初始化虎误差缓存,第一列为是否有效的标志位,第二列为实际的误差E的值。
            self.K = np.mat(np.zeros((self.m,self.m)))        #初始化核K
            for i in range(self.m):                            #计算所有数据的核K
                self.K[:,i] = kernelTrans(self.X, self.X[i,:], kTup)
    
    def kernelTrans(X, A, kTup):
        """
        通过核函数将数据转换更高维的空间
        Parameters:
            X - 数据矩阵
            A - 单个数据的向量
            kTup - 包含核函数信息的元组
        Returns:
            K - 计算的核K
        """
        m,n = np.shape(X)
        K = np.mat(np.zeros((m,1)))
        if kTup[0] == 'lin': K = X * A.T                       #线性核函数,只进行内积。
        elif kTup[0] == 'rbf':                                 #高斯核函数,根据高斯核函数公式进行计算
            for j in range(m):
                deltaRow = X[j,:] - A
                K[j] = deltaRow*deltaRow.T
            K = np.exp(K/(-1*kTup[1]**2))                     #计算高斯核K
        else: raise NameError('核函数无法识别')
        return K                                             #返回计算的核K
    
    def loadDataSet(fileName):
        """
        读取数据
        Parameters:
            fileName - 文件名
        Returns:
            dataMat - 数据矩阵
            labelMat - 数据标签
        """
        dataMat = []; labelMat = []
        fr = open(fileName)
        for line in fr.readlines():                                     #逐行读取,滤除空格等
            lineArr = line.strip().split('\t')
            dataMat.append([float(lineArr[0]), float(lineArr[1])])      #添加数据
            labelMat.append(float(lineArr[2]))                          #添加标签
        return dataMat,labelMat
    
    def calcEk(oS, k):
        """
        计算误差
        Parameters:
            oS - 数据结构
            k - 标号为k的数据
        Returns:
            Ek - 标号为k的数据误差
        """
        fXk = float(np.multiply(oS.alphas,oS.labelMat).T*oS.K[:,k] + oS.b)
        Ek = fXk - float(oS.labelMat[k])
        return Ek
    
    def selectJrand(i, m):
        """
        函数说明:随机选择alpha_j的索引值
    
        Parameters:
            i - alpha_i的索引值
            m - alpha参数个数
        Returns:
            j - alpha_j的索引值
        """
        j = i                                 #选择一个不等于i的j
        while (j == i):
            j = int(random.uniform(0, m))
        return j
    
    def selectJ(i, oS, Ei):
        """
        内循环启发方式2
        Parameters:
            i - 标号为i的数据的索引值
            oS - 数据结构
            Ei - 标号为i的数据误差
        Returns:
            j, maxK - 标号为j或maxK的数据的索引值
            Ej - 标号为j的数据误差
        """
        maxK = -1; maxDeltaE = 0; Ej = 0                         #初始化
        oS.eCache[i] = [1,Ei]                                      #根据Ei更新误差缓存
        validEcacheList = np.nonzero(oS.eCache[:,0].A)[0]        #返回误差不为0的数据的索引值
        if (len(validEcacheList)) > 1:                            #有不为0的误差
            for k in validEcacheList:                           #遍历,找到最大的Ek
                if k == i: continue                             #不计算i,浪费时间
                Ek = calcEk(oS, k)                                #计算Ek
                deltaE = abs(Ei - Ek)                            #计算|Ei-Ek|
                if (deltaE > maxDeltaE):                        #找到maxDeltaE
                    maxK = k; maxDeltaE = deltaE; Ej = Ek
            return maxK, Ej                                        #返回maxK,Ej
        else:                                                   #没有不为0的误差
            j = selectJrand(i, oS.m)                            #随机选择alpha_j的索引值
            Ej = calcEk(oS, j)                                    #计算Ej
        return j, Ej                                             #j,Ej
    
    def updateEk(oS, k):
        """
        计算Ek,并更新误差缓存
        Parameters:
            oS - 数据结构
            k - 标号为k的数据的索引值
        Returns:
            无
        """
        Ek = calcEk(oS, k)                                        #计算Ek
        oS.eCache[k] = [1,Ek]                                    #更新误差缓存
    
    def clipAlpha(aj,H,L):
        """
        修剪alpha_j
        Parameters:
            aj - alpha_j的值
            H - alpha上限
            L - alpha下限
        Returns:
            aj - 修剪后的alpah_j的值
        """
        if aj > H:
            aj = H
        if L > aj:
            aj = L
        return aj
    
    def innerL(i, oS):
        """
        优化的SMO算法
        Parameters:
            i - 标号为i的数据的索引值
            oS - 数据结构
        Returns:
            1 - 有任意一对alpha值发生变化
            0 - 没有任意一对alpha值发生变化或变化太小
        """
        #步骤1:计算误差Ei
        Ei = calcEk(oS, i)
        #优化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)):
            #使用内循环启发方式2选择alpha_j,并计算Ej
            j,Ej = selectJ(i, oS, Ei)
            #保存更新前的aplpha值,使用深拷贝
            alphaIold = oS.alphas[i].copy(); alphaJold = oS.alphas[j].copy();
            #步骤2:计算上下界L和H
            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
            #步骤3:计算eta
            eta = 2.0 * oS.K[i,j] - oS.K[i,i] - oS.K[j,j]
            if eta >= 0:
                print("eta>=0")
                return 0
            #步骤4:更新alpha_j
            oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej)/eta
            #步骤5:修剪alpha_j
            oS.alphas[j] = clipAlpha(oS.alphas[j],H,L)
            #更新Ej至误差缓存
            updateEk(oS, j)
            if (abs(oS.alphas[j] - alphaJold) < 0.00001):
                print("alpha_j变化太小")
                return 0
            #步骤6:更新alpha_i
            oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])
            #更新Ei至误差缓存
            updateEk(oS, i)
            #步骤7:更新b_1和b_2
            b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,i] - oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[i,j]
            b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,j]- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[j,j]
            #步骤8:根据b_1和b_2更新b
            if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1
            elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2
            else: oS.b = (b1 + b2)/2.0
            return 1
        else:
            return 0
    
    def smoP(dataMatIn, classLabels, C, toler, maxIter, kTup = ('lin',0)):
        """
        完整的线性SMO算法
        Parameters:
            dataMatIn - 数据矩阵
            classLabels - 数据标签
            C - 松弛变量
            toler - 容错率
            maxIter - 最大迭代次数
            kTup - 包含核函数信息的元组
        Returns:
            oS.b - SMO算法计算的b
            oS.alphas - SMO算法计算的alphas
        """
        oS = optStruct(np.mat(dataMatIn), np.mat(classLabels).transpose(), C, toler, kTup)                #初始化数据结构
        iter = 0                                                                                         #初始化当前迭代次数
        entireSet = True; alphaPairsChanged = 0
        while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):                            #遍历整个数据集都alpha也没有更新或者超过最大迭代次数,则退出循环
            alphaPairsChanged = 0
            if entireSet:                                                                                #遍历整个数据集                           
                for i in range(oS.m):       
                    alphaPairsChanged += innerL(i,oS)                                                    #使用优化的SMO算法
                    print("全样本遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter,i,alphaPairsChanged))
                iter += 1
            else:                                                                                         #遍历非边界值
                nonBoundIs = np.nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]                        #遍历不在边界0和C的alpha
                for i in nonBoundIs:
                    alphaPairsChanged += innerL(i,oS)
                    print("非边界遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter,i,alphaPairsChanged))
                iter += 1
            if entireSet:                                                                                #遍历一次后改为非边界遍历
                entireSet = False
            elif (alphaPairsChanged == 0):                                                                #如果alpha没有更新,计算全样本遍历
                entireSet = True 
            print("迭代次数: %d" % iter)
        return oS.b,oS.alphas                                                                             #返回SMO算法计算的b和alphas
    
    def testRbf(k1 = 1.3):
        """
        测试函数
        Parameters:
            k1 - 使用高斯核函数的时候表示到达率
        Returns:
            无
        """
        dataArr,labelArr = loadDataSet('testSetRBF.txt')                        #加载训练集
        b,alphas = smoP(dataArr, labelArr, 200, 0.0001, 100, ('rbf', k1))        #根据训练集计算b和alphas
        datMat = np.mat(dataArr); labelMat = np.mat(labelArr).transpose()
        svInd = np.nonzero(alphas.A > 0)[0]                                        #获得支持向量
        sVs = datMat[svInd]                                                     
        labelSV = labelMat[svInd];
        print("支持向量个数:%d" % np.shape(sVs)[0])
        m,n = np.shape(datMat)
        errorCount = 0
        for i in range(m):
            kernelEval = kernelTrans(sVs,datMat[i,:],('rbf', k1))                #计算各个点的核
            predict = kernelEval.T * np.multiply(labelSV,alphas[svInd]) + b     #根据支持向量的点,计算超平面,返回预测结果
            if np.sign(predict) != np.sign(labelArr[i]): errorCount += 1        #返回数组中各元素的正负符号,用1和-1表示,并统计错误个数
        print("训练集错误率: %.2f%%" % ((float(errorCount)/m)*100))             #打印错误率
        dataArr,labelArr = loadDataSet('testSetRBF2.txt')                         #加载测试集
        errorCount = 0
        datMat = np.mat(dataArr); labelMat = np.mat(labelArr).transpose()         
        m,n = np.shape(datMat)
        for i in range(m):
            kernelEval = kernelTrans(sVs,datMat[i,:],('rbf', k1))                 #计算各个点的核           
            predict=kernelEval.T * np.multiply(labelSV,alphas[svInd]) + b         #根据支持向量的点,计算超平面,返回预测结果
            if np.sign(predict) != np.sign(labelArr[i]): errorCount += 1        #返回数组中各元素的正负符号,用1和-1表示,并统计错误个数
        print("测试集错误率: %.2f%%" % ((float(errorCount)/m)*100))             #打印错误率
    
    def showDataSet(dataMat, labelMat):
        """
        数据可视化
        Parameters:
            dataMat - 数据矩阵
            labelMat - 数据标签
        Returns:
            无
        """
        data_plus = []                                  #正样本
        data_minus = []                                 #负样本
        for i in range(len(dataMat)):
            if labelMat[i] > 0:
                data_plus.append(dataMat[i])
            else:
                data_minus.append(dataMat[i])
        data_plus_np = np.array(data_plus)              #转换为numpy矩阵
        data_minus_np = np.array(data_minus)            #转换为numpy矩阵
        plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1])   #正样本散点图
        plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1]) #负样本散点图
        plt.show()
    
    if __name__ == '__main__':
        testRbf()
    

    运行结果如下图所示:

    可以看到,训练集错误率为1%,测试集错误率都是4%,训练耗时1.7s。可以尝试更换不同的K1参数以观察测试错误率、训练错误率、支持向量个数随k1的变化情况。你会发现K1过大,会出现过拟合的情况,即训练集错误率低,但是测试集错误率高。


    五 klearn构建SVM分类器

    在第一篇文章中,我们使用了kNN进行手写数字识别。它的缺点是存储空间大,因为要保留所有的训练样本,如果你的老板让你节约这个内存空间,并达到相同的识别效果,甚至更好。那这个时候,我们就要可以使用SVM了,因为它只需要保留支持向量即可,而且能获得可比的效果。

    使用的数据集还是kNN用到的数据集(testDigits和trainingDigits):https://github.com/Jack-Cherish/Machine-Learning/tree/master/kNN/3.%E6%95%B0%E5%AD%97%E8%AF%86%E5%88%AB

    如果对这个数据集不了解的,可以先看看我的第一篇文章:

    CSDN:http://blog.csdn.net/c406495762/article/details/75172850
    知乎:https://zhuanlan.zhihu.com/p/28656126

    首先,我们先使用自己用python写的代码进行训练。创建文件svm-digits.py文件,编写代码如下:

    # -*-coding:utf-8 -*-
    import matplotlib.pyplot as plt
    import numpy as np
    import random
    
    """
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-10-03
    """
    
    class optStruct:
        """
        数据结构,维护所有需要操作的值
        Parameters:
            dataMatIn - 数据矩阵
            classLabels - 数据标签
            C - 松弛变量
            toler - 容错率
            kTup - 包含核函数信息的元组,第一个参数存放核函数类别,第二个参数存放必要的核函数需要用到的参数
        """
        def __init__(self, dataMatIn, classLabels, C, toler, kTup):
            self.X = dataMatIn                                #数据矩阵
            self.labelMat = classLabels                        #数据标签
            self.C = C                                         #松弛变量
            self.tol = toler                                 #容错率
            self.m = np.shape(dataMatIn)[0]                 #数据矩阵行数
            self.alphas = np.mat(np.zeros((self.m,1)))         #根据矩阵行数初始化alpha参数为0   
            self.b = 0                                         #初始化b参数为0
            self.eCache = np.mat(np.zeros((self.m,2)))         #根据矩阵行数初始化虎误差缓存,第一列为是否有效的标志位,第二列为实际的误差E的值。
            self.K = np.mat(np.zeros((self.m,self.m)))        #初始化核K
            for i in range(self.m):                            #计算所有数据的核K
                self.K[:,i] = kernelTrans(self.X, self.X[i,:], kTup)
    
    def kernelTrans(X, A, kTup):
        """
        通过核函数将数据转换更高维的空间
        Parameters:
            X - 数据矩阵
            A - 单个数据的向量
            kTup - 包含核函数信息的元组
        Returns:
            K - 计算的核K
        """
        m,n = np.shape(X)
        K = np.mat(np.zeros((m,1)))
        if kTup[0] == 'lin': K = X * A.T                       #线性核函数,只进行内积。
        elif kTup[0] == 'rbf':                                 #高斯核函数,根据高斯核函数公式进行计算
            for j in range(m):
                deltaRow = X[j,:] - A
                K[j] = deltaRow*deltaRow.T
            K = np.exp(K/(-1*kTup[1]**2))                     #计算高斯核K
        else: raise NameError('核函数无法识别')
        return K                                             #返回计算的核K
    
    def loadDataSet(fileName):
        """
        读取数据
        Parameters:
            fileName - 文件名
        Returns:
            dataMat - 数据矩阵
            labelMat - 数据标签
        """
        dataMat = []; labelMat = []
        fr = open(fileName)
        for line in fr.readlines():                                     #逐行读取,滤除空格等
            lineArr = line.strip().split('\t')
            dataMat.append([float(lineArr[0]), float(lineArr[1])])      #添加数据
            labelMat.append(float(lineArr[2]))                          #添加标签
        return dataMat,labelMat
    
    def calcEk(oS, k):
        """
        计算误差
        Parameters:
            oS - 数据结构
            k - 标号为k的数据
        Returns:
            Ek - 标号为k的数据误差
        """
        fXk = float(np.multiply(oS.alphas,oS.labelMat).T*oS.K[:,k] + oS.b)
        Ek = fXk - float(oS.labelMat[k])
        return Ek
    
    def selectJrand(i, m):
        """
        函数说明:随机选择alpha_j的索引值
    
        Parameters:
            i - alpha_i的索引值
            m - alpha参数个数
        Returns:
            j - alpha_j的索引值
        """
        j = i                                 #选择一个不等于i的j
        while (j == i):
            j = int(random.uniform(0, m))
        return j
    
    def selectJ(i, oS, Ei):
        """
        内循环启发方式2
        Parameters:
            i - 标号为i的数据的索引值
            oS - 数据结构
            Ei - 标号为i的数据误差
        Returns:
            j, maxK - 标号为j或maxK的数据的索引值
            Ej - 标号为j的数据误差
        """
        maxK = -1; maxDeltaE = 0; Ej = 0                         #初始化
        oS.eCache[i] = [1,Ei]                                      #根据Ei更新误差缓存
        validEcacheList = np.nonzero(oS.eCache[:,0].A)[0]        #返回误差不为0的数据的索引值
        if (len(validEcacheList)) > 1:                            #有不为0的误差
            for k in validEcacheList:                           #遍历,找到最大的Ek
                if k == i: continue                             #不计算i,浪费时间
                Ek = calcEk(oS, k)                                #计算Ek
                deltaE = abs(Ei - Ek)                            #计算|Ei-Ek|
                if (deltaE > maxDeltaE):                        #找到maxDeltaE
                    maxK = k; maxDeltaE = deltaE; Ej = Ek
            return maxK, Ej                                        #返回maxK,Ej
        else:                                                   #没有不为0的误差
            j = selectJrand(i, oS.m)                            #随机选择alpha_j的索引值
            Ej = calcEk(oS, j)                                    #计算Ej
        return j, Ej                                             #j,Ej
    
    def updateEk(oS, k):
        """
        计算Ek,并更新误差缓存
        Parameters:
            oS - 数据结构
            k - 标号为k的数据的索引值
        Returns:
            无
        """
        Ek = calcEk(oS, k)                                        #计算Ek
        oS.eCache[k] = [1,Ek]                                    #更新误差缓存
    
    
    def clipAlpha(aj,H,L):
        """
        修剪alpha_j
        Parameters:
            aj - alpha_j的值
            H - alpha上限
            L - alpha下限
        Returns:
            aj - 修剪后的alpah_j的值
        """
        if aj > H:
            aj = H
        if L > aj:
            aj = L
        return aj
    
    def innerL(i, oS):
        """
        优化的SMO算法
        Parameters:
            i - 标号为i的数据的索引值
            oS - 数据结构
        Returns:
            1 - 有任意一对alpha值发生变化
            0 - 没有任意一对alpha值发生变化或变化太小
        """
        #步骤1:计算误差Ei
        Ei = calcEk(oS, i)
        #优化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)):
            #使用内循环启发方式2选择alpha_j,并计算Ej
            j,Ej = selectJ(i, oS, Ei)
            #保存更新前的aplpha值,使用深拷贝
            alphaIold = oS.alphas[i].copy(); alphaJold = oS.alphas[j].copy();
            #步骤2:计算上下界L和H
            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
            #步骤3:计算eta
            eta = 2.0 * oS.K[i,j] - oS.K[i,i] - oS.K[j,j]
            if eta >= 0:
                print("eta>=0")
                return 0
            #步骤4:更新alpha_j
            oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej)/eta
            #步骤5:修剪alpha_j
            oS.alphas[j] = clipAlpha(oS.alphas[j],H,L)
            #更新Ej至误差缓存
            updateEk(oS, j)
            if (abs(oS.alphas[j] - alphaJold) < 0.00001):
                print("alpha_j变化太小")
                return 0
            #步骤6:更新alpha_i
            oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])
            #更新Ei至误差缓存
            updateEk(oS, i)
            #步骤7:更新b_1和b_2
            b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,i] - oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[i,j]
            b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,j]- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[j,j]
            #步骤8:根据b_1和b_2更新b
            if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1
            elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2
            else: oS.b = (b1 + b2)/2.0
            return 1
        else:
            return 0
    
    def smoP(dataMatIn, classLabels, C, toler, maxIter, kTup = ('lin',0)):
        """
        完整的线性SMO算法
        Parameters:
            dataMatIn - 数据矩阵
            classLabels - 数据标签
            C - 松弛变量
            toler - 容错率
            maxIter - 最大迭代次数
            kTup - 包含核函数信息的元组
        Returns:
            oS.b - SMO算法计算的b
            oS.alphas - SMO算法计算的alphas
        """
        oS = optStruct(np.mat(dataMatIn), np.mat(classLabels).transpose(), C, toler, kTup)                #初始化数据结构
        iter = 0                                                                                         #初始化当前迭代次数
        entireSet = True; alphaPairsChanged = 0
        while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):                            #遍历整个数据集都alpha也没有更新或者超过最大迭代次数,则退出循环
            alphaPairsChanged = 0
            if entireSet:                                                                                #遍历整个数据集                           
                for i in range(oS.m):       
                    alphaPairsChanged += innerL(i,oS)                                                    #使用优化的SMO算法
                    print("全样本遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter,i,alphaPairsChanged))
                iter += 1
            else:                                                                                         #遍历非边界值
                nonBoundIs = np.nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]                        #遍历不在边界0和C的alpha
                for i in nonBoundIs:
                    alphaPairsChanged += innerL(i,oS)
                    print("非边界遍历:第%d次迭代 样本:%d, alpha优化次数:%d" % (iter,i,alphaPairsChanged))
                iter += 1
            if entireSet:                                                                                #遍历一次后改为非边界遍历
                entireSet = False
            elif (alphaPairsChanged == 0):                                                                #如果alpha没有更新,计算全样本遍历
                entireSet = True 
            print("迭代次数: %d" % iter)
        return oS.b,oS.alphas                                                                             #返回SMO算法计算的b和alphas
    
    
    def img2vector(filename):
        """
        将32x32的二进制图像转换为1x1024向量。
        Parameters:
            filename - 文件名
        Returns:
            returnVect - 返回的二进制图像的1x1024向量
        """
        returnVect = np.zeros((1,1024))
        fr = open(filename)
        for i in range(32):
            lineStr = fr.readline()
            for j in range(32):
                returnVect[0,32*i+j] = int(lineStr[j])
        return returnVect
    
    def loadImages(dirName):
        """
        加载图片
        Parameters:
            dirName - 文件夹的名字
        Returns:
            trainingMat - 数据矩阵
            hwLabels - 数据标签
        """
        from os import listdir
        hwLabels = []
        trainingFileList = listdir(dirName)           
        m = len(trainingFileList)
        trainingMat = np.zeros((m,1024))
        for i in range(m):
            fileNameStr = trainingFileList[i]
            fileStr = fileNameStr.split('.')[0]     
            classNumStr = int(fileStr.split('_')[0])
            if classNumStr == 9: hwLabels.append(-1)
            else: hwLabels.append(1)
            trainingMat[i,:] = img2vector('%s/%s' % (dirName, fileNameStr))
        return trainingMat, hwLabels   
    
    def testDigits(kTup=('rbf', 10)):
        """
        测试函数
        Parameters:
            kTup - 包含核函数信息的元组
        Returns:
            无
        """
        dataArr,labelArr = loadImages('trainingDigits')
        b,alphas = smoP(dataArr, labelArr, 200, 0.0001, 10, kTup)
        datMat = np.mat(dataArr); labelMat = np.mat(labelArr).transpose()
        svInd = np.nonzero(alphas.A>0)[0]
        sVs=datMat[svInd]
        labelSV = labelMat[svInd];
        print("支持向量个数:%d" % np.shape(sVs)[0])
        m,n = np.shape(datMat)
        errorCount = 0
        for i in range(m):
            kernelEval = kernelTrans(sVs,datMat[i,:],kTup)
            predict=kernelEval.T * np.multiply(labelSV,alphas[svInd]) + b
            if np.sign(predict) != np.sign(labelArr[i]): errorCount += 1
        print("训练集错误率: %.2f%%" % (float(errorCount)/m))
        dataArr,labelArr = loadImages('testDigits')
        errorCount = 0
        datMat = np.mat(dataArr); labelMat = np.mat(labelArr).transpose()
        m,n = np.shape(datMat)
        for i in range(m):
            kernelEval = kernelTrans(sVs,datMat[i,:],kTup)
            predict=kernelEval.T * np.multiply(labelSV,alphas[svInd]) + b
            if np.sign(predict) != np.sign(labelArr[i]): errorCount += 1   
        print("测试集错误率: %.2f%%" % (float(errorCount)/m))
    
    if __name__ == '__main__':
        testDigits()
    

    SMO算法实现部分跟上文是一样的,我们新创建了img2vector()、loadImages()、testDigits()函数,它们分别用于二进制图形转换、图片加载、训练SVM分类器。我们自己的SVM分类器是个二类分类器,所以在设置标签的时候,将9作为负类,其余的0-8作为正类,进行训练。这是一种’ovr’思想,即one vs rest,就是对一个类别和剩余所有的类别进行分类。如果想实现10个数字的识别,一个简单的方法是,训练出10个分类器。这里简单起见,只训练了一个用于分类9和其余所有数字的分类器,运行结果如下:

    可以看到,虽然我们进行了所谓的**“优化”**,但是训练仍然很耗时,迭代10次,花费了307.4s。因为我们没有多进程、没有设置自动的终止条件,总之一句话,需要优化的地方太多了。尽管如此,我们训练后得到的结果还是不错的,可以看到训练集错误率为0,测试集错误率也仅为0.01%。

    接下来,就是讲解本文的重头戏:sklearn.svm.SVC。

    1 Sklearn.svm.SVC

    官方英文文档手册:http://scikit-learn.org/stable/modules/generated/sklearn.svm.SVC.html

    sklearn.svm模块提供了很多模型供我们使用,本文使用的是svm.SVC,它是基于libsvm实现的。

    让我们先看下SVC这个函数,一共有14个参数:

    参数说明如下:

    • C:惩罚项,float类型,可选参数,默认为1.0,C越大,即对分错样本的惩罚程度越大,因此在训练样本中准确率越高,但是泛化能力降低,也就是对测试数据的分类准确率降低。相反,减小C的话,容许训练样本中有一些误分类错误样本,泛化能力强。对于训练样本带有噪声的情况,一般采用后者,把训练样本集中错误分类的样本作为噪声。
    • kernel:核函数类型,str类型,默认为’rbf’。可选参数为:
      • ‘linear’:线性核函数
      • ‘poly’:多项式核函数
      • ‘rbf’:径像核函数/高斯核
      • ‘sigmod’:sigmod核函数
      • ‘precomputed’:核矩阵
      • precomputed表示自己提前计算好核函数矩阵,这时候算法内部就不再用核函数去计算核矩阵,而是直接用你给的核矩阵,核矩阵需要为n*n的。
    • degree:多项式核函数的阶数,int类型,可选参数,默认为3。这个参数只对多项式核函数有用,是指多项式核函数的阶数n,如果给的核函数参数是其他核函数,则会自动忽略该参数。
    • gamma:核函数系数,float类型,可选参数,默认为auto。只对’rbf’ ,‘poly’ ,'sigmod’有效。如果gamma为auto,代表其值为样本特征数的倒数,即1/n_features。
    • coef0:核函数中的独立项,float类型,可选参数,默认为0.0。只有对’poly’ 和,'sigmod’核函数有用,是指其中的参数c。
    • probability:是否启用概率估计,bool类型,可选参数,默认为False,这必须在调用fit()之前启用,并且会fit()方法速度变慢。
    • shrinking:是否采用启发式收缩方式,bool类型,可选参数,默认为True。
    • tol:svm停止训练的误差精度,float类型,可选参数,默认为1e^-3。
    • cache_size:内存大小,float类型,可选参数,默认为200。指定训练所需要的内存,以MB为单位,默认为200MB。
    • class_weight:类别权重,dict类型或str类型,可选参数,默认为None。给每个类别分别设置不同的惩罚参数C,如果没有给,则会给所有类别都给C=1,即前面参数指出的参数C。如果给定参数’balance’,则使用y的值自动调整与输入数据中的类频率成反比的权重。
    • verbose:是否启用详细输出,bool类型,默认为False,此设置利用libsvm中的每个进程运行时设置,如果启用,可能无法在多线程上下文中正常工作。一般情况都设为False,不用管它。
    • max_iter:最大迭代次数,int类型,默认为-1,表示不限制。
    • decision_function_shape:决策函数类型,可选参数’ovo’和’ovr’,默认为’ovr’。'ovo’表示one vs one,'ovr’表示one vs rest。
    • random_state:数据洗牌时的种子值,int类型,可选参数,默认为None。伪随机数发生器的种子,在混洗数据时用于概率估计。

    其实,只要自己写了SMO算法,每个参数的意思,大概都是能明白的。

    2 编写代码

    SVC很是强大,我们不用理解算法实现的具体细节,不用理解算法的优化方法。同时,它也满足我们的多分类需求。创建文件svm-svc.py文件,编写代码如下:

    # -*- coding: UTF-8 -*-
    import numpy as np
    import operator
    from os import listdir
    from sklearn.svm import SVC
    
    """
    Author:
        Jack Cui
    Blog:
        http://blog.csdn.net/c406495762
    Zhihu:
        https://www.zhihu.com/people/Jack--Cui/
    Modify:
        2017-10-04
    """
    
    def img2vector(filename):
        """
        将32x32的二进制图像转换为1x1024向量。
        Parameters:
            filename - 文件名
        Returns:
            returnVect - 返回的二进制图像的1x1024向量
        """
        #创建1x1024零向量
        returnVect = np.zeros((1, 1024))
        #打开文件
        fr = open(filename)
        #按行读取
        for i in range(32):
            #读一行数据
            lineStr = fr.readline()
            #每一行的前32个元素依次添加到returnVect中
            for j in range(32):
                returnVect[0, 32*i+j] = int(lineStr[j])
        #返回转换后的1x1024向量
        return returnVect
    
    def handwritingClassTest():
        """
        手写数字分类测试
        Parameters:
            无
        Returns:
            无
        """
        #测试集的Labels
        hwLabels = []
        #返回trainingDigits目录下的文件名
        trainingFileList = listdir('trainingDigits')
        #返回文件夹下文件的个数
        m = len(trainingFileList)
        #初始化训练的Mat矩阵,测试集
        trainingMat = np.zeros((m, 1024))
        #从文件名中解析出训练集的类别
        for i in range(m):
            #获得文件的名字
            fileNameStr = trainingFileList[i]
            #获得分类的数字
            classNumber = int(fileNameStr.split('_')[0])
            #将获得的类别添加到hwLabels中
            hwLabels.append(classNumber)
            #将每一个文件的1x1024数据存储到trainingMat矩阵中
            trainingMat[i,:] = img2vector('trainingDigits/%s' % (fileNameStr))
        clf = SVC(C=200,kernel='rbf')
        clf.fit(trainingMat,hwLabels)
        #返回testDigits目录下的文件列表
        testFileList = listdir('testDigits')
        #错误检测计数
        errorCount = 0.0
        #测试数据的数量
        mTest = len(testFileList)
        #从文件中解析出测试集的类别并进行分类测试
        for i in range(mTest):
            #获得文件的名字
            fileNameStr = testFileList[i]
            #获得分类的数字
            classNumber = int(fileNameStr.split('_')[0])
            #获得测试集的1x1024向量,用于训练
            vectorUnderTest = img2vector('testDigits/%s' % (fileNameStr))
            #获得预测结果
            # classifierResult = classify0(vectorUnderTest, trainingMat, hwLabels, 3)
            classifierResult = clf.predict(vectorUnderTest)
            print("分类返回结果为%d\t真实结果为%d" % (classifierResult, classNumber))
            if(classifierResult != classNumber):
                errorCount += 1.0
        print("总共错了%d个数据\n错误率为%f%%" % (errorCount, errorCount/mTest * 100))
    
    if __name__ == '__main__':
        handwritingClassTest()
    

    代码和kNN的实现是差不多的,就是换了个分类器而已。运行结果如下:


    六 总结

    1 SVM的优缺点

    优点

    • 可用于线性/非线性分类,也可以用于回归,泛化错误率低,也就是说具有良好的学习能力,且学到的结果具有很好的推广性。
    • 可以解决小样本情况下的机器学习问题,可以解决高维问题,可以避免神经网络结构选择和局部极小点问题。
    • SVM是最好的现成的分类器,现成是指不加修改可直接使用。并且能够得到较低的错误率,SVM可以对训练集之外的数据点做很好的分类决策。

    缺点

    • 对参数调节和和函数的选择敏感。

    2 其他

    • 至此,关于SVM的文章已经写完,还有一些理论和细节,可能会在今后的文章提及。
    • 下篇文章将讲解AdaBoost,欢迎各位的捧场!
    • 如有问题,请留言。如有错误,还望指正,谢谢!

    PS: 如果觉得本篇本章对您有所帮助,欢迎关注、评论、赞!

    参考资料:

    • [1] SVM多维空间线性可分的理解:https://www.zhihu.com/question/27210162/answer/44815488
    • [2] 《机器学习实战》第六章
    展开全文
  • 使用线性SVM和HOG功能对照片进行行人分类 HOG简介 HOG是“定向直方图”的首字母缩写。 这是一种称为特征描述符的算法,可帮助在计算机视觉和图像处理模型中进行对象检测。 HOG是一种特征描述符,可对图像局部区域中...
  • 线性SVM,非线性SVM计算量太大。 OpenCV支持HOG+SVM和HOG+Cascade,都采用滑动窗口技术。 2.2 HOG+Adaboost 由于HOG + SVM的方案计算量太大,为了提高速度,后面有研究者参考了VJ[6]在人脸检测中的分类器设计思路,...

    目录

    1. 行人检测算法研究综述

    2. 基于机器学习的方法:人工特征+分类器

    2.1 HOG+SVM

    2.2 HOG+Adaboost

    2.3 ICF+AdaBoost

    2.4 DPM+ latent SVM

    3. HOG+SVM环境配置

    3.1 数据集INRIADATA

    3.2 算法原理

    3.2.1 HOG梯度方向直方图 Histogram of Oriented Gradients

    3.2.2 使用线性SVM进行训练

    3.3.3 代码

    4. 环境配置:Matlab安装libsvm

    5. 行人检测方面的数据集

     1. HOG+SVM使用的行人识别数据集

    2. 红外行人数据集:KAIST Multispectral Pedestrain Detection Benchmark

    3. 南方科技大学的红外数据集:SCUT_FIR_Pedestrian_Dataset

    6. 在自己的数据集上分析和优化

    1. 如何学习怎么查找模型问题所在?

    1.1 SVM深入理解

    1.3 SVM中解决泛化问题

    2. Matlab svm训练的模型如何在C++中使用

    7. 分析和优化模型

    当前状态:

    2.其他的一些疑问 

    3. 红外数据集

    Dataset 01: OSU Thermal Pedestrian Database

    Dataset 03: OSU Color-Thermal Database

    8. 指标方面

     1. 图像指标对应 

    2. 正负样本的问题

    2.1 样本的数量

    2.2 负样本

    2.3 应该把脸部的局部部位作为负样本吗?如眼睛或鼻子

    2.4 该怎么收集负样本信息?

    2.5 在步骤4中。

    2.6 探测器的最后增强步骤

    2.7 real numbers

    3. 负样本的concerns

    3.1 不同环境下的负样本

    3.2 不同尺度(距离)下的负样本

    4、目标检测——降低误检测率:负样本制作、增强、训练

    为什么要训练负样本?

    如何收集负样本?

    怎么消除误检测?

    5、评价指标

    6、大规模训练结果



    1. 行人检测算法研究综述

    • 参考资料https://zhuanlan.zhihu.com/p/51438953
    • 要解决的问题:找出图像或视频帧中所有的行人,包括位置和大小,一般用矩形框表示。
    • 主要难题: 外观差异大,遮挡问题,背景复杂,检测速度。
    • 行人数据库: INRIA 数据库、Caltech 数据库和 TUD 行人数据库。

    目标检测现在主要分为以下几个方向:

    • (1)基于运动检测的算法——背景建模。只利用像素级的信息,没有利用图像中更高层的语义信息。利用背景建模方法,提取出前景运动的目标,在目标区域内进行特征提取,然后利用分类器进行分类,判断是否包含行人;
      • 只能检测运动目标,不处理静止目标
      • 受光照变化、阴影影响大
      • 目标伪装色造成漏检和断裂
      • 受恶劣天气,雨雪,以及树叶摇晃等干扰物的影响
      • 多目标粘连,重叠,则无法处理。
    • (2)基于机器学习的方法(现阶段行人检测算法的主流),也是目前行人检测最常用的方法,根据大量的样本构建行人检测分类器。提取的特征主要有目标的灰度、边缘、纹理、颜色、梯度直方图等信息。分类器主要包括神经网络、SVM、adaboost以及现在被计算机视觉视为宠儿的深度学习。
      • 人工特征(括颜色,边缘,纹理等)+分类器(有神经网络,线性SVM,AdaBoost,随机森林等计算机视觉领域常用的算法)

    统计机器学习目前存在的难点

        • 行人的姿态、服饰各不相同、复杂的背景、不同的行人尺度以及不同的光照环境。
        • 提取的特征在特征空间中的分布不够紧凑;
        • 分类器的性能受训练样本的影响较大;
        • 离线训练时的负样本无法涵盖所有真实应用场景的情况;

    目前的行人检测基本上都是基于法国研究人员Dalal在2005的CVPR发表的HOG+SVM的行人检测算法(Histograms of Oriented Gradients for Human Detection, Navneet Dalel,Bill Triggs, CVPR2005)。HOG+SVM作为经典算法也集成到opencv里面去了,可以直接调用实现行人检测。

    为了解决速度问题可以采用背景差分法的统计学习行人检测,前提是背景建模的方法足够有效(即效果好速度快),目前获得比较好的检测效果的方法通常采用多特征融合的方法以及级联分类器。(常用的特征有Harry-like、Hog特征、LBP特征、Edgelet特征、CSS特征、COV特征、积分通道特征以及CENTRIST特征)。

    https://blog.csdn.net/SIGAI_CSDN/article/details/80693322的内容提要如下:

    2. 基于机器学习的方法:人工特征+分类器

    2.1 HOG+SVM

    人体有自身的外观特征,我们可以手工设计出特征,然后用这种特征来训练分类器用于区分行人和背景。这些特征包括颜色,边缘,纹理等机器学习中常用的特征,采用的分类器有神经网络,SVM,AdaBoost,随机森林等计算机视觉领域常用的算法。由于是检测问题,因此一般采用滑动窗口的技术,这在SIGAI之前的公众号文章“人脸检测算法综述”,“基于深度学习的目标检测算法综述”中已经介绍过了。

    行人检测第一个有里程碑意义的成果是Navneet Dalal在2005的CVPR中提出的基于HOG + SVM的行人检测算法[5]。Navneet Dalal是行人检测中之前经常使用的INRIA数据集的缔造者。 

    梯度方向直方图(HOG)是一种边缘特征,它利用了边缘的朝向和强度信息,后来被广泛应用于车辆检测,车牌检测等视觉目标检测问题。HOG的做法是固定大小的图像先计算梯度然后进行网格划分计算每个点处的梯度朝向和强度,然后形成网格内的所有像素的梯度方向分布直方图,最后汇总起来,形成整个直方图特征

    这一特征很好的描述了行人的形状、外观信息,比Haar特征更为强大,另外,该特征对光照变化和小量的空间平移不敏感。下图为用HOG特征进行行人检测的流程:

    得到候选区域的HOG特征后,需要利用分类器对该区域进行分类,确定是行人还是背景区域。在实现时,使用了线性支持向量机,这是因为采用非线性核的支持向量机在预测时的计算量太大,与支持向量的个数成正比。如果读者对这一问题感兴趣,可以阅读SIGAI之前关于SVM的文章。

    目前OpenCV中的行人检测算法支持HOG+SVM以及HOG+Cascade两种,二者都采用了滑动窗口技术,用固定大小的窗口扫描整个图像,然后对每一个窗口进行前景和背景的二分类。为了检测不同大小的行人,还需要对图像进行缩放。

    • HOG计算梯度方向和强度->梯度直方图,比Haar特征更为强大,且对光照变化和小量空间平移不敏感。
    • 得到HOG特征后,使用分类器进行分类。线性SVM,非线性SVM计算量太大。
    • OpenCV支持HOG+SVM和HOG+Cascade,都采用滑动窗口技术。

    2.2 HOG+Adaboost

    由于HOG + SVM的方案计算量太大,为了提高速度,后面有研究者参考了VJ[6]在人脸检测中的分类器设计思路,AdaBoost分类器级联的策略应用到了人体检测中,只是将Haar特征替换成HOG特征,因为Haar特征过于简单,无法描述人体这种复杂形状的目标。下图为基于级联Cascade分类器的检测流程:

    图中每一级中的分类器都是利用AdaBoost算法学习到的一个强分类器,处于前面的几个强分类器由于在分类器训练的时候会优先选择弱分类器,可以把最好的几个弱分类器进行集成,所有只需要很少的几个就可以达到预期效果,计算会非常简单,速度很快,大部分背景窗口很快会被排除掉,剩下很少一部分候选区域或通过后续的几级分类器进行判别,最终整体的检测速度有了很大的提升,相同条件下的预测时间只有基于SVM方法的十分之一。

    2.3 ICF+AdaBoost

    HOG特征只关注了物体的边缘和形状信息,对目标的表观信息并没有有效记录,所以很难处理遮挡问题,而且由于梯度的性质,该特征对噪点敏感。针对这些问题后面有人提出了积分通道特征(ICF[7],积分通道特征包括10个通道:

    6 个方向的梯度直方图,3 个LUV 颜色通道和1 梯度幅值,见下图,这些通道可以高效计算并且捕获输入图像不同的信息。

    在这篇文章里,AdaBoost分类器采用了soft cascade的级联方式。为了检测不同大小的行人,作者并没有进行图像缩放然后用固定大小的分类器扫描,而是训练了几个典型尺度大小的分类器,对于其他尺度大小的行人,采用这些典型尺度分类器的预测结果进行插值来逼近,这样就不用对图像进行缩放。因为近处的行人和远处的行人在外观上有很大的差异,因此这样做比直接对图像进行缩放精度更高。这一思想在后面的文章中还得到了借鉴。通过用GPU加速,这一算法达到了实时,并且有很高的精度,是当时的巅峰之作。

    2.4 DPM+ latent SVM

    行人检测中的一大难题是遮挡问题,为了解决这一问题,出现了采用部件检测的方法,把人体分为头肩,躯干,四肢等部分,对这些部分分别进行检测,然后将结果组合起来,使用的典型特征依然是HOG,采用的分类器有SVM和AdaBoost针对密集和遮挡场景下的行人检测算法可以阅读文献[15]。

    DPM(Deformable Parts Models)算法在SIGAI在之前的文章“基于深度学习的目标检测算法综述”已经提到过。这是是一种基于组件的检测算法,DPM检测中使用的特征是HOG,针对目标物不同部位的组建进行独立建模。DPM中根模型和部分模型的作用,根模型(Root-Filter)主要是对物体潜在区域进行定位,获取可能存在物体的位置,但是是否真的存在我们期望的物体,还需要结合组件模型(Part-Filter)进行计算后进一步确认,DPM的算法流程如下:

    DPM算法在人体检测中取得取得了很好的效果,主要得益于以下几个原因:

    1. 基于方向梯度直方图(HOG)的低级特征(具有较强的描述能力)
    2. 基于可变形组件模型的高效匹配算法
    3. 采用了鉴别能力很强的latent-SVM分类器

    DPM算法同时存在明显的局限性,首先DPM特征计算复杂计算速度慢(论文[8]中针对DPM提出了多个加速的策略,有兴趣的读者可以参考);其次,人工特征对于旋转、拉伸、视角变化的物体检测效果差。这些弊端很大程度上限制了算法的应用场景,这一点也是基于人工特征+分类器的通病。

    采用经典机器学习的算法虽然取得了不错的成绩,但依然存在下面的问题

    1.对于外观,视角,姿态各异的行人检测精度还是不高

    2.提取的特征在特征空间中的分布不够紧凑

    3.分类器的性能受训练样本的影响较大

    4.离线训练时的负样本无法涵盖所有真实应用场景的情况

    基于机器学习的更多方法以参考综述文章[10][18][19]。文献[10]对常见的16种行人检测算法进行了简单描述,并在6个公开测试库上进行测试,给出了各种方法的优缺点及适用情况。

    文献[18]提出了Caltech数据集,相比之前的数据集,它的规模大了2个数量级。作者在这个数据集上比较了当时的主要算法,分析了一些失败的的原因,为后续的研究指出了方向。

    文献[19]也比较了10年以来的行人检测算法,总结了各种改进措施,并支持了以后的研究方向。

    3. HOG+SVM环境配置

    3.1 数据集INRIADATA

    来自HOG+SVM的论文,论文中使用的图片是归一化之后的。

    训练和测试的时候有几种选择方式:

    用normalized_images目录下的图片做训练,或者用original_images目录下的图片+annotations获取行人区域做训练;测试则都在original_images/test/pos上测试。

    这个数据集来自于 github搜索关键字:HOG+SVM ,或者github搜索关键字:HOG+SVM行人检测,后者star最多的的 https://github.com/icsfy/Pedestrian_Detection,找不到原作者了。反正就是作者已经分好了如下的数据集:

    3.2 算法原理

    资料:https://blog.csdn.net/hongbin_xu/article/details/79810544 HOG特征检测学习笔记,MATLAB code求HOG。

    3.2.1 HOG梯度方向直方图 Histogram of Oriented Gradients

    • 特征描述子:特征向量作为SVM的输入很重要,选取怎样的特征决定了分类的准确率。
    • 哪些特征是有用的呢?哪些特征冗余?特征需要有很好的区分能力。
    • HOG特征描述子,选用梯度方向的分布作为特征。梯度是描述边缘信息的。

    如何计算HOG梯度方向直方图呢?

    • 图像需要保持固定的长宽比,如1:2。从原始图像得到小的切片。如得到64*128的切片
    • 再把64*128的切片分成8*8的Cell。用特征描述子描述图像的一小块,用更细微的方式刻画了原图像。
    • 梯度方向图+梯度强度图,方向在[0,180],经验表明:无符号的梯度比有符号梯度效果更好。
    • 计算梯度方向直方图:每个bin是按照梯度方向选出来的。一般是9个bin
    • 总结:8*8的cell用9*1的直方图表示,4个cell组成一个block,block为16*16。所以一个block有4*9*1=36*1的列向量表示HOG特征。那么对于一个64*128的图像而言,其特征的维度为:7*15*36=3780。

    梯度方向直方图Histogram of Oriented Gradients (HOG) 这篇文章详细描述了HOG特征是如何提取和计算的,超级棒!

    关键的几张图示如下:

    这里写图片描述

    这里写图片描述

    这里写图片描述

    这里写图片描述

    其他资料:

    https://blog.csdn.net/carson2005/article/details/8316835   CSDN 266

    https://blog.csdn.net/hongbin_xu/article/details/79845290  CSDN 3081 code

    https://blog.csdn.net/qq_37753409/article/details/79047063    CSDN C++ code

    https://www.cnblogs.com/heleifz/p/train-hog-svm-detector.html

    • Vlfeat有一个HOG的实现,比OpenCV的实现要好。
    • 按照这种策略训练出来的检测器,在 INRIA 上达到了 Precision:50%,Recall:71% 的性能,超过了 OpenCV 自带的 HOG Detector(其实没有可比性,因为我用了更好的 HOG 特征),如果您有更好的方案,请告诉我,谢谢。

    https://blog.csdn.net/Young__Fan/article/details/85000200   HOG+SVM参考资料大合集

    3.2.2 使用线性SVM进行训练

    提取特征+分类器,提取了好的特征,那么用什么样的分类器都能得到不错的结果。

    为什么使用线性SVM?

    另外一个问题什么时候选择SVM或逻辑回归?如图:
    1. 特征维度n很大:
    使用逻辑回归和线性SVM.
    2. 特征维度n小,样本数量m中等:
    使用高斯核SVM。
    3. 特征维度n小,且样本数量m巨大:
    - 可以创建新的特征
    - 然后使用逻辑回归和无核SVM

    什么是线性SVM?

    左图就是线性SVM,就是不使用核函数,将\theta^TX送入SVM,形成一个大间距分类器。后者使用了高斯核函数,可以得到复杂的曲线边界。

    image

    • 核函数:为了得到非线性的拟合

    假设函数hθ(x)=θ0+θ1f1+θ2f2+θ3f3+⋯(用f代替x的参数)这个函数为新的假设函数。

    我们设置f函数,衡量标记点和原先样本点的相似性。 为fi=similarity(x^{(i))},l^{(i)})=exp(-\frac{\left \| x^{(i))}-l^{(i)}\right \|^2}{2\sigma^2})

    即设计一些landmarks li,这些landmarks对于得到非线性的边界十分重要,计算xi和li的相似度fi,并将fi作为特征输入:

    hθ(x)=θ0+θ1f1+θ2f2+θ3f3+⋯,就可以得到非线性的分类边界!

    • 如何得到参数θ呢?

      • 输入X,y,C等参数,经过svmTrain后得到模型model,得到参数w和b,由于是线性的,所以y=(w^Tx+b)

      • 则:w^Tx+b>1为正样本,w^Tx+b<-1为负样本。

    model = svmTrain(X, y, C, @linearKernel, 1e-3, 20);
    
    function visualizeBoundaryLinear(X, y, model)
    %VISUALIZEBOUNDARYLINEAR plots a linear decision boundary learned by the
    %SVM
    %   VISUALIZEBOUNDARYLINEAR(X, y, model) plots a linear decision boundary 
    %   learned by the SVM and overlays the data on it
    
    w = model.w;
    b = model.b;
    xp = linspace(min(X(:,1)), max(X(:,1)), 100);
    yp = - (w(1)*xp + b)/w(2);
    plotData(X, y);
    hold on;
    plot(xp, yp, '-b'); 
    hold off
    
    end

    3.3.3 代码

    4. 环境配置:Matlab安装libsvm

     资料:


    • 1. 下载libsvm,下载地址为:https://www.csie.ntu.edu.tw/~cjlin/libsvm/找到这个download LIBSVM 下载
    • 2. 下载得到的libsvm-3.24文件夹复制到E:\Program Files\matlab r2018\toolbox下
    • 3. 在matlab的主页栏找到 设置路径 ,添加 E:\Program Files\matlab r2018\toolbox\libsvm-3.24\matlab 到 matlab的路径中

    • 4. 在MATLAB命令窗口输入mex –setup C++ ,

    • 5. matlab目录切换到E:\Program Files\matlab r2018\toolbox\libsvm-3.24\matlab 文件下,并执行make

    • 5.解决办法:——之前没有解决这个问题,好像是权限导致的,使用管理员方式打开matlab
    • 6. 不明所以,直接找matlab的 安装MinGW-w64,转向了matlab file exchange,登录matlab账号进行下载,得到
    • 直接将这个文件拖到matlab的命令行,然后就是漫长的等待安装下载

    参考 Matlab选择mingw编译器 进行安装。

    寻找新办法

    https://blog.csdn.net/SKY_yiyi_9/article/details/88140283 如何在Matlab2018a中配置MinGW-w64 C/C++ 编译器

    无奈放弃,让别人跑的,生成.mex64给我,哎!

    • 猜测我是vs2013,别人的时vs2015,哎!
    • 然后要添加文件夹,添加libsvm的很多文件夹到当前路径下面
    • 然后使用的时候为了与matlab自带的函数冲突,使用libsvmtrain等。使用的时候也用libsvmtrian即可。解决问题

    5. 行人检测方面的数据集

    https://blog.csdn.net/haronchou/article/details/106120804 其他数据集

     1. HOG+SVM使用的行人识别数据集

    • (1) INRIA Person Dataset(INRIA行人数据库)——可见光数据集,样本大小128*64

    2. 红外行人数据集:KAIST Multispectral Pedestrain Detection Benchmark

    • github地址:https://github.com/SoonminHwang/rgbt-ped-detection 
    • 这个是多光谱图像:对齐的RGB+热像仪图像,640×480像素,手动注释所有人,骑自行车的人。注释包括边界框(如Caltech Pedestrian Dataset)之间的时间对应。
    • 许多人也在努力提高我们的基准行人检测性能,如:
    • 其他研究人员也采用多模态的方法,Also, another researches to employ multi-modality are presented.

      • Image-to-image translation [Arxiv '17]
      • Calibrations

    3. 南方科技大学的红外数据集:SCUT_FIR_Pedestrian_Dataset

    • 执行之前,先运行一下starup,将所有的文件夹加入当前路径,才能方便的调用其他函数
    •  extract_img_anno_scut
    • dbExtract_scut 从数据库seq提取images,同时提取annotations到Txt files
    • 从pth的annatation里面load ground truth,即annatations文件夹要放在pth的路径下面
    • 然后将提取到的Images放在tDir的images文件夹下面
    • seq文件要放在pth的videos文件夹下面。
    • 然后运行 extract_img_anno_scut函数即可

    得到的结果如下:

    import cv2
    import os
    import sys
    import glob
    from PIL import Image
    
    # VEDAI 图像存储位置
    src_img_dir = "F:\\pedestrain detection benchmark\\SCUT_FIR_Pedestrian_Dataset-master\\train02\\images\\"
    # VEDAI 图像的 ground truth 的 txt 文件存放位置
    src_txt_dir = "F:\\pedestrain detection benchmark\\SCUT_FIR_Pedestrian_Dataset-master\\train02\\annotations\\"
    # 所有的图像名称
    img_Lists = glob.glob(src_img_dir + '/*.jpg')
    
    img_basenames = []  # e.g. 100.jpg
    for item in img_Lists:
        img_basenames.append(os.path.basename(item))
    
    img_names = []  # e.g. 100
    for item in img_basenames:
        temp1, temp2 = os.path.splitext(item)
        img_names.append(temp1)
    
    for img in img_names:
        # open the crospronding txt file
        with open(src_txt_dir + '/' + img + '.txt',"r") as f:
            line_count = 0
            for line in f.readlines():
                if line_count == 0:
                    line_count += 1
                    continue
                line = line.strip('\n')
                spt = line.split(' ')
                cat_name = str(spt[0])
                x_min = int(spt[1])
                y_min = int(spt[2])
                bbox_wid = int(spt[3])
                bbox_hei = int(spt[4])
    
                if cat_name == "people" or cat_name == "walk person":
                    if bbox_hei >= 10 and bbox_hei <= 30:
                        bbox_hei_new = 2*bbox_wid
                        if y_min+bbox_hei > 576:
                            y_min = 576-bbox_hei_new
                        else:
                            y_min = y_min-(bbox_hei_new-bbox_hei)/2
    
                        # 确保尺寸为1:2,扩充为原来的1.5倍
                        '''
                        bbox_wid_new = int(bbox_wid)
                        bbox_hei_new = int(2 * bbox_wid_new)
                        # 防止越界
                        if y_min >= 0.25*bbox_hei and y_min + bbox_hei_new< 576:
                            y_min = y_min - 0.25*bbox_hei
                        elif y_min < 0.25*bbox_hei:
                            y_min = 0
                        else:
                            y_min = 576 - bbox_hei_new
    
                        if x_min >= 0.25 * bbox_wid and x_min + bbox_wid_new < 720:
                            x_min = x_min - 0.25 * bbox_wid
                        elif y_min < 0.25 * bbox_wid:
                            x_min = 0
                        else:
                            x_min = 720 - bbox_wid_new
                        '''
    
                        # 读入图像
                        img_data = cv2.imread((src_img_dir + '/' + img + '.jpg'), cv2.IMREAD_UNCHANGED)
                        cropped = img_data[int(y_min):int(y_min+bbox_hei_new),int(x_min):int(x_min+bbox_wid)]  # y为高度,x为宽度
                        cv2.imwrite("G:\\scut_pedestrain_crop\\train\\" + "test_" + str(img) + ".jpg", cropped)
                        # cv2.rectangle(img_data, (x_min, y_min), (x_min+bbox_wid, y_min+bbox_hei), (255, 0, 0), 2)
                        # cv2.imshow('label', img_data)
                        # cv2.waitKey(0)
    
    
            '''
        # gt = open(src_txt_dir + '/gt_' + img + '.txt').read().splitlines()
    
        # write the region of image on xml file
        for img_each_label in gt:
            spt = img_each_label.split(' ')  # 这里如果txt里面是以逗号‘,’隔开的,那么就改为spt = img_each_label.split(',')。
    
    
            xml_file.write('        <name>' + str(spt[4]) + '</name>\n')
            xml_file.write('        <pose>Unspecified</pose>\n')
            xml_file.write('        <truncated>0</truncated>\n')
            xml_file.write('        <difficult>0</difficult>\n')
            xml_file.write('        <bndbox>\n')
            xml_file.write('            <xmin>' + str(spt[0]) + '</xmin>\n')
            xml_file.write('            <ymin>' + str(spt[1]) + '</ymin>\n')
            xml_file.write('            <xmax>' + str(spt[2]) + '</xmax>\n')
            xml_file.write('            <ymax>' + str(spt[3]) + '</ymax>\n')
            xml_file.write('        </bndbox>\n')
            xml_file.write('    </object>\n')
    
    
    # 获取路径下的所有txt格式的文件名
    def file_name(file_dir):
        File_Name=[]
        img_dir = "F:\\pedestrain detection benchmark\\SCUT_FIR_Pedestrian_Dataset-master\\train02\\images\\"
        for files in os.listdir(file_dir):
            if os.path.splitext(files)[1] == '.txt':
                img_dir += files[:-3]+"jpg"
                # File_Name.append(files)
                # 以原格式读入图像,原图像为灰度图
                img = cv2.imread(img_dir,cv2.IMREAD_UNCHANGED)
                # 需要读入标签信息
    
    
    
        return File_Name
    txt_file_name=file_name("F:\\pedestrain detection benchmark\\SCUT_FIR_Pedestrian_Dataset-master\\train02\\annotations\\")
    print("txt_file_name",txt_file_name)
            '''

    6. 在自己的数据集上分析和优化

    背景:由于HOG+SVM行人检测在自己的数据集上出现了问题,如下:

    • 正样本召回率95%,负样本召回率59%——导致Recall=95%,Precision=59%,F1=0.6
    • 需要对模型进行分析,到底是哪里出了问题,导致了现在的情况?
    • 初步认为借助:学习曲线等一系列诊断工具可以分析出是特征不够?不具有区分性?还是样本不够?还是样本的结构不太好?
      • (1)现象:模型对样本没有任何区分性,即对负样本的预测上没有任何趋势
      • (2)原因:未知
      • (3)怎么办?:暂时准备给出一些clean干净的负样本数据,降低负样本识别的难度。——临时措施
        • 目的1:仅仅是为了测试当前模型对于easy task场景能否胜任
        • 目的2:若easy场景可以适用,则写入工程中测试初步的效果,然后出第一版
      • (4)下一步优化目的:要加入midum task,主要是midum的负样本,如何提高负样本的识别率。
      • (5)下一步优化方向:肯定是要优化HOG,主要是分析问题出在特征的哪个环节,后者负样本数据集结构不合理的哪个环节,以及数据尺寸问题的哪个环节这部分复杂又重要,特征的选择、分析、优化决定了一切!

    后记:坑爹的原因是读图的路径设置的有问题,导致读入的仅仅是文件名,没有加入绝对路径,导致读的图是错误的,才出现了上述问题。我也是醉了!!!以后要避免这种低级错误。

    当然,模型的分析和优化技能也是相当重要的,是后期最重要的技能。所以还是要掌握的!

    2020/05/15 update:

    • 当前问题:样本太少太少。项目数据量在200+左右,且需要识别400m左右的人,<200m的行人。

    1. 如何学习怎么查找模型问题所在?

    • 参考学习别人的资料:《机器学习实战》chap6-SVM支持向量机的代码学习,学习如何使用python的numpy等计算包
    • 参考《机器学习实战》chap6如何计算SVM的训练误差和测试误差的,然后计算学习曲线等诊断模型
    • 模型诊断的手段和方法有待学习!!!!

    《机器学习实战》中英文电子书+源码下载:https://redstonewill.com/1884/

    1.1 SVM深入理解

    资料:Python3《机器学习实战》学习笔记(八):支持向量机原理篇之手撕线性SVM,SVM原理详解,公式推导等

    • 线性可分:棍子;非线性可分:kerneling;trick: Optimization
    • 哪个分类器更好?:拥有更大分类间隔的分类器更鲁棒!有左右两个极限位置:越过该位置就会错分
    • 具有最大间隔的决策面就是SVM寻找的最优解。
    • 两侧虚线所穿过的样本点,就是SVM中的支持样本点,称为”支持向量

    最优化过程的数学建模:

    • (1)目标函数——分类间隔最宽;希望什么东西的什么指标达到最好

    • (2)优化对象——通过移动决策面。期望通过改变哪些因素来使得目标函数达到最优

    对上述两个对象进行建模,“决策面方程”w^TX+b=0,分类间隔方程:|w^TX+b|=1

    目的是找出一个分类效果好的超平面作为分类器,分类间隔越大,这个超平面的分类效果越好。所以变成了求解分类间隔d最大的分类器。

    • (3)约束条件:虽然有了目标函数的数学形式,为了求解最大的d,仍有如下问题:即我们的要优化的变量d受到下面条件的限制和约束:
      • 如何判断超平面将样本点正确分类?
      • 向求距离d的最大值,首先要找到支持向量上的点,怎么在众多的点中选出支持向量上的点呢?

    针对问题1,假设决策面在间隔区域的中轴线上,且支持向量的样本点到决策面的距离为d,约束公式如下:

    第三个公式同时处以了d,也就是说样本点到决策面的距离>d,才预测为1 or -1.

    (1)(2)(3)

    • (4)线性SVM优化问题描述

    目的是使得d最大化(最大间隔的分类器),d的表达式如下。支持向量又满足如下式子:

    我们只关心支持向量上的点,因此我们的目标函数简化为:d=\frac{1}{\left \| w \right \|},那么求max{d}变成了求min{||w||},进一步为min{||w||^2}。等价于支持向量所在的分类间隔要最大!所以d的分子才能为1,这里做了一个概念上的等价替换!

    那么就有如下的目标函数和约束条件:分类间距最大,同时满足分类的约束条件:

    上述带有不等式约束的优化式子常用的方法是KKT条件。把所有的等式、不等式约束与目标函数f(x)写成一个式子,即拉格朗日乘数法。通过KKT条件,可以求出最优值。KKT条件是可以取得最优值的必要条件。

    随后,人们又发现,使用拉格朗日获得的函数,使用求导的方法求解依然困难。进而,需要对问题再进行一次转换,即使用一个数学技巧:拉格朗日对偶

    所以,显而易见的是,我们在拉格朗日优化我们的问题这个道路上,需要进行下面二个步骤

    • 将有约束的原始目标函数转换为无约束的新构造的拉格朗日目标函数
    • 使用拉格朗日对偶性,将不易求解的优化问题转化为易求解的优化

    将上面的有约束的目标函数转为如下的无约束的拉格朗日目标函数为:

    其中αi是拉格朗日乘子,αi大于等于0,是我们构造新目标函数时引入的系数变量(我们自己设置)。现在我们令:

    \theta(w)=max_{\alpha_i\geqslant 0}(\frac{1}{2}\left \| w \right \|^2-\sum_{i=1}^{n}\alpha_i(y_i(w^Tx_i+b)-1))

    每一个样本都作为约束条件加上去!

    这部分的内容要详细推到,一直到SMO怎么求解的,关于SVM的深入分析和理解。暂时放在这里。我当前需要了解如何优化模型,如何分析特征,如何评价和诊断模型。

    机器学习——支持向量机SVM(Support Vector Machine)(下)

    SVM训练的模型的复杂度由支持向量的个数决定的,而不是数据的维度决定。SVM不太容易过拟合。

    SVM训练出的模型完全依赖于支持向量。

    一个SVM如果训练得到的支持向量个数比较少,则泛化能力越强。

    1.3 SVM中解决泛化问题

    机器学习模型的泛化能力不足,有什么改进思路?

    2. Matlab svm训练的模型如何在C++中使用

    参考网址:http://www.jeepxie.net/article/646115.html

    • 失败,放弃!决定将代码移植到C++中​​​​​​​

    7. 分析和优化模型

    当前状态:

    • 算法使用HOG+SVM对检测的结果进行了二分类。算法选择上:运动目标建模形成切片+HOG-SVM二分类
    • 训练集:正样本为400m左右的行人251张,负样本100张为检测里面的其他部分,包括误检等
    • 验证集:92张正样本,109张负样本
      • 当前唯一实现是:训练集和验证集是绝对的不重合,不存在重叠元素。且来自同一分布。分别从相同场景下提取的前景目标。
      • 最大问题是:数据量过小,明显的训练集过拟合了,模型泛化能力不行。
      • 最大问题:验证集的设计太小,且不合理。

    • 吴恩达深度学习中表明:数据量<1W的,Train/Dev/Test的比例在60%、20%、20%,即6K训练集,2k验证集和2k测试集
    • Train和Dev,Test的数据源于同一分布。即将所有的数据重新洗牌,放入验证集和测试集。
    • 基于大量数据的简单模型 > 基于少量数据的复杂模型
    • 更多的数据 胜于 聪明的算法
    • 好的数据 胜于 多的数据

    • 数据集数据量不够是一个比较大的问题,有哪些办法可以解决呢?​​
      • (1)多采集一些数据量
      • (2)去公开的红外数据集爬取一些可用的数据集
      • (3)利用别的数据集来训练,用自己的数据集测试。也许迁移学习是可以的哦???当然保持训练集和验证集分布一致是最好的了。
    • 如果已经获得了足够的数据,又该怎么办呢?进行模型优化四部曲
      • (1)训练集指标:过拟合 or 欠拟合?
      • (2)验证集指标:泛化能力不好?增大训练集的正则化?增大训练集的样本数?
      • (3)测试集指标:泛化能力不好,验证集数据增加,或者...
      • (4)Real Life表现:回去改变验证集或代价函数。

    2.其他的一些疑问 

    • (1)把验证集里面不通过的hard example直接丢进训练集再训练一遍?
      • 嘿,这不是自欺欺人嘛。选择验证集的目的就是为了测试算法效果,你这样的话还用验证集有什么用呢?不就没用了嘛。所以大多数比赛的测试集都是不公开的啊
    • (2)更多的数据 > 聪明的算法,好的数据 > 多的数据
      • 对于红外数据而言,数据冗余是否是一个问题因为本来目标就小,姿态什么的多样性就比较少!一个人换来换去就那几个姿势和形态,都放进去吗?会不会数据冗余?
        • 这个应该就属于是否是好的数据了,如果冗余的数据不多的话,那么就会是好的数据!
      • 那怎么能够具有区分性呢?现在是靠hog去提取行人的特征,提取到的特征到底是怎么样的,也没有去看?

    3. 红外数据集

    https://blog.csdn.net/hunnzi/article/details/103445059

    (1)KAIST——用于行人检测任务

    • 双光系统,长波红外+可见光。长波红外相机:320*256像素,FLIR-A35。可见光相机:640*480,帧率20fps。
    • 75%的行人输入Medium尺度,在[40,128]pixel;10.38%的行人<20pixel;13.67%的行人>128pixel
    • 有四种标签:person(好分辨);people(难分辨);cyclist(骑车的人);person?(人工不能分辨)。
    • 共计95,328个可见光-红外图像对,总数为103,128个标注,1182个人。

    获取地址:https://blog.csdn.net/Lcd_2018_7_18/article/details/103064150 

    1. 红外数据集

    http://vcipl-okstate.org/pbvs/bench/ OTCBVS数据集里面有用的部分

    Dataset 01: OSU Thermal Pedestrian Database

    Topic of Interest:
    Person detection in thermal imagery. 热成像行人检测 75mm镜头,8楼的房顶拍摄,10个序列共284张图,8-bit灰度图,360*240像素,有ground-truth的txt以及matlab程序可以画bounding box,不包含遮挡的行人,50%可见的行人才进入ground-truth中

    Sensor Details:
    Raytheon 300D thermal sensor core
    75 mm lens 
    Camera mounted on rooftop of 8-story building
    Gain/focus on manual control

    Data Details:
    Pedestrian intersection on the Ohio State University campus
    Number of sequences = 10
    Total number of images = 284
    Format of images = 8-bit grayscale bitmap
    Image size = 360 x 240 pixels
    Sampling rate = non-uniform, less than 30Hz
    Environmental information for each sequence provided in subdirectories
    Ground truth provided in subdirectories as list of bounding boxes (with approximately same aspect ratio) around people.
    For the ground truth data, we selected only those people that were at least 50% visible in the image (i.e., highly occluded people were not selected).

    Dataset 03: OSU Color-Thermal Database

    彩色和热成像的双光系统,双光融合的目标检测。红外23mm镜头。6个视频序列共17089张图像,8-bit灰度图,320*240pixe的红外图像。还有跟踪的result的文件。

    Fusion of color and thermal imagery,
    Fusion-based object detection in color and thermal imagery

    Sensor Details:
    Thermal Sensor: Raytheon PalmIR 250D, 25 mm lens
    Color Sensor: Sony TRV87 Handycam

    Cameras mounted adjacent to each other on tripod at two locations approximately 3 stories above ground
    Gain/focus on manual control

    Data Details:
    Busy pathway intersections on the Ohio State University campus
    Number of color/thermal sequences = 6 (3 at each location)
    Total number of images = 17089
    Format of images = Thermal: 8-bit grayscale bitmap, Color: 24-bit color bitmap
    Image size = 320 x 240 pixels
    Sampling rate = approx. 30Hz
    Color/Thermal images registered using homography with manually-selected points
    Files containing tracking results on the dataset are provided by Alex Leykin

    8. 指标方面

     1. 图像指标对应 

    为了明确图像指标上的对应,给出如下图示。

    • flir的行人<30pixel的与我们的工程数据<30pixel的不一样
    • 同样height都是30pixel,flir的人看着要近一些。应该是flir的分辨率比1280的要小些。
    • flir只有<30pixel的数据才有意义,因为对应着300m的行人
    • 20pixel height以下的数据就很难识别了,具有轮廓信息还是非常重要的。
    • 300m-30pixel, 400m-20pixel, 500m-16pixel 且降低误识别才是最重要的!
    • 先做<30pixel height的行人识别,重点在于降低误检率!验证集很重要!
    % 行人高度、镜头、距离对应得到的Pixel
    
    height = 1.6;%m
    distance = 700;%m
    
    f = 55;%55mm镜头
    size_p = 12*0.001;%12um像元尺寸12um*10e-3=0.012mm
    
    pixel = f*height / (size_p * distance);
    fprintf("距离=%d m \t pixel_height=%f\n",distance,pixel);

     

    理论计算距离与pixel的对应 

    300m                      400m                                  500m     


    2. 正负样本的问题

    以下来自:https://stackoverflow.com/questions/25598335/collect-negative-samples-of-adaboost-algorithm-for-face-detection

    • 负样本的选取与问题场景有关。不能是所研究问题毫不相关的场景图片,这样的负样本没有意义
    • 如果不知道应用环境,那就尽量采集多种环境下的图像。
    • 负样本之间不能有太多的相关性
    • 负样本的尺寸应该跟正样本对齐
    • 如果老对一些负样本误识别,可以将这样的负样本丢进训练集再训练一次
    • 负样本的尺寸要一样大。所有的正样本和负样本都要resize到相同的尺寸,resize to exact resolution,same format.
    • 负样本每个stage以两倍的数量递增!负样本为:1k,2k,4k,8k,16k,...,512k,这就是10个stage,累计就是1,024k=102.4W的样本。

    分类问题中的正负样本问题,如人脸识别中的例子。负样本的选取与问题场景相关,具体而言,如果你要进行教室中学生的人脸识别,那么负样本就是教室的窗子、墙等等,也就是说,不能是与你要研究的问题毫不相关的乱七八糟的场景图片,这样的负样本并没有意义,还有一个比较好的网址是(http://www.doc.ic.ac.uk/~sgc/teaching/pre2012/v231/lecture10.html)

    2.1 样本的数量

    训练一个检测器时,每个stage应该要有几千张正负样本。通常一个检测器会有10~20个stage。每个stage会依次递减负样本的数量为上个stage的一半。所以你大概需要3k-10k的正样本,5,000,000~100million的负样本(500万-10亿)。

    2.2 负样本

    负样本与环境相关,如果你要进行教室中学生的人脸识别,那么负样本就是教室的窗子、墙等等,也就是说,不能是与你要研究的问题毫不相关的乱七八糟的场景图片,这样的负样本并没有意义。如果你不知道你对应用环境,那就尽量采取多种可能的不同的环境图像。

    2.3 应该把脸部的局部部位作为负样本吗?如眼睛或鼻子

    你可以把眼睛和鼻子作为负样本,但是只将这些作为负样本远远不够。一个探测器的real strength真正韧劲will come from negative images which represent the typical background of the faces。一个探测器的真正韧劲来源于负样本,这些负样本代表了人脸的典型背景

    2.4 该怎么收集负样本信息?

    你其实并不需要很多负样本,你可以采取1000张图像,然后从中产生10,000,000张负样本。假设你有一张汽车的图片,像素为1k×1k pixels。假设你的人脸识别检测器工作在20×20pixel进行检测,如OpenCV就是这样的。那么你可以将1k×1k的汽车图片剪切成20×20的切片。你也可以获得2500个切片(50×50)。

    2.5 在步骤4中。

    我必须提醒你注意,负样本之间不能有太多相关性,所以我不建议只从一张图片中截取和产生1百万个负样本。create a library of 1000 images 或从google随机下载数据,确定所有的负样本都不包含人脸。从每个图像中crop截取10,000个负样本,这样得到了大约10million的负样本。训练你的检测器。下一步你可以将每一张图像cut到~50,000张部分重叠的pieces,这样就将负样本扩大到50millions。这样你就会有很好的结果了。

    2.6 探测器的最后增强步骤

    当你最后已经有一个很好的探测器之后,那么在更多的图像上运行试一下。它会产生误识别。收集所有这些误识别的检测图像,然后将他们加到你的负样本中。然后将检测器再训练一遍。这样的迭代次数越多,你的检测器会更鲁棒。

    2.7 real numbers

    今天最好的face detectors(如Facebooks)使用了hundreds of millions几百个million(1亿)的正样本,billions10亿负样本。正样本不仅仅是正面的脸,还有各个方向的脸,各种表情的脸,各种年龄的脸,各种性别的脸,各种种族的脸。带着眼镜/帽子/太阳镜/化妆的脸等。你不必要跟最好的比,所以miss掉一些faces时不要get angry。

    3. 负样本的concerns

    3.1 不同环境下的负样本

    • 如何考虑到不同环境下的使用情况?不同对比度下的~
    • 还是根据用户场景,定制一套流程,进行负样本收集等?

    3.2 不同尺度(距离)下的负样本

    <100m,200m,300m,400m下的负样本问题,不同距离的负样本的尺度也是不一样的,比如车?比如误检测?

    4、目标检测——降低误检测率:负样本制作、增强、训练

    目标检测——降低误检测率:负样本制作、增强、训练

    如何从负样本的角度去提高模型的识别精度呢?解释如下:

    为什么要训练负样本?

    • 训练负样本的目的是为了降低误检测率、误识别率,提高网络模型的泛化能力。通俗地讲就是告诉检测器,这些“不是你要检测的目标”。
    • 负样本的核心目的是:增强泛化能力!

    如何收集负样本?

    可以通过下面两种方式收集负样本:

    • 采用本任务场景的不包含目标物体的背景图像,例如你的目标是识别某园区内的行人,那么所有本园区内不包含行人的图片都视作负样本。不要使用不属于本任务场景的负图像,因为其对检测器的性能影响不大。
    • 测试图像中被识别错误的目标所在区域。(通常对原图像进行裁剪,使得裁剪下来的图像只包含误识别的物体,而不包含目标)
    • 正负样本必须放在一起训练,不能单独训练负样本,否则经过训练,网络会把所有的图像都识别为背景。
    • 正负样本的比例最好为1:1到1:2左右,数量差距不能太悬殊,特别是正样本数量本来就不太多的情况下。

    怎么消除误检测?

      把使用正样本训练好的模型拿来进行测试,此时会得到一些被错误识别的图片。把这些图片收集起来作为负样本加入到正样本集(如果图片中同时包含误识别物体和目标,可以将图像裁剪,裁剪后的图像包含误识别物体而不包含目标并尽量覆盖原图大部分区域,然后再将其分辨率resize回原图大小),组成新的训练集,送入模型进行训练。

           如果负样本的来源只有误识别的图片,那么由于误识别的图片往往占少数,可以利用图像增强(如高斯滤波、对比度增强、锐化、平滑等)的方法扩充负图像数量至和正样本数量相同,并组合在一起。将这样得到的训练集送入模型进行训练,经过若干个epoch当Loss收敛到稳定值时,再次测试原来的出现误识别的图像你会发现误识别现象基本消失了,并且类似原来误识别的场景将会被正确识别。

    5、评价指标

    虚警概率:false positive / (false positive + true positive),就是1-precision=false alarm。

    6、大规模训练结果

    • (1)多采集项目的数据,制作样本切片,达到1k-2k
    • (2)特别是负样本的制作,为正样本的2倍,考虑环境和其他目标。剔除<10pixel的数据。
    • (3)然后再考虑是否进行特征和调参
    展开全文
  • SVM 适合中小型数据样本、非线性、高维的分类问题。它将实例的特征向量映射为空间中的一些点。如: 而SVM要做的事情就是找到那么一条线, “最好地” 区分这两类点,以后有了新的点,这条线也能做出很好的分类。...
  • 支持向量机原理篇之手撕线性SVM

    千次阅读 2018-07-07 16:00:09
    本篇文章参考了诸多大牛的文章写成的,对于什么是SVM做出了生动的阐述,同时也进行了线性SVM的理论推导,以及最后的编程实践,公式较多,还需静下心来一点一点推导。 本文出现的所有代码,均可在我的github上下载,...
  • 多类潜在局部线性 SVM 版权所有 (c) 2013 Idiap 研究所, 作者: < > Idiap研究所, Center du Parc,邮政信箱 592, 马可尼街 19, 1920 年,瑞士马蒂尼 电话:+41 27 721 77 57 传真:+41 27 721 77 12 ...
  • 为了增强局部线性嵌入(LLE)算法对人脸识别中特征的分类性能,将最小生成树算法思想引入,提出一种邻域参数动态变化的新的局部线性嵌入算法。该算法采用单链聚类算法以及对其进一步优化自动确定数据点邻域,改善了...
  • LLE局部线性嵌入算法

    千次阅读 2016-07-09 14:50:25
    LLE局部线性嵌入算法 导语:很久没发博文了,今日抽个小空,整理下上个学期做过的东西,写成博文,供给初学者参考: 一、LLE算法概念 LLE算法,即局部线性嵌入算法,是一种非线性降维算法,它利用线性重构的局部...
  • 线性SVM线性SVM SMO算法 李航老师讲东西向来不拖泥带水,但本节也几乎用了40页的篇幅,足见涉及内容之多。7月底硬着头皮啃了下来,边看边在笔记上摘抄标注,花了2天时间对SVM有了基本认识。一个月过去了,如今...
  • 开发了局部线性嵌入(LLE)算法,以直接从高光谱散射图像数据中提取特征。 应用偏最小二乘判别分析(PLSDA)和支持向量机(SVM)来开发基于LLE,均值LLE和均值光谱算法的分类模型。 基于LLE算法的模型实现了80.4%...
  • 4、线性分类: SVM, Softmax

    千次阅读 2018-07-11 16:09:26
    4、线性分类 上一课最我们介绍了图像分类的问题,任务是从一个固定的类别集合中选一个分配给待识别的图像。最后,我们描述了k-近邻(KNN)分类器,它通过将待标记的图像与训练集中已经标记好的图像进行比较来标记...
  • SVM

    万次阅读 多人点赞 2018-03-05 22:51:51
    SVM是用来解决二分类问题的有监督学习算法,在引入了核方法之后SVM也可以用来解决非线性问题。 一般SVM有下面三种: 硬间隔支持向量机(线性可分支持向量机):当训练数据线性可分时,可通过硬间隔最大化学得一个...
  • 目录极大似然估计最小二乘法线性回归概念公式逻辑回归概述优点代价函数总结线性回归和LR对比SVM概述参数核函数常用核函数类型核函数选择多分类lr svm对比相同点不同点 极大似然估计 极大似然估计,通俗理解来说,...
  • svm

    千次阅读 2018-11-28 23:41:10
    (1) 非线性映射是SVM方法的理论基础,SVM利用内积核函数代替向高维空间的非线性映射; (2) 对特征空间划分的最优超平面是SVM的目标,最大化分类边际的思想是SVM方法的核心; (3) 支持向量是SVM的训练结果,在SVM分类...
  • 前面博客有讲到,样本如果不是线性的可以通过多项式扩展,映射到多维空间来拟合。如此之外,还可以做一个局部加权线性回归(Locally Weighted Linear Regression,LWLR)。 ...
  • 简单的线性回归有可能出现欠拟合的现象,这是由于数据可能不是
  • 监督学习 | SVM线性支持向量机原理 1. 非线性支持向量机 对解线性分类问题,线性分类支持向量机是一种非常有效的方法。但是,有时分类问题是非线性的,这时可以使用非线性支持向量机(non-linear support...
  • 0 引 言  在传感器非线性校正领域,国内外许多学者提出多种方法,并得到广泛...但是,该方法也存在一定的局限性,主要表现在:1)神经网络存在局部极小和过学习问题,易影响网络的泛化能力,因此,对样本的数量和质
  • 支持向量机(SVM)之线性分类

    万次阅读 2018-04-20 22:07:56
      支持向量机(Support Vector Machine, SVM)是曾经打败神经网络的分类方法,从90年代后期开始在很多领域均有举足轻重的应用,近年来,由于深度学习的兴起,SVM的风光开始衰退,但是其仍然不失为一种经典的分类方法...
  • 本文主要介绍图像的线性分类器。首先介绍感知机,并由感知机引出评价函数与代价函数两个概念;比较了SVM 与 Softmax 两类分类器的区别。
  • cs231n-(2)线性分类器:SVM和Softmax

    千次阅读 2016-09-03 21:17:41
    讲解了线性分类器,要理解score function和loss function,以及svm的评分和softmax的概率
  • 常用的分类器包括SVM、KNN、贝叶斯、线性回归、逻辑回归、决策树、随机森林、xgboost、GBDT、boosting、神经网络NN。 常见的降维方法包括TF-IDF、主题模型LDA、主成分分析PCA等等 常用的回归技术:线性回归、多项...
  • 基本线性回归和局部加权线性回归、岭回归、前向逐步回归 in Python
  • 线性SVM可以理解为通过核函数将非线性特征映射为线性特征,然后再求线性分界面(只是方便理解,实际上并不是线性)。实际运算时直接在低维空间中直接计算,而不是先映射到高维空间再进行内积运算 然而,LR算法里...
  • 通过这些知识我们不难发现每个算法都是有很多功能的,这些功能能够更好地帮助大家理解机器学习的相关知识,在这篇文章中我们给大家介绍一下关于SVM线性回归的优缺点。 首先我们给大家介绍一下SVM支持向量机。其实...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 15,551
精华内容 6,220
关键字:

局部线性svm