function_svmtrain

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, including % the support vectors, 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(.,RBF_SIGMA,SIGMA) allows you to specify the scaling % factor, sigma, in the radial basis function kernel. % % 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 2 error(Bioinfo:svmtrain:TooManyGroups,. SVMTRAIN only supports classification into two groups.nGROUP contains %d different groups.,ngroups) end % convert to 1, -1. groupIndex = 1 - (2* (groupIndex-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,. boxconstraint,rbf_sigma,autoscale, smo_opts; for j=1:2:numoptargs pname = optargsj; pval = optargsj+1; k = find(strncmpi(pname, okargs,numel(pname); 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); 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; useSigma = true; case 5 kfun = poly_kernel; usePoly = true; case 6 kfun = mlp_kernel; useMLP = true; otherwise error(Bioinfo:svmtrain:UnknownKernelFunction,. Unknown Kernel Function %s.,pval); end elseif isa (pval, function_handle) kfun = pval; else error(Bioinfo:svmtrain:BadKernelFunction,. The kernel function input does not appear to be a function handlenor valid function name.) end case 2 % method if strncmpi(pval,qp,2) if isempty(which(quadprog) warning(Bioinfo:svmtrain:NoOptim,. The Optimization Toolbox is required to use the quadratic programming method.) optimMethod = SMO; else optimMethod = QP; end elseif strncmpi(pval,ls,numel(pval) optimMethod = LS; elseif strncmpi(pval, smo, numel(pval) optimMethod = SMO; else error(Bioinfo:svmtrain:UnknownMethod,. Unknown method option %s. Valid methods are SMO, QP and LS,pval); end case 3 % display plotflag = opttf(pval,okargsk,mfilename); if plotflag = true if size(training,2) = 2 plotflag = true; else plotflag = false; 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 = optimset(qp_opts,pval); elseif iscell(pval) qp_opts = optimset(qp_opts,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 0 warning(Bioinfo:svmtrain:MLPBiasNotNegative,. The bias for MLP kernel should be negative.) end kfunargs = pval(1),pval(2); setMLP = true; case 8 % box constraint: it can be a positive numeric scalar % or a numeric vector of the same length as there are % data points if isscalar(pval) n2 = length(find(groupIndex=-1); c1 = 0.5 * pval * nPoints / n1; c2 = 0.5 * pval * nPoints / n2; boxconstraint(groupIndex=1) = c1; boxconstraint(groupIndex=-1) = c2; elseif isvector(pval) end boxconstraint = pval; else error(Bioinfo:svmtrain:BoxConstraintNotScalar,. The box constraint must be a numeric scalar or vector 0.); end % If boxconstraint = Inf then convergence will not % happen so fix the value to 1/sqrt(eps). boxconstraint = min(boxconstraint,repmat(1/sqrt(eps(class(boxconstraint),. size(boxconstraint); case 9 % rbf sigma if isscalar(pval) setSigma = true; else error(Bioinfo:svmtrain:RBFSigmaNotScalar,. Sigma must be a numeric scalar.); end case 10 % autoscale autoScale = opttf(pval,okargsk,mfilename); case 11 % smo_opts smo_opts = pval; end end end end if setPoly end if setMLP end if setSigma end % plot the data if requested if plotflag hAxis,hLines = svmplotdata(training,groupIndex); legend(hLines,cellstr(groupString); end % autoscale data if required, we cant use the zscore function here, % because we need the shift and scale values. scaleData = ; if autoScale scaleData.shift = - mean(training); stdVals = std(training); scaleData.scaleFactor = 1./stdVals; % leave zero-variance data unscaled: scaleData.scaleFactor(isfinite(scaleData.scaleFactor) = 1; % shift and scale columns of data matrix: for c = 1:size(training, 2) training(:,c) = scaleData.scaleFactor(c) * . (training(:,c) + scaleData.shift(c); end end if strcmpi(optimMethod, SMO) % if we have a kernel that takes extra arguments we must define a new % kernel function handle to be passed to seqminopt if isempty(kfunargs) tmp_kfun = (x,y) feval(kfun, x,y, kfunargs:); else tmp_kfun = kfun; end alpha bias = seqminopt(training, groupIndex, . boxconstraint, tmp_kfun, smo_opts); svIndex = find(alpha sqrt(eps); sv = training(svIndex,:); alphaHat = groupIndex(svIndex).*alpha(svIndex); else % QP and LS both need the kernel matrix: % calculate kernel function and add additional term required % for two-norm soft margin try kx = feval(kfun,training,training,kfunargs:); % ensure function is symmetric kx = (kx+kx)/2 + diag(1./boxconstraint); catch theException error(Bioinfo:svmtrain:KernelFunctionError,. Error calculating the kernel function:n%sn, theException.message); end % create Hessian H =(groupIndex * groupIndex).*kx); if strcmpi(optimMethod, QP) % X=QUADPROG(H,f,A,b,Aeq,beq,LB,UB,X0,opts) alpha, fval, exitflag, . output, lambda = quadprog(H,-ones(nPoints,1),. groupIndex,0,zeros(nPoints,1), . Inf *ones(nPoints,1), ones(nPoints,1), . qp_opts); %#ok if exitflag sqrt(eps); sv = training(svIndex,:); % calculate the parameters of the separating line from the support % vectors. alphaHat = groupIndex(svIndex).*alpha(svIndex); % Calculate the bias by applying the indicator function to the support % vector with largest alpha. maxAlpha,maxPos = max(alpha); bias = groupIndex(maxPos) - sum(alphaHat.*kx(svIndex,maxPos); % an alternative method is to average the values over all support vectors % bias = mean(groupIndex(sv) - sum(alphaHat(:,ones(1,numSVs).*kx(sv,sv); % An alternative way to calculate support vectors is to look for zeros of % th