pso优化bp网络机械故障诊断(MATLAB代码)

在混合算法中,需要优化的对象(粒子)是 BP 神经网络的权值和值。首先应把要优化的神经网络的全部权值和闽值构成一个实数数组,并赋予它们 [0,1之间的随机数。然后,按照选定的网络结构,用前向算法计算出对应于每组输入样本的神经网络输出。这里BP网络的激活函数都选为sigmoid 函数然后用改进PSO算法搜索出最优位置,使如下的均方误差指标(适应度函数值)达到最小:

式中,n为样本个数;c 为网络神经元的输出个数;tkp为第p个样本的第个理想输出值;Y.p为第p个样本的第  个实际输出值。
这样适应度函数值达到最小时搜索到的便是 BP 网络的最佳权值和闻值。下面是混合算法的实现步骤:
(1)确定适应度权值w、最大允许迭代步数 t、搜索范围[-Xmax,Xmax]、最大速度，并根据网络规模确定粒子数。 

通过改进的粒子群算法对网络权值进行训练,改善了使用多层感知器的BP算中存在的网络学习收敛速度慢,以及容易陷入局部极小等问题,提高了网络收敛速度,可以有效地用于设备故障诊断,对设备故障进行特征识别和诊断分析。 

文件夹如下

数据集介绍

data5.m

训练集特征

训练集的标签 

测试集的特征 

 主程序（main.m文件）

%% 该代码为基于PSO和BP网络的预测
%% 初始化
clc
clear
warning off
%读取数据
load data5 input output test

%节点个数
inputnum=9;
hiddennum=13;
outputnum=5;

%训练数据和预测数据
input_train=input';
input_test=test';
output_train=output';
output_test=output';
inputn=input_train;
outputn=output_train;
%构建网络
dx=[0,1;0,1;0,1;0,1;0,1;0,1;0,1;0,1;0,1];
net=newff(dx,[hiddennum,outputnum],{'tansig' 'logsig'},'trainlm');

% 参数初始化
%粒子群算法中的两个参数
c1 = 1.49445;
c2 = 1.49445;

maxgen=100;   % 进化次数  
sizepop=20;   %种群规模

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

for i=1:sizepop
    pop(i,:)=5*rands(1,200);
    V(i,:)=rands(1,200);
    fitness(i)=fun(pop(i,:),inputnum,hiddennum,outputnum,net,inputn,outputn);
end
% 个体极值和群体极值
[bestfitness bestindex]=min(fitness);
zbest=pop(bestindex,:);   %全局最佳
gbest=pop;    %个体最佳
fitnessgbest=fitness;   %个体最佳适应度
fitnesszbest=bestfitness;   %全局最佳适应度

%% 迭代寻优
for i=1:maxgen    
    for j=1:sizepop   
        %速度更新
        V(j,:) = V(j,:) + c1*rand*(gbest(j,:) - pop(j,:)) + c2*rand*(zbest - pop(j,:));
        V(j,find(V(j,:)>Vmax))=Vmax;
        V(j,find(V(j,:)<Vmin))=Vmin;
        
        %种群更新
        pop(j,:)=pop(j,:)+0.2*V(j,:);
        pop(j,find(pop(j,:)>popmax))=popmax;
        pop(j,find(pop(j,:)<popmin))=popmin; 
        
        %自适应变异
        pos=unidrnd(21);
        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(i)=fitnesszbest;            
end

%% 结果分析
plot(yy)
title(['适应度曲线  ' '终止代数＝' num2str(maxgen)]);
xlabel('进化代数');ylabel('适应度');
x=zbest;

%% 把最优初始阀值权值赋予网络预测
% %用pso优化的BP网络进行值预测
w1=x(1:inputnum*hiddennum);
B1=x(inputnum*hiddennum+1:inputnum*hiddennum+hiddennum);
w2=x(inputnum*hiddennum+hiddennum+1:inputnum*hiddennum+hiddennum+hiddennum*outputnum);
B2=x(inputnum*hiddennum+hiddennum+hiddennum*outputnum+1:inputnum*hiddennum+hiddennum+hiddennum*outputnum+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);

%% BP网络训练
%网络进化参数
net.trainParam.epochs=1000;
net.trainParam.lr=0.001;
net.trainParam.goal=1e-6;
net.trainParam.showwindow = 0;

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

%% BP网络预测
an=sim(net,input_test)

效果图

测试集的预测标签

完整代码可以关注

function error = fun(x,inputnum,hiddennum,outputnum,net,inputn,outputn)
% 该函数用来计算适应度
%x          input     个体
%inputnum   input     输入层节点数
%outputnum  input     隐含层节点数
%net        input     网络
%inputn     input     训练输入数据
%outputn    input     训练输出数据

%error      output    个体适应度

%提取
w1=x(1:inputnum*hiddennum);
B1=x(inputnum*hiddennum+1:inputnum*hiddennum+hiddennum);
w2=x(inputnum*hiddennum+hiddennum+1:inputnum*hiddennum+hiddennum+hiddennum*outputnum);
B2=x(inputnum*hiddennum+hiddennum+hiddennum*outputnum+1:inputnum*hiddennum+hiddennum+hiddennum*outputnum+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);
%可以关注https://mbd.pub/o/bread/ZJyTmZZt