tic % 计时
%% 清空环境导入数据
clear
clc
close all
load aa %
%% PSO-SVR
train_input=aa(1:120,1:4); %输入训练
train_output=aa(1:120,5); %输出训练
test_input=aa(121:150,1:4); %输入测试
test_input=test_input';
test_output=aa(121:150,5); %输出测试
% test_output=rand(100,1);
test_output=test_output';
output_test1=aa(121:150,5);
[input_train,rule1]=mapminmax(train_input',0,1);
[output_train,rule2]=mapminmax(train_output',0,1);
input_test=mapminmax('apply',test_input,rule1);
output_test=mapminmax('apply',test_output,rule2);
%% 参数设置
% 粒子位置和速度更新参数
c1=1.5;
c2=1.5;
maxgen=5; % 迭代次数
sizepop=5; % 种群规模（粒子数）
% 目标函数信息
dim=2; % 待优化参数个数
xub=[100,100]; % 参数取值上界
xlb=[0.01,0.01]; % 参数取值下界
vub=[5,5]; % 速度上界
vlb=[-5,-5]; % 速度下界
%% 初始化
% 粒子位置初始化
xRange=repmat((xub-xlb),[sizepop,1]);
xLower=repmat(xlb,[sizepop,1]);
pop=rand(sizepop,dim).*xRange+xLower;
% 粒子速度初始化
vRange=repmat((vub-vlb),[sizepop,1]);
vLower=repmat(vlb,[sizepop,1]);
V=rand(sizepop,dim).*vRange+vLower;
% 计算初始化的目标函数值（适应度函数值）
fitness=ones(sizepop,1);
for k=1:sizepop
fitness(k)=fobj(pop(k,:),input_train,output_train,input_test,output_test);
end
%% 寻找初始极值
[bestfitness,bestindex]=min(fitness);
zbest=pop(bestindex,:); % 全局最优解的位置
gbest=pop; % 个体位置
fitnessgbest=fitness; % 个体适应度值
fitnesszbest=bestfitness; % 全局最优适应度值
%% 迭代寻优
yy=ones(sizepop,1); % 用于保存每次迭代的最优目标函数值
for k=1:maxgen
% 粒子位置和速度更新
for m=1:sizepop
% 粒子速度更新
sol=V(m,:)+c1*rand*(gbest(m,:)-pop(m,:))+c2*rand*(zbest-pop(m,:));
% 确保粒子速度取值范围不越界
ind=find(sol<vlb);
sol(ind)=vlb(ind);
ind=find(sol>vub);
sol(ind)=vub(ind);
V(m,:)=sol;
% 粒子位置更新
sol=pop(m,:)+0.5*V(m,:);
% 确保粒子位置取值范围不越界
ind=find(sol<xlb);
sol(ind)=xlb(ind);
ind=find(sol>xub);
sol(ind)=xub(ind);
pop(m,:)=sol;
% 更新粒子适应度值
fitness(m)=fobj(pop(m,:),input_train,output_train,input_test,output_test);
end
% 个体极值及位置和群体极值及位置更新
for m=1:sizepop
% 个体极值及其位置更新
if fitness(m)<fitnessgbest(m)
gbest(m,:)=pop(m,:);
fitnessgbest(m)=fitness(m);
end
% 群体极值及其位置更新
if fitness(m)<fitnesszbest
zbest=pop(m,:);
fitnesszbest=fitness(m);
end
end
% 记录每代最优值
yy(k)=fitnesszbest;
end
%% 打印参数选择结果
bestobjfun=fitnesszbest;
bestc=zbest(1);
bestg=zbest(2);
disp('打印参数选择结果');
str=sprintf('Best c = %g，Best g = %g',bestc,bestg);
disp(str)
%% 利用回归预测分析最佳的参数进行SVM网络训练
cmd_cs_svr=['-s 3 -t 2',' -c ',num2str(bestc),' -g ',num2str(bestg)];
model_cs_svr=libsvmtrain(output_train',input_train',cmd_cs_svr); % SVM模型训练
%% SVM网络回归预测
[output_test_pre,acc,~]=libsvmpredict(output_test',input_test',model_cs_svr); % SVM模型预测及其精度
test_pre=mapminmax('reverse',output_test_pre',rule2);
test_pre = test_pre';
err_pre=output_test1-test_pre;
figure('Name','测试数据残差图')
set(gcf,'unit','centimeters','position',[0.5,5,30,5])
plot(err_pre,'*-');
figure('Name','原始-预测图')
plot(test_pre,'*r-');
hold on;
plot(output_test1,'bo-');
legend('预测','原始')
% set(gcf,'unit','centimeters','position',[0.5,13,30,5]) %set指令设置图片的大小位置 0.5位距离屏幕左侧0.5厘米，13为距离屏幕下方13厘米，30为图片长度，5位图片宽度
result=[output_test1,test_pre]
MAE1=mymae(output_test1,test_pre)
MSE1=mymse(output_test1,test_pre)
MAPE1=mymape(output_test1,test_pre)
%% 显示程序运行时间
toc
