%%%%%基于模拟退火粒子群优化函数——轮盘赌策略迭代 %%%%惯性权重二次缩减
clc clear %%%% %导入训练数据和测试数据 train_input=xlsread('Sample Data2',['AA2:AG111']); train_out=xlsread('Sample Data2',['Z2:Z111']); test_input=xlsread('Sample Data2',['AA112:AH141']); test_out=xlsread('Sample Data2',['Z112:Z141']); %%归一化 [train_data ,pstrain0] = mapminmax(train_input',0,1); [test_data] = mapminmax('apply',test_input',pstrain0); [train_result,pstrain1] = mapminmax(train_out',0,1); [test_result] = mapminmax('apply',test_out',pstrain1); %%取转置 train_data = train_data'; train_result=train_result'; test_data = test_data'; %%参数初始化 N=80; %%%种群规模 c1=1.5; c2=2.0; c1f=0.5;c1i=2.5; %c1f、c1i表示对应的最终和初始值 c2f=2.5;c2i=0.5; %%%%认知因子 wmax=1.2; wmin=0.8; %%%%惯性权重的取值范围 xmax=[1000 100]; xmin=[0.1 0.01];%%%%%%位置范围 M=400; %%%%%迭代次数 %MaxC=10; %%%%%混沌搜索步数 D=2; %%%%%问题维数 vmax=[500 50]; vmin=[-500 -50]; %%%%速度范围 lamda=0.8; %%%%退火常数 eps=10^(-8); %机器的浮点运算误差限,当某个数小于它时则记为0 %%定义LSSVM相关参数 type='function estimation'; %回归 %type=’c’,表示classification 分类,type=’f’,表示回归。 kernel='RBF_kernel'; %LSSVM的核函数选取为径向基内核 proprecess='proprecess'; %proprecess表明数据已归一化处理,original表明数据没有进行归一化处理 format long; %% %------初始化种群的个体------------
for i=1:N
for j=1:D
x(i,j)=(xmax(1,j)-xmin(1,j))*rand(1,1)+xmin(1,j); %随机初始化位置
v(i,j)=vmax(1,j)*rands(1,1); %随机初始化速度
end end %------先计算各个粒子的适应度,并初始化个体适应度P(i)和群体最优适应度Pg----------------------
for i=1:N %%%%%计算各个粒子的适应度p(i) gam=x(i,1);sig2=x(i,2); model=initlssvm(train_data,train_result,type,gam,sig2,kernel,proprecess); model=trainlssvm(model); %求出训练集和测试集的预测值 [train_predict_y,zt,model]=simlssvm(model,train_data); [test_predict_y,zt,model]=simlssvm(model,test_data); %训练集和测试集的预测数据反归一化 train_predict=mapminmax('reverse',train_predict_y',pstrain1);%训练集预测值 test_predict=mapminmax('reverse',test_predict_y',pstrain1); %测试集预测值 %%%%预测数据四舍五入保留一位有效数字 train_predict=roundn(train_predict,-1); test_predict=roundn(test_predict,-1); %计算均方差 result=[train_out;test_out]; predict=[train_predict';test_predict']; allmse=(sum(((result-predict)./result).^2))/length(result); fitness(i)=allmse; %以所有数据的预测值计算的均方差为适应度值
p(i)=fitness(i);
y(i,:)=x(i,:);
end
pg = x(N,:); %Pg为全局最优
for i=1:(N-1) %%%%%计算群体粒子的最优适应度pg gam=pg(1);sig2=pg(2); model=initlssvm(train_data,train_result,type,gam,sig2,kernel,proprecess); model=trainlssvm(model); %求出训练集和测试集的预测值 [train_predict_y,zt,model]=simlssvm(model,train_data); [test_predict_y,zt,model]=simlssvm(model,test_data); %训练集和测试集的预测数据反归一化 train_predict=mapminmax('reverse',train_predict_y',pstrain1);%训练集预测值 test_predict=mapminmax('reverse',test_predict_y',pstrain1); %测试集预测值 %%%%预测数据四舍五入保留一位有效数字 train_predict=roundn(train_predict,-1); test_predict=roundn(test_predict,-1); %计算均方差 result=[train_out;test_out]; predict=[train_predict';test_predict']; allmse=(sum(((result-predict)./result).^2))/length(result); global_fitness=allmse; %以所有数据的预测值计算的均方差为适应度值
if p(i)
pg=x(i,:);
end
end
%% %------进入主要循环,按照公式依次迭代------------
T=global_fitness/log(5); %%%设置初始温度
for t=1:M groupFit=global_fitness %%%%群体最优
%%%根据轮盘赌策略公式确定当前温度下各粒子的适配值 for i=1:N Tfit(i) = exp( - (p(i) - groupFit)/T); end SumTfit = sum(Tfit); Tfit = Tfit/SumTfit; pBet = rand(); for i=1:N ComFit(i)=sum(Tfit(1:i)); if pBetvmax(1,j) v(i,j)=vmax(1,j); end if v(i,j) v(i,j)=vmin(1,j); end if x(i,j)>xmax(1,j) x(i,j)=xmax(1,j); end if x(i,j) x(i,j)=xmin(1,j); end end %%%%%%%计算更新后粒子的适应度 gam=x(i,1);sig2=x(i,2); model=initlssvm(train_data,train_result,type,gam,sig2,kernel,proprecess); model=trainlssvm(model); %求出训练集和测试集的预测值 [train_predict_y,zt,model]=simlssvm(model,train_data); [test_predict_y,zt,model]=simlssvm(model,test_data); %训练集和测试集的预测数据反归一化 train_predict=mapminmax('reverse',train_predict_y',pstrain1);%训练集预测值 test_predict=mapminmax('reverse',test_predict_y',pstrain1); %测试集预测值 %%%%预测数据四舍五入保留一位有效数字 train_predict=roundn(train_predict,-1); test_predict=roundn(test_predict,-1); %计算均方差 result=[train_out;test_out]; predict=[train_predict';test_predict']; allmse=(sum(((result-predict)./result).^2))/length(result); fitness(i)=allmse; %以所有数据的预测值计算的均方差为适应度值
local_fit(i)=fitness(i); if local_fit(i)
p(i)=local_fit(i); y(i,:)=x(i,:); end if p(i) pg=y(i,:); end end T=T*lamda;%%%%%%降温 yy(t)=groupFit; end xm=pg' fv=groupFit hold on plot(yy,'g') hold off %title(['适应度曲线 ' '终止代数=' num2str(M)]); %xlabel('进化代数');ylabel('适应度');