数学建模常用代码
SVM分类器
1.命令函数部分：
clear;%清屏
clc;
X =load('data.txt');
n = length(X);%总样本数量
y = X(:,4);%类别标志
X = X(:,1:3);
TOL = 0.0001;%精度要求
C = 1;%参数，对损失函数的权重
b = 0;%初始设置截距b
Wold = 0;%未更新a时的W(a)
Wnew = 0;%更新a后的W(a)
for i = 1 : 50%设置类别标志为1或者-1y(i) = -1;
end
a = zeros(n,1);%参数a
for i = 1 : n%随机初始化a,a属于[0,C]a(i) = 0.2;
end%为简化计算，减少重复计算进行的计算
K = ones(n,n);
for i = 1 :n%求出K矩阵，便于之后的计算for j = 1 : nK(i,j) = k(X(i,:),X(j,:));end
end
sum = zeros(n,1);%中间变量，便于之后的计算，sum(k)=sigma a(i)*y(i)*K(k,i);
for k = 1 : nfor i = 1 : nsum(k) = sum(k) + a(i) * y(i) * K(i,k);end
endwhile 1%迭代过程%启发式选点
n1 = 1;%初始化，n1,n2代表选择的2个点
n2 = 2;
%n1按照第一个违反KKT条件的点选择
while n1 <= nif y(n1) * (sum(n1) + b) == 1 && a(n1) >= C && a(n1) <= 0break;endif y(n1) * (sum(n1) + b) > 1 && a(n1) ~= 0break;endif y(n1) * (sum(n1) + b) < 1 && a(n1) ~=Cbreak;endn1 = n1 + 1;
end
%n2按照最大化|E1-E2|的原则选取
E1 = 0;
E2 = 0;
maxDiff = 0;%假设的最大误差
E1 = sum(n1) + b - y(n1);%n1的误差
for i = 1 : ntempSum = sum(i) + b - y(i);if abs(E1 - tempSum)> maxDiffmaxDiff = abs(E1 - tempSum);n2 = i;E2 = tempSum;end
end%以下进行更新
a1old = a(n1);
a2old = a(n2);
KK = K(n1,n1) + K(n2,n2) - 2*K(n1,n2);
a2new = a2old + y(n2) *(E1 - E2) / KK;%计算新的a2
%a2必须满足约束条件
S = y(n1) * y(n2);
if S == -1U = max(0,a2old - a1old);V = min(C,C - a1old + a2old);
elseU = max(0,a1old + a2old - C);V = min(C,a1old + a2old);
end
if a2new > Va2new = V;
end
if a2new < Ua2new = U;
end
a1new = a1old + S * (a2old - a2new);%计算新的a1
a(n1) = a1new;%更新a
a(n2) = a2new;%更新部分值
sum = zeros(n,1);
for k = 1 : nfor i = 1 : nsum(k) = sum(k) + a(i) * y(i) * K(i,k);end
end
Wold = Wnew;
Wnew = 0;%更新a后的W(a)
tempSum = 0;%临时变量
for i = 1 : nfor j = 1 : ntempSum= tempSum + y(i )*y(j)*a(i)*a(j)*K(i,j);endWnew= Wnew+ a(i);
end
Wnew= Wnew - 0.5 * tempSum;
%以下更新b：通过找到某一个支持向量来计算
support = 1;%支持向量坐标初始化
while abs(a(support))< 1e-4 && support <= nsupport = support + 1;
end
b = 1 / y(support) - sum(support);
%判断停止条件
if abs(Wnew/ Wold - 1 ) <= TOLbreak;
end
end
%输出结果：包括原分类，辨别函数计算结果，svm分类结果
for i = 1 : nfprintf('第%d点:原标号 ',i);if i <= 50fprintf('-1');elsefprintf(' 1');endfprintf(' 判别函数值%f 分类结果',sum(i) + b);if abs(sum(i) + b - 1) < 0.5fprintf('1\n');else if abs(sum(i) + b + 1) < 0.5fprintf('-1\n');elsefprintf('归类错误\n');endend
end2.名为f的功能函数部分:
function y = k(x1,x2)y = exp(-0.5*norm(x1 - x2).^2);
endK-means算法代码
function [Idx, Center] = K_means(X, xstart)
% K-means聚类
% Idx是数据点属于哪个类的标记，Center是每个类的中心位置
% X是全部二维数据点，xstart是类的初始中心位置len = length(X); %X中的数据点个数
Idx = zeros(len, 1); %每个数据点的Id，即属于哪个类C1 = xstart(1,:); %第1类的中心位置
C2 = xstart(2,:); %第2类的中心位置
C3 = xstart(3,:); %第3类的中心位置for i_for = 1:100%为避免循环运行时间过长，通常设置一个循环次数%或相邻两次聚类中心位置调整幅度小于某阈值则停止%更新数据点属于哪个类for i = 1:lenx_temp = X(i,:); %提取出单个数据点d1 = norm(x_temp - C1); %与第1个类的距离d2 = norm(x_temp - C2); %与第2个类的距离d3 = norm(x_temp - C3); %与第3个类的距离d = [d1;d2;d3];[~, id] = min(d); %离哪个类最近则属于那个类Idx(i) = id;end%更新类的中心位置L1 = X(Idx == 1,:); %属于第1类的数据点L2 = X(Idx == 2,:); %属于第2类的数据点L3 = X(Idx == 3,:); %属于第3类的数据点C1 = mean(L1); %更新第1类的中心位置C2 = mean(L2); %更新第2类的中心位置C3 = mean(L3); %更新第3类的中心位置
endCenter = [C1; C2; C3]; %类的中心位置%演示数据
%% 1 random sample
%随机生成三组数据
a = rand(30,2) * 2;
b = rand(30,2) * 5;
c = rand(30,2) * 10;
figure(1);
subplot(2,2,1);
plot(a(:,1), a(:,2), 'r.'); hold on
plot(b(:,1), b(:,2), 'g*');
plot(c(:,1), c(:,2), 'bx'); hold off
grid on;
title('raw data');%% 2 K-means cluster
X = [a; b; c]; %需要聚类的数据点
xstart = [2 2; 5 5; 8 8]; %初始聚类中心
subplot(2,2,2);
plot(X(:,1), X(:,2), 'kx'); hold on
plot(xstart(:,1), xstart(:,2), 'r*'); hold off
grid on;
title('raw data center');[Idx, Center] = K_means(X, xstart);
subplot(2,2,4);
plot(X(Idx==1,1), X(Idx==1,2), 'kx'); hold on
plot(X(Idx==2,1), X(Idx==2,2), 'gx');
plot(X(Idx==3,1), X(Idx==3,2), 'bx');
plot(Center(:,1), Center(:,2), 'r*'); hold off
grid on;
title('K-means cluster result');disp('xstart = ');
disp(xstart);
disp('Center = ');
disp(Center);