>>edit svmtrain
>>edit svmclassify
>>edit svmpredict
function [svm_struct, svIndex] = svmtrain(training, groupnames, varargin)
%SVMTRAIN trains a support vector machine classifier
%
% SVMStruct = SVMTRAIN(TRAINING,GROUP) trains a support vector machine
% classifier using data TRAINING taken from two groups given by GROUP.
% SVMStruct contains information about the trained classifier that is
% used by SVMCLASSIFY for classification. GROUP is a column vector of
% values of the same length as TRAINING that defines two groups. Each
% element of GROUP specifies the group the corresponding row of TRAINING
% belongs to. GROUP can be a numeric vector, a string array, or a cell
% array of strings. SVMTRAIN treats NaNs or empty strings in GROUP as
% missing values and ignores the corresponding rows of TRAINING.
%
% SVMTRAIN(...,'KERNEL_FUNCTION',KFUN) allows you to specify the kernel
% function KFUN used to map the training data into kernel space. The
% default kernel function is the dot product. KFUN can be one of the
% following strings or a function handle:
%
% 'linear' Linear kernel or dot product
% 'quadratic' Quadratic kernel
% 'polynomial' Polynomial kernel (default order 3)
% 'rbf' Gaussian Radial Basis Function kernel
% 'mlp' Multilayer Perceptron kernel (default scale 1)
% function A kernel function specified using @,
% for example @KFUN, or an anonymous function
%
% A kernel function must be of the form
%
% function K = KFUN(U, V)
%
% The returned value, K, is a matrix of size M-by-N, where U and V have M
% and N rows respectively. If KFUN is parameterized, you can use
% anonymous functions to capture the problem-dependent parameters. For
% example, suppose that your kernel function is
%
% function k = kfun(u,v,p1,p2)
% k = tanh(p1*(u*v')+p2)
%
% You can set values for p1 and p2 and then use an anonymous function:
% @(u,v) kfun(u,v,p1,p2).
%
% SVMTRAIN(...,'POLYORDER',ORDER) allows you to specify the order of a
% polynomial kernel. The default order is 3.
%
% SVMTRAIN(...,'MLP_PARAMS',[P1 P2]) allows you to specify the
% parameters of the Multilayer Perceptron (mlp) kernel. The mlp kernel
% requires two parameters, P1 and P2, where K = tanh(P1*U*V' + P2) and P1
% >0 and P2 <0. Default values are P1 = 1 and P2 = -1.
%
% SVMTRAIN(...,'METHOD',METHOD) allows you to specify the method used
% to find the separating hyperplane. Options are
%
% 'QP' Use quadratic programming (requires the Optimization Toolbox)
% 'LS' Use least-squares method
%
% If you have the Optimization Toolbox, then the QP method is the default
% method. If not, the only available method is LS.
%
% SVMTRAIN(...,'QUADPROG_OPTS',OPTIONS) allows you to pass an OPTIONS
% structure created using OPTIMSET to the QUADPROG function when using
% the 'QP' method. See help optimset for more details.
%
% SVMTRAIN(...,'SHOWPLOT',true), when used with two-dimensional data,
% creates a plot of the grouped data and plots the separating line for
% the classifier.
%
% Example:
% % Load the data and select features for classification
% load fisheriris
% data = [meas(:,1), meas(:,2)]
% % Extract the Setosa class
% groups = ismember(species,'setosa')
% % Randomly select training and test sets
% [train, test] = crossvalind('holdOut',groups)
% cp = classperf(groups)
% % Use a linear support vector machine classifier
% svmStruct = svmtrain(data(train,:),groups(train),'showplot',true)
% classes = svmclassify(svmStruct,data(test,:),'showplot',true)
% % See how well the classifier performed
% classperf(cp,classes,test)
% cp.CorrectRate
%
% See also CLASSIFY, KNNCLASSIFY, QUADPROG, SVMCLASSIFY.
% Copyright 2004 The MathWorks, Inc.
% $Revision: 1.1.12.1 $ $Date: 2004/12/24 20:43:35 $
% References:
% [1] Kecman, V, Learning and Soft Computing,
% MIT Press, Cambridge, MA. 2001.
% [2] Suykens, J.A.K., Van Gestel, T., De Brabanter, J., De Moor, B.,
% Vandewalle, J., Least Squares Support Vector Machines,
% World Scientific, Singapore, 2002.
% [3] Scholkopf, B., Smola, A.J., Learning with Kernels,
% MIT Press, Cambridge, MA. 2002.
%
% SVMTRAIN(...,'KFUNARGS',ARGS) allows you to pass additional
% arguments to kernel functions.
% set defaults
plotflag = false
qp_opts = []
kfunargs = {}
setPoly = falseusePoly = false
setMLP = falseuseMLP = false
if ~isempty(which('quadprog'))
useQuadprog = true
else
useQuadprog = false
end
% set default kernel function
kfun = @linear_kernel
% check inputs
if nargin <2
error(nargchk(2,Inf,nargin))
end
numoptargs = nargin -2
optargs = varargin
% grp2idx sorts a numeric grouping var ascending, and a string grouping
% var by order of first occurrence
[g,groupString] = grp2idx(groupnames)
% check group is a vector -- though char input is special...
if ~isvector(groupnames) &&~ischar(groupnames)
error('Bioinfo:svmtrain:GroupNotVector',...
'Group must be a vector.')
end
% make sure that the data is correctly oriented.
if size(groupnames,1) == 1
groupnames = groupnames'
end
% make sure data is the right size
n = length(groupnames)
if size(training,1) ~= n
if size(training,2) == n
training = training'
else
error('Bioinfo:svmtrain:DataGroupSizeMismatch',...
'GROUP and TRAINING must have the same number of rows.')
end
end
% NaNs are treated as unknown classes and are removed from the training
% data
nans = find(isnan(g))
if length(nans) >0
training(nans,:) = []
g(nans) = []
end
ngroups = length(groupString)
if ngroups >2
error('Bioinfo:svmtrain:TooManyGroups',...
'SVMTRAIN only supports classification into two groups.\nGROUP contains %d different groups.',ngroups)
end
% convert to 1, -1.
g = 1 - (2* (g-1))
% handle optional arguments
if numoptargs >= 1
if rem(numoptargs,2)== 1
error('Bioinfo:svmtrain:IncorrectNumberOfArguments',...
'Incorrect number of arguments to %s.',mfilename)
end
okargs = {'kernel_function','method','showplot','kfunargs','quadprog_opts','polyorder','mlp_params'}
for j=1:2:numoptargs
pname = optargs{j}
pval = optargs{j+1}
k = strmatch(lower(pname), okargs)%#ok
if isempty(k)
error('Bioinfo:svmtrain:UnknownParameterName',...
'Unknown parameter name: %s.',pname)
elseif length(k)>1
error('Bioinfo:svmtrain:AmbiguousParameterName',...
'Ambiguous parameter name: %s.',pname)
else
switch(k)
case 1 % kernel_function
if ischar(pval)
okfuns = {'linear','quadratic',...
'radial','rbf','polynomial','mlp'}
funNum = strmatch(lower(pval), okfuns)%#ok
if isempty(funNum)
funNum = 0
end
switch funNum %maybe make this less strict in the future
case 1
kfun = @linear_kernel
case 2
kfun = @quadratic_kernel
case {3,4}
kfun = @rbf_kernel
case 5
kfun = @poly_kernel
usePoly = true
case 6
kfun = @mlp_kernel
useMLP = true
otherwise
error('Bioinfo:svmtrain:UnknownKernelFunction',...
'Unknown Kernel Function %s.',kfun)
end
elseif isa (pval, 'function_handle')
kfun = pval
else
error('Bioinfo:svmtrain:BadKernelFunction',...
'The kernel function input does not appear to be a function handle\nor valid function name.')
end
case 2 % method
if strncmpi(pval,'qp',2)
useQuadprog = true
if isempty(which('quadprog'))
warning('Bioinfo:svmtrain:NoOptim',...
'The Optimization Toolbox is required to use the quadratic programming method.')
useQuadprog = false
end
elseif strncmpi(pval,'ls',2)
useQuadprog = false
else
error('Bioinfo:svmtrain:UnknownMethod',...
'Unknown method option %s. Valid methods are ''QP'' and ''LS''',pval)
end
case 3 % display
if pval ~= 0
if size(training,2) == 2
plotflag = true
else
warning('Bioinfo:svmtrain:OnlyPlot2D',...
'The display option can only plot 2D training data.')
end
end
case 4 % kfunargs
if iscell(pval)
kfunargs = pval
else
kfunargs = {pval}
end
case 5 % quadprog_opts
if isstruct(pval)
qp_opts = pval
elseif iscell(pval)
qp_opts = optimset(pval{:})
else
error('Bioinfo:svmtrain:BadQuadprogOpts',...
'QUADPROG_OPTS must be an opts structure.')
end
case 6 % polyorder
if ~isscalar(pval) || ~isnumeric(pval)
error('Bioinfo:svmtrain:BadPolyOrder',...
'POLYORDER must be a scalar value.')
end
if pval ~=floor(pval) || pval <1
error('Bioinfo:svmtrain:PolyOrderNotInt',...
'The order of the polynomial kernel must be a positive integer.')
end
kfunargs = {pval}
setPoly = true
case 7 % mlpparams
if numel(pval)~=2
error('Bioinfo:svmtrain:BadMLPParams',...
'MLP_PARAMS must be a two element array.')
end
if ~isscalar(pval(1)) || ~isscalar(pval(2))
error('Bioinfo:svmtrain:MLPParamsNotScalar',...
'The parameters of the multi-layer perceptron kernel must be scalar.')
end
kfunargs = {pval(1),pval(2)}
setMLP = true
end
end
end
end
if setPoly &&~usePoly
warning('Bioinfo:svmtrain:PolyOrderNotPolyKernel',...
'You specified a polynomial order but not a polynomial kernel')
end
if setMLP &&~useMLP
warning('Bioinfo:svmtrain:MLPParamNotMLPKernel',...
'You specified MLP parameters but not an MLP kernel')
end
% plot the data if requested
if plotflag
[hAxis,hLines] = svmplotdata(training,g)
legend(hLines,cellstr(groupString))
end
% calculate kernel function
try
kx = feval(kfun,training,training,kfunargs{:})
% ensure function is symmetric
kx = (kx+kx')/2
catch
error('Bioinfo:svmtrain:UnknownKernelFunction',...
'Error calculating the kernel function:\n%s\n', lasterr)
end
% create Hessian
% add small constant eye to force stability
H =((g*g').*kx) + sqrt(eps(class(training)))*eye(n)
if useQuadprog
% The large scale solver cannot handle this type of problem, so turn it
% off.
qp_opts = optimset(qp_opts,'LargeScale','Off')
% X=QUADPROG(H,f,A,b,Aeq,beq,LB,UB,X0,opts)
alpha = quadprog(H,-ones(n,1),[],[],...
g',0,zeros(n,1),inf *ones(n,1),zeros(n,1),qp_opts)
% The support vectors are the non-zeros of alpha
svIndex = find(alpha >sqrt(eps))
sv = training(svIndex,:)
% calculate the parameters of the separating line from the support
% vectors.
alphaHat = g(svIndex).*alpha(svIndex)
% Calculate the bias by applying the indicator function to the support
% vector with largest alpha.
[maxAlpha,maxPos] = max(alpha)%#ok
bias = g(maxPos) - sum(alphaHat.*kx(svIndex,maxPos))
% an alternative method is to average the values over all support vectors
% bias = mean(g(sv)' - sum(alphaHat(:,ones(1,numSVs)).*kx(sv,sv)))
% An alternative way to calculate support vectors is to look for zeros of
% the Lagrangians (fifth output from QUADPROG).
%
% [alpha,fval,output,exitflag,t] = quadprog(H,-ones(n,1),[],[],...
% g',0,zeros(n,1),inf *ones(n,1),zeros(n,1),opts)
%
% sv = t.lower <sqrt(eps) &t.upper <sqrt(eps)
else % Least-Squares
% now build up compound matrix for solver
A = [0 g'g,H]
b = [0ones(size(g))]
x = A\b
% calculate the parameters of the separating line from the support
% vectors.
sv = training
bias = x(1)
alphaHat = g.*x(2:end)
end
svm_struct.SupportVectors = sv
svm_struct.Alpha = alphaHat
svm_struct.Bias = bias
svm_struct.KernelFunction = kfun
svm_struct.KernelFunctionArgs = kfunargs
svm_struct.GroupNames = groupnames
svm_struct.FigureHandles = []
if plotflag
hSV = svmplotsvs(hAxis,svm_struct)
svm_struct.FigureHandles = {hAxis,hLines,hSV}
end
支持向量机1 简介
支持向量机基本上是最好的有监督学习算法了。最开始接触SVM是去年暑假的时候,老师要求交《统计学习理论》的报告,那时去网上下了一份入门教程,里面讲的很通俗,当时只是大致了解了一些相关概念。这次斯坦福提供的学习材料,让我重新学习了一些SVM知识。我看很多正统的讲法都是从VC 维理论和结构风险最小原理出发,然后引出SVM什么的,还有些资料上来就讲分类超平面什么的。这份材料从前几节讲的logistic回归出发,引出了SVM,既揭示了模型间的联系,也让人觉得过渡更自然。
2 重新审视logistic回归
Logistic回归目的是从特征学习出一个0/1分类模型,而这个模型是将特性的线性组合作为自变量,由于自变量的取值范围是负无穷到正无穷。因此,使用logistic函数(或称作sigmoid函数)将自变量映射到(0,1)上,映射后的值被认为是属于y=1的概率。
形式化表示就是
假设函数
其中x是n维特征向量,函数g就是logistic函数。
的图像是
可以看到,将无穷映射到了(0,1)。
而假设函数就是特征属于y=1的概率。
当我们要判别一个新来的特征属于哪个类时,只需求 ,若大于0.5就是y=1的类,反之属于y=0类。
再审视一下 ,发现 只和 有关, >0,那么 ,g(z)只不过是用来映射,真实的类别决定权还在 。还有当 时, =1,反之 =0。如果我们只从 出发,希望模型达到的目标无非就是让训练数据中y=1的特征 ,而是y=0的特征 。Logistic回归就是要学习得到 ,使得正例的特征远大于0,负例的特征远小于0,强调在全部训练实例上达到这个目标。
图形化表示如下:
中间那条线是 ,logistic回顾强调所有点尽可能地远离中间那条线。学习出的结果也就中间那条线。考虑上面3个点A、B和C。从图中我们可以确定A是×类别的,然而C我们是不太确定的,B还算能够确定。这样我们可以得出结论,我们更应该关心靠近中间分割线的点,让他们尽可能地远离中间线,而不是在所有点上达到最优。因为那样的话,要使得一部分点靠近中间线来换取另外一部分点更加远离中间线。我想这就是支持向量机的思路和logistic回归的不同点,一个考虑局部(不关心已经确定远离的点),一个考虑全局(已经远离的点可能通过调整中间线使其能够更加远离)。这是我的个人直观理解。
3 形式化表示
我们这次使用的结果标签是y=-1,y=1,替换在logistic回归中使用的y=0和y=1。同时将 替换成w和b。以前的 ,其中认为 。现在我们替换 为b,后面替换 为 (即 )。这样,我们让 ,进一步 。也就是说除了y由y=0变为y=-1,只是标记不同外,与logistic回归的形式化表示没区别。再明确下假设函数
上一节提到过我们只需考虑 的正负问题,而不用关心g(z),因此我们这里将g(z)做一个简化,将其简单映射到y=-1和y=1上。映射关系如下:
4 函数间隔(functional margin)和几何间隔(geometric margin)
给定一个训练样本 ,x是特征,y是结果标签。i表示第i个样本。我们定义函数间隔如下:
可想而知,当 时,在我们的g(z)定义中, , 的值实际上就是 。反之亦然。为了使函数间隔最大(更大的信心确定该例是正例还是反例),当 时, 应该是个大正数,反之是个大负数。因此函数间隔代表了我们认为特征是正例还是反例的确信度。
继续考虑w和b,如果同时加大w和b,比如在 前面乘个系数比如2,那么所有点的函数间隔都会增大二倍,这个对求解问题来说不应该有影响,因为我们要求解的是 ,同时扩大w和b对结果是无影响的。这样,我们为了限制w和b,可能需要加入归一化条件,毕竟求解的目标是确定唯一一个w和b,而不是多组线性相关的向量。这个归一化一会再考虑。
刚刚我们定义的函数间隔是针对某一个样本的,现在我们定义全局样本上的函数间隔
说白了就是在训练样本上分类正例和负例确信度最小那个函数间隔。
接下来定义几何间隔,先看图
假设我们有了B点所在的 分割面。任何其他一点,比如A到该面的距离以 表示,假设B就是A在分割面上的投影。我们知道向量BA的方向是 (分割面的梯度),单位向量是 。A点是 ,所以B点是x= (利用初中的几何知识),带入 得,
进一步得到
实际上就是点到平面距离。
再换种更加优雅的写法:
当 时,不就是函数间隔吗?是的,前面提到的函数间隔归一化结果就是几何间隔。他们为什么会一样呢?因为函数间隔是我们定义的,在定义的时候就有几何间隔的色彩。同样,同时扩大w和b,w扩大几倍, 就扩大几倍,结果无影响。同样定义全局的几何间隔
5 最优间隔分类器(optimal margin classifier)
回想前面我们提到我们的目标是寻找一个超平面,使得离超平面比较近的点能有更大的间距。也就是我们不考虑所有的点都必须远离超平面,我们关心求得的超平面能够让所有点中离它最近的点具有最大间距。形象的说,我们将上面的图看作是一张纸,我们要找一条折线,按照这条折线折叠后,离折线最近的点的间距比其他折线都要大。形式化表示为:
这里用 =1规约w,使得 是几何间隔。
到此,我们已经将模型定义出来了。如果求得了w和b,那么来一个特征x,我们就能够分类了,称为最优间隔分类器。接下的问题就是如何求解w和b的问题了。
由于 不是凸函数,我们想先处理转化一下,考虑几何间隔和函数间隔的关系, ,我们改写一下上面的式子:
这时候其实我们求的最大值仍然是几何间隔,只不过此时的w不受 的约束了。然而这个时候目标函数仍然不是凸函数,没法直接代入优化软件里计算。我们还要改写。前面说到同时扩大w和b对结果没有影响,但我们最后要求的仍然是w和b的确定值,不是他们的一组倍数值,因此,我们需要对 做一些限制,以保证我们解是唯一的。这里为了简便我们取 。这样的意义是将全局的函数间隔定义为1,也即是将离超平面最近的点的距离定义为 。由于求 的最大值相当于求 的最小值,因此改写后结果为:
这下好了,只有线性约束了,而且是个典型的二次规划问题(目标函数是自变量的二次函数)。代入优化软件可解。
到这里发现,这个讲义虽然没有像其他讲义一样先画好图,画好分类超平面,在图上标示出间隔那么直观,但每一步推导有理有据,依靠思路的流畅性来推导出目标函数和约束。
接下来介绍的是手工求解的方法了,一种更优的求解方法。
6 拉格朗日对偶(Lagrange duality)
先抛开上面的二次规划问题,先来看看存在等式约束的极值问题求法,比如下面的最优化问题:
目标函数是f(w),下面是等式约束。通常解法是引入拉格朗日算子,这里使用 来表示算子,得到拉格朗日公式为
L是等式约束的个数。
然后分别对w和 求偏导,使得偏导数等于0,然后解出w和 。至于为什么引入拉格朗日算子可以求出极值,原因是f(w)的dw变化方向受其他不等式的约束,dw的变化方向与f(w)的梯度垂直时才能获得极值,而且在极值处,f(w)的梯度与其他等式梯度的线性组合平行,因此他们之间存在线性关系。(参考《最优化与KKT条件》)
然后我们探讨有不等式约束的极值问题求法,问题如下:
我们定义一般化的拉格朗日公式
这里的 和 都是拉格朗日算子。如果按这个公式求解,会出现问题,因为我们求解的是最小值,而这里的 已经不是0了,我们可以将 调整成很大的正值,来使最后的函数结果是负无穷。因此我们需要排除这种情况,我们定义下面的函数:
这里的P代表primal。假设 或者 ,那么我们总是可以调整 和 来使得 有最大值为正无穷。而只有g和h满足约束时, 为f(w)。这个函数的精妙之处在于 ,而且求极大值。
因此我们可以写作
这样我们原来要求的min f(w)可以转换成求 了。
我们使用 来表示 。如果直接求解,首先面对的是两个参数,而 也是不等式约束,然后再在w上求最小值。这个过程不容易做,那么怎么办呢?
我们先考虑另外一个问题
D的意思是对偶, 将问题转化为先求拉格朗日关于w的最小值,将 和 看作是固定值。之后在 求最大值的话:
这个问题是原问题的对偶问题,相对于原问题只是更换了min和max的顺序,而一般更换顺序的结果是Max Min(X) <= MinMax(X)。然而在这里两者相等。用 来表示对偶问题如下:
下面解释在什么条件下两者会等价。假设f和g都是凸函数,h是仿射的(affine, )。并且存在w使得对于所有的i, 。在这种假设下,一定存在 使得 是原问题的解, 是对偶问题的解。还有 另外, 满足库恩-塔克条件(Karush-Kuhn-Tucker, KKT condition),该条件如下:
所以如果 满足了库恩-塔克条件,那么他们就是原问题和对偶问题的解。让我们再次审视公式(5),这个条件称作是KKT dual complementarity条件。这个条件隐含了如果 ,那么 。也就是说, 时,w处于可行域的边界上,这时才是起作用的约束。而其他位于可行域内部( 的)点都是不起作用的约束,其 。这个KKT双重补足条件会用来解释支持向量和SMO的收敛测试。
这部分内容思路比较凌乱,还需要先研究下《非线性规划》中的约束极值问题,再回头看看。KKT的总体思想是将极值会在可行域边界上取得,也就是不等式为0或等式约束里取得,而最优下降方向一般是这些等式的线性组合,其中每个元素要么是不等式为0的约束,要么是等式约束。对于在可行域边界内的点,对最优解不起作用,因此前面的系数为0。
7 最优间隔分类器(optimal margin classifier)
重新回到SVM的优化问题:
我们将约束条件改写为:
从KKT条件得知只有函数间隔是1(离超平面最近的点)的线性约束式前面的系数 ,也就是说这些约束式 ,对于其他的不在线上的点( ),极值不会在他们所在的范围内取得,因此前面的系数 .注意每一个约束式实际就是一个训练样本。
看下面的图:
实线是最大间隔超平面,假设×号的是正例,圆圈的是负例。在虚线上的点就是函数间隔是1的点,那么他们前面的系数 ,其他点都是 。这三个点称作支持向量。构造拉格朗日函数如下:
注意到这里只有 没有 是因为原问题中没有等式约束,只有不等式约束。
下面我们按照对偶问题的求解步骤来一步步进行,
首先求解 的最小值,对于固定的 , 的最小值只与w和b有关。对w和b分别求偏导数。
并得到
将上式带回到拉格朗日函数中得到,此时得到的是该函数的最小值(目标函数是凸函数)
代入后,化简过程如下:
最后得到
由于最后一项是0,因此简化为
这里我们将向量内积 表示为
此时的拉格朗日函数只包含了变量 。然而我们求出了 才能得到w和b。
接着是极大化的过程 ,
前面提到过对偶问题和原问题满足的几个条件,首先由于目标函数和线性约束都是凸函数,而且这里不存在等式约束h。存在w使得对于所有的i, 。因此,一定存在 使得 是原问题的解, 是对偶问题的解。在这里,求 就是求 了。
如果求出了 ,根据 即可求出w(也是 ,原问题的解)。然后
即可求出b。即离超平面最近的正的函数间隔要等于离超平面最近的负的函数间隔。
关于上面的对偶问题如何求解,将留给下一篇中的SMO算法来阐明。
这里考虑另外一个问题,由于前面求解中得到
我们通篇考虑问题的出发点是 ,根据求解得到的 ,我们代入前式得到
也就是说,以前新来的要分类的样本首先根据w和b做一次线性运算,然后看求的结果是大于0还是小于0,来判断正例还是负例。现在有了 ,我们不需要求出w,只需将新来的样本和训练数据中的所有样本做内积和即可。那有人会说,与前面所有的样本都做运算是不是太耗时了?其实不然,我们从KKT条件中得到,只有支持向量的 ,其他情况 。因此,我们只需求新来的样本和支持向量的内积,然后运算即可。这种写法为下面要提到的核函数(kernel)做了很好的铺垫。这是上篇,先写这么多了。
7 核函数(Kernels)
考虑我们最初在“线性回归”中提出的问题,特征是房子的面积x,这里的x是实数,结果y是房子的价格。假设我们从样本点的分布中看到x和y符合3次曲线,那么我们希望使用x的三次多项式来逼近这些样本点。那么首先需要将特征x扩展到三维 ,然后寻找特征和结果之间的模型。我们将这种特征变换称作特征映射(feature mapping)。映射函数称作 ,在这个例子中
我们希望将得到的特征映射后的特征应用于SVM分类,而不是最初的特征。这样,我们需要将前面 公式中的内积从 ,映射到 。
至于为什么需要映射后的特征而不是最初的特征来参与计算,上面提到的(为了更好地拟合)是其中一个原因,另外的一个重要原因是样例可能存在线性不可分的情况,而将特征映射到高维空间后,往往就可分了。(在《数据挖掘导论》Pang-Ning Tan等人著的《支持向量机》那一章有个很好的例子说明)
将核函数形式化定义,如果原始特征内积是 ,映射后为 ,那么定义核函数(Kernel)为
到这里,我们可以得出结论,如果要实现该节开头的效果,只需先计算 ,然后计算 即可,然而这种计算方式是非常低效的。比如最初的特征是n维的,我们将其映射到 维,然后再计算,这样需要 的时间。那么我们能不能想办法减少计算时间呢?
先看一个例子,假设x和z都是n维的,
展开后,得
这个时候发现我们可以只计算原始特征x和z内积的平方(时间复杂度是O(n)),就等价与计算映射后特征的内积。也就是说我们不需要花 时间了。
现在看一下映射函数(n=3时),根据上面的公式,得到
也就是说核函数 只能在选择这样的 作为映射函数时才能够等价于映射后特征的内积。
再看一个核函数
对应的映射函数(n=3时)是
更一般地,核函数 对应的映射后特征维度为 。(这个我一直没有理解)。
由于计算的是内积,我们可以想到IR中的余弦相似度,如果x和z向量夹角越小,那么核函数值越大,反之,越小。因此,核函数值是 和 的相似度。
再看另外一个核函数
这时,如果x和z很相近( ),那么核函数值为1,如果x和z相差很大( ),那么核函数值约等于0。由于这个函数类似于高斯分布,因此称为高斯核函数,也叫做径向基函数(Radial Basis Function 简称RBF)。它能够把原始特征映射到无穷维。
既然高斯核函数能够比较x和z的相似度,并映射到0到1,回想logistic回归,sigmoid函数可以,因此还有sigmoid核函数等等。
下面有张图说明在低维线性不可分时,映射到高维后就可分了,使用高斯核函数。
来自Eric Xing的slides
注意,使用核函数后,怎么分类新来的样本呢?线性的时候我们使用SVM学习出w和b,新来样本x的话,我们使用 来判断,如果值大于等于1,那么是正类,小于等于是负类。在两者之间,认为无法确定。如果使用了核函数后, 就变成了 ,是否先要找到 ,然后再预测?答案肯定不是了,找 很麻烦,回想我们之前说过的
只需将 替换成 ,然后值的判断同上。
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)