基于粒子群优化的BP神经网络算法

大家好，我是带我去滑雪！

基于粒子群优化的BP神经网络算法（Particle Swarm Optimization Backpropagation Neural Network，PSO-BPNN）是一种利用粒子群优化算法优化BP神经网络的算法。它将BP神经网络的权重和偏置值作为粒子群的位置，并利用PSO算法来搜索最优解。该算法具体许多优点，例如：1、能够在搜索过程中不陷入局部最优解，而是更有可能找到全局最优解或者接近全局最优解的解；2、可以通过利用粒子群的速度来调整BP神经网络的权重和偏置值，使得网络能够更快地优化；3、可以通过直接根据网络性能来更新权重和偏置值，避免了梯度消失问题；4、由于采用了随机性初始化的粒子群，可以减轻对初始值的依赖，使得算法更加稳定和鲁棒；5、PSO-BP算法的粒子群可以并行运行，每个粒子都可以独立进行计算和更新。这使得算法能够有效利用并行计算的优势，加速了训练过程。

PSO-BP神经网络算法的基本步骤：

初始化粒子群的位置和速度。每个粒子的位置表示BP神经网络的权重和偏置值，速度表示更新的步长；
对于每个粒子，根据当前位置计算网络的输出，并计算误差（例如均方误差）；
根据每个粒子的误差评估其适应度，适应度可以使用误差的倒数或其他评价指标来表示；
更新每个粒子的最优位置和最优适应度。如果当前适应度优于历史最优适应度，则更新最优位置；
更新每个粒子的速度和位置。速度的更新考虑了个体经验和群体共享经验，以及随机项，使得粒子能够在搜索空间中进行探索和利用；
重复步骤2到步骤5，直到满足停止条件（例如达到最大迭代次数或达到期望的网络性能）。

下面开始代码实战：

clc;

clear;

tic

close all;

load('basket.mat')#导入数据

P = trains(:,1:end-1) ;

T = trains(:,end) ;

P_test = tests(:,1:end-1) ;

T_test = tests(:,end) ;

cur_season = pred ;

inputnum=size(P,2);

hiddennum=2*inputnum+1;

outputnum=size(T,2);

w1num=inputnum*hiddennum;                               w2num=outputnum*hiddennum;

N=w1num+hiddennum+w2num+outputnum;

nVar=N;

VarSize=[1,nVar];

VarMin=-0.5;

VarMax=0.5;

MaxIt=200;

nPop=359;

w=1;

wdamp=0.99;

c1=1.5;

c2=2.0;

VelMax=0.1*(VarMax-VarMin);

VelMin=-VelMax;

empty_particle.Position=[];

empty_particle.Cost=[];

empty_particle.Velocity=[];

empty_particle.Best.Position=[];

empty_particle.Best.Cost=[];

particle=repmat(empty_particle,nPop,1);

GlobalBest.Cost=inf;

for i=1:nPop

    particle(i).Position=unifrnd(VarMin,VarMax,VarSize);

    particle(i).Velocity=zeros(VarSize);

    particle(i).Cost=BpFunction(particle(i).Position,P,T,hiddennum,P_test,T_test);

    particle(i).Best.Position=particle(i).Position;

    particle(i).Best.Cost=particle(i).Cost;

    if particle(i).Best.Cost<GlobalBest.Cost

        GlobalBest=particle(i).Best;

    end

end

BestCost=zeros(MaxIt,1);

for it=1:MaxIt

    for i=1:nPop

        particle(i).Velocity = w*particle(i).Velocity ...

            +c1*rand(VarSize).*(particle(i).Best.Position-particle(i).Position) ...

            +c2*rand(VarSize).*(GlobalBest.Position-particle(i).Position);

        particle(i).Velocity = max(particle(i).Velocity,VelMin);

        particle(i).Velocity = min(particle(i).Velocity,VelMax);

        particle(i).Position = particle(i).Position + particle(i).Velocity;

        IsOutside=(particle(i).Position<VarMin | particle(i).Position>VarMax);

        particle(i).Velocity(IsOutside)=-particle(i).Velocity(IsOutside);

        particle(i).Position = max(particle(i).Position,VarMin);

        particle(i).Position = min(particle(i).Position,VarMax);

        particle(i).Cost=BpFunction(particle(i).Position,P,T,hiddennum,P_test,T_test);

        if particle(i).Cost<particle(i).Best.Cost

            particle(i).Best.Position=particle(i).Position;

            particle(i).Best.Cost=particle(i).Cost;

            if particle(i).Best.Cost<GlobalBest.Cost

                GlobalBest=particle(i).Best;

            end

        end

    end

    BestCost(it)=GlobalBest.Cost;

    disp(['Iteration ' num2str(it) ': Best Cost = ' num2str(BestCost(it))]);

    w=w*wdamp;

end

BestSol=GlobalBest;

%% Results

figure;

%plot(BestCost,'LineWidth',7);

semilogy(BestCost,'LineWidth',7);

xlabel('Number of iterations')

ylabel('The variation of the error')

title('Evolution')

grid on;

fprintf([' Optimal initial weight and threshold:\n=',num2str(BestSol.Position),'\n Minimum error=',num2str(BestSol.Cost),'\n'])

cur_test=zeros(size(cur_season,1),1);

[~,bestCur_sim]=BpFunction(BestSol.Position,P,T,hiddennum,cur_season,cur_test);

prob=softmax(bestCur_sim);

disp(['1次尝试通过游戏的概率为 ',num2str(prob(1))]);

disp(['2次尝试通过游戏的概率',num2str(prob(2))]);

disp(['3次尝试通过游戏的概率',num2str(prob(3))]);

disp(['4次尝试通过游戏的概率',num2str(prob(4))]);

disp(['5次尝试通过游戏的概率',num2str(prob(5))]);

disp(['6次尝试通过游戏的概率',num2str(prob(6))]);

disp(['7次以上尝试通过游戏的概率',num2str(prob(7))]);

toc

function [err,T_sim]=BpFunction(x,P,T,hiddennum,P_test,T_test)

inputnum=size(P,7);

outputnum=size(T,7);

[p_train,ps_train]=mapminmax(P',0,1);

p_test=mapminmax('apply',P_test',ps_train);

[t_train,ps_output]=mapminmax(T',0,1);

net=newff(p_train,t_train,hiddennum);

net.trainParam.epochs=1000;

net.trainParam.goal=1e-3;

net.trainParam.lr=0.01;

net.trainParam.showwindow=false;                                    w1num=inputnum*hiddennum;                                           w2num=outputnum*hiddennum;                                        W1=x(1:w1num);                                                      B1=x(w1num+1:w1num+hiddennum);                                      W2=x(w1num+hiddennum+1:w1num+hiddennum+w2num);                    B2=x(w1num+hiddennum+w2num+1:w1num+hiddennum+w2num+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 = train(net,p_train,t_train);

t_sim = sim(net,p_test);

T_sim1 = mapminmax('reverse',t_sim,ps_output);

T_sim=T_sim1';

err=norm(T_sim-T_test);

index0= T_sim<0;

index1= T_sim>1;

penalty=1000*abs(sum(T_sim(index0)))+1000*sum(T_sim(index1)-1);

err=err+penalty;

end