使用pso优化BP神经网络算法
|%% 该代码为基于PSO和BP网络的预测
%% 清空环境
clc
clear

%读取数据
load data input_train input_test output_train output_test

%节点个数
inputnum=8;
hiddennum=17;
outputnum=8;

%训练数据和预测数据
input_train=input_train’;
input_test=input_test’;
output_train=output_train’;
output_test=output_test’;
%样本输入输出数据归一化
[inputn,inputps]=mapminmax(input_train);
[outputn,outputps]=mapminmax(output_train);

%构建网络
net=newff(inputn,outputn,hiddennum);

%%新添加
%初始化编码长度，以免后面出现超维错误gl
numsum=inputnumhiddennum+hiddennum+hiddennumoutputnum+outputnum; %编码长度
% 参数初始化
%粒子群算法中的两个参数
c1 = 1.49445;
c2 = 1.49445;

maxgen=2; % 进化次数 100
sizepop=30; %种群规模
wmax=0.9;
wmin=0.4;

Vmax=1;
Vmin=-1;
popmax=5;
popmin=-5;

%% 产生初始粒子和速度
for i=1:sizepop
%随机产生一个种群
% pop(i,:)=normrnd(0,1,numsum,20);
% vov(i,:)=normrnd(0,1,numsum,20);
pop(i,:)=5*rands(1,numsum); %初始种群numsum
vov(i,:)=rands(1,numsum); %初始化速度
%计算适应度

end
x=pop(i,:);
fitness(i)=fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn); %染色体的适应度
% 个体极值和群体极值
[bestfitness bestindex]=min(fitness);
zbest=pop(bestindex,:); %全局最佳
gbest=pop; %个体最佳
fitnessgbest=fitness; %个体最佳适应度值
fitnesszbest=bestfitness; %全局最佳适应度值

%% 迭代寻优
for i=1:maxgen
i;
%粒子位置和速度更新
for j=1:sizepop
w=wmax-(wmax-wmin)j/maxgen;
%速度更新
vov(j,:) = wvov(j,:) + c1rand(gbest(j,:) – pop(j,:)) + c2rand(zbest – pop(j,:));
vov(j,find(vov(j,:)>Vmax))=Vmax;
vov(j,find(vov(j,:)<Vmin))=Vmin;

    %种群更新
    pop(j,:)=pop(j,:)+0.5*vov(j,:);
    pop(j,find(pop(j,:)>popmax))=popmax;
    pop(j,find(pop(j,:)<popmin))=popmin;
    
    %引入变异算子，重新初始化粒子

% if rand>0.9
% k=ceil(21*rand);
% pop(j,k)=rand;
% end

    pos=unidrnd(297);  %gl 此处例程是pos=unidrnd(21); ;把21换成了numsum。
    if rand>0.95
        pop(j,pos)=5*rands(1,1);
    end
    %新粒子适应度值
    fitness(j)=fun(pop(j,:),inputnum,hiddennum,outputnum,net,inputn,outputn);
end
%%个体极值和群体极值更新
for j=1:sizepop
%个体最优更新
if fitness(j) < fitnessgbest(j)
    gbest(j,:) = pop(j,:);
    fitnessgbest(j) = fitness(j);
end

%群体最优更新 
if fitness(j) < fitnesszbest
    zbest = pop(j,:);
    fitnesszbest = fitness(j);
end

end
%%每代最优值记录到yy数组中
yy(i)=fitnesszbest;    

end

%% 结果分析
plot(yy)
title([‘适应度曲线 ’ ‘终止代数＝’ num2str(maxgen)],‘fontsize’,12);
xlabel(‘进化代数’,‘fontsize’,12);ylabel(‘适应度’,‘fontsize’,12);

x=zbest;
%% 把最优初始阀值权值赋予网络预测
% %用遗传算法优化的BP网络进行值预测
w1=x(1:inputnumhiddennum);
B1=x(inputnumhiddennum+1:inputnumhiddennum+hiddennum);
w2=x(inputnumhiddennum+hiddennum+1:inputnumhiddennum+hiddennum+hiddennumoutputnum);
B2=x(inputnumhiddennum+hiddennum+hiddennumoutputnum+1:inputnumhiddennum+hiddennum+hiddennumoutputnum+outputnum);

net.iw{1,1}=reshape(w1,hiddennum,inputnum);
net.lw{2,1}=reshape(w2,outputnum,hiddennum);
net.b{1}=reshape(B1,hiddennum,1);
net.b{2}=reshape(B2,outputnum,1);
% net.b{2}=B2;

%% BP网络训练
%网络进化参数
net.trainParam.epochs=100;
net.trainParam.lr=0.1;
net.trainParam.goal=0.00001;

%网络训练
[net,tr]=train(net,inputn,outputn);

%% BP网络预测
%数据归一化
inputn_test=mapminmax(‘apply’,input_test,inputps);
an=sim(net,inputn_test);
test_simu=mapminmax(‘reverse’,an,outputps);
error=test_simu-output_test;
figure(2)
plot(error)
title(‘仿真预测误差’,‘fontsize’,12);
xlabel(‘仿真次数’,‘fontsize’,12);ylabel(‘误差百分值’,‘fontsize’,12);
子程序：function
| function error = fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn)
% [m,n]=size(x);
%该函数用来计算适应度值
%x input 个体
%inputnum input 输入层节点数
%outputnum input 隐含层节点数
%net input 网络
%inputn input 训练输入数据
%outputn input 训练输出数据

%error output 个体适应度值

%提取 BP神经网络初始权值和阈值，x为个体

w1=x(1:inputnumhiddennum);
B1=x(inputnumhiddennum+1:inputnumhiddennum+hiddennum);
w2=x(inputnumhiddennum+hiddennum+1:inputnumhiddennum+hiddennum+hiddennumoutputnum);
B2=x(inputnumhiddennum+hiddennum+hiddennumoutputnum+1:inputnumhiddennum+hiddennum+hiddennumoutputnum+outputnum);

%网络权值赋值
net.iw{1,1}=reshape(w1,hiddennum,inputnum);
net.lw{2,1}=reshape(w2,outputnum,hiddennum);
net.b{1}=reshape(B1,hiddennum,1);
net.b{2}=reshape(B2,outputnum,1);
% net.b{2}=B2;
%BP神经网络构建
net=newff(inputn,outputn,hiddennum);
net.trainParam.epochs=100;
net.trainParam.lr=0.1;
net.trainParam.goal=0.00001;
net.trainParam.show=100;
net.trainParam.showWindow=0;

%BP神经网络训练
net=train(net,inputn,outputn);

%网络训练
an=sim(net,inputn);

error=sum(abs(an-outputn));