改进的麻雀算法优化最大相关峭度解卷积（SCSSA-MCKD），实现早期微弱故障诊断，MATLAB代码实现

01 引言

由于一些设备的早期故障产生的冲击十分微弱，易被系统噪声干扰，如何有效地对设备的原始故障信号进行降噪并增强信号中微弱冲击成分，是进行该类部件早期故障诊断的关键。

最大相关峭度解卷积（MCKD）通过解卷积运算突出被噪声淹没的连续冲击脉冲，提高原始信号的相关峭度值，非常适用于提取微弱故障信号的连续瞬态冲击。已有很多文献证明MCKD可以用于早期微弱故障诊断。

02 流程介绍

本篇文章的思路大家可以参考这篇文献，这篇文献的被引量和下载次数都是相当高的：张俊,张建群,钟敏等.基于PSO-VMD-MCKD方法的风机轴承微弱故障诊断[J].振动.测试与诊断,2020,40(02):287-296+418.DOI:10.16450/j.cnki.issn.1004-6801.2020.02.011.

流程介绍：

①模拟早期故障信号，并添加强烈的高斯白噪声来模拟实际工况设备被环境噪声所覆盖的微弱故障；

②找最佳IMF分量：采用VMD或者EMD等方式对故障模拟信号分解，计算各IMF分量的包络谱峰值因子，包络谱峰值因子最大的IMF分量即为最佳IMF分量；

③采用对最大相关峭度解卷积（MCKD）的三个参数进行寻优，得到增强的IMF分量；

④对IMF分量进行包络谱分析，得到故障特征频率。

03 结果分析

①模拟早期微弱故障信号

信号的原理公式：

信号结果图：

（1）冲击信号，该信号是不加高斯白噪声之前的冲击信号。

（2）故障模拟信号，该信号是加高斯白噪声之后的冲击信号，可以看到，周期冲击信号则完全被噪声淹没。

（3）故障模拟信号的包络谱，从这个图中难以发现突出频率，无法辨别故障特征。

②找最佳IMF分量

采用VMD对故障信号进行分解，得到若干个IMF分量，然后计算IMF分量的包络谱峰值因子，包络谱峰值因子最大的IMF分量即为最佳IMF分量；

将VMD模态分解数设置为K = 10，这个数字是根据经验自己设置的。(大家如果想通过优化VMD参数的方式得到模态分解数K也是可以的，可以看之前发的有关优化VMD参数的文章，这篇文章就不再讲解优化VMD了)

得到各IMF分量的包络谱峰值因子图如下：

可以看到最佳的IMF分量是IMF7。然后接下来的操作都是对这个IMF7展开的。

对IMF7进行包络谱分析，

想获得更多VMD相关分析的可以看这篇文章：VMD分解，matlab代码，包络线，包络谱，中心频率，峭度值，能量熵，样本熵，模糊熵，排列熵，多尺度排列熵，西储大学数据集为例

③改进的麻雀算法优化MCKD参数

原文献中只对MCKD的[L,T]两个参数进行了优化，作者在这里进行了改进，对[L,T,M]三个参数进行了优化，以确保找到最佳的MCKD参数。

适应度函数依旧选择包络谱峰值因子。改进的麻雀算法优化MCKD适应度曲线如下。

处理过程如下：

得到的最佳参数分别是[L,T,M] = [242,99,1];

MCKD处理后的包络谱如下：

可以看到，与优化前的IMF分量包络谱比较，得到了明显的增强。解卷积后的包络谱中故障特征频率fi 及其２倍频、３倍频、4倍频的谱线均清晰可见，表明特征频率被准确提取。

04代码

%% 
clear
clc  
close all
fobj=@mckdcost;       %包络谱峰值因子
%% 选取数据
load data.mat  %注意这里替换为自己的数据即可，数据形式为n行*1列，列数必须为1。
data = best_imf;
%% 设置参数
lb = [100 85 1];    %mckd下限
ub = [1000 142 7];  %mckd上限
dim = 3;            % 优化变量数目
Max_iter=20;       % 最大迭代数目
SearchAgents_no=20;       %种群规模
fs = 12800;   %采样频率
%% 调用SCSSA函数
[fMin , bestX, Convergence_curve ] = SCSSA(SearchAgents_no,Max_iter,lb,ub,dim,fobj,data,fs);


%% 画适应度函数曲线图，并输出最佳参数
figure
plot(-Convergence_curve,'Color',[0.9 0.5 0.1],'Marker','>','LineStyle','--','linewidth',1);


title('Objective space')
xlabel('Iteration');
ylabel('Best score obtained so far');
legend('SCSSA优化MCKD')
display(['The best solution obtained by SCSSA is : ', num2str(fix(bestX))]);  %输出最佳位置
display(['The best optimal value of the objective funciton found by SCSSA is : ', num2str(-fMin)]);  %输出最佳适应度值


%% 将最佳的MCKD参数回带，求出结果


filterSize = fix(bestX(1));
termIter = 30;
T = fix(bestX(2));
M = fix(bestX(3));
plotMode = 0;
%--------------- Run actual mckd code:数据进行mckd分解---------------------------
[y_final,~,~] = mckd(data,filterSize,termIter,T,M,plotMode);

代码目录如下：