基于VMD-PE-MCKD低速重载滚动轴承故障特征提取

毛 欢，魏志刚，刘迎松，韦雅宁，王宏元，陆 强

（1.安徽工业大学 机械工程学院，安徽 马鞍山 243002；2.马钢（集团）股份有限公司，安徽 马鞍山 244500）

滚动轴承是旋转机械的重要组成部分，其工作状况直接影响整个机械系统的稳定与安全[1]。低速重载滚动轴承运行工况恶劣且复杂，由于其运行转速低，故障特征频率低，故难以被提取。低速重载轴承故障信号具有非线性和非平稳性的特点[2]，其故障特征的提取一直以来是研究的重难点。

时频分析方法是一种常用的可以有效处理非平稳信号的方法，例如小波分析、短时傅里叶变换、经验模态分解等。小波分析作为典型的时频分析方法应用广泛，但其缺点是分析结果的精确性依赖小波基函数的选择。经验模态分解（Empirical Mode Decomposition，EMD）作为一种自适应分解的方法，避免了函数选择的问题，但其存在端点效应、模态混叠的问题。集合经验模态分解（Ensemble Empirical Mode Decomposition，EEMD）改善了EMD的模态混叠问题，但其分解效果取决于添加白噪声的选择。自适应噪声的完备集合经验模态分解（Complete Ensemble Empirical Mode Decomposition with Adaptive Noise，CEEMDAN）分解信号时可自适应添加白噪声，可以克服模态混叠，使分解更具完整性。变分模态分解[3](Variational Mode Decomposition，VMD)是Dragomiretskiy 提出的一种非递归自适应信号处理方法，克服了EMD 的不足问题。夏均忠等[4]运用VMD 方法进行轴承故障特征提取，证明VMD 可以克服EMD 模态混叠问题；施杰等[5]将VMD 与深度迁移学习结合实现了轴承故障特征提取与智能诊断。最大相关峭度解卷积[6]（Maximum Correlation Kurtosis Deconvolution，MCKD）是基于最小熵解卷积而提出的新的解卷积方法。该方法以相关峭度为评价指标，可以增强信号的冲击特性，适用于微弱故障信号的诊断。张永鑫等[7]利用MCKD进行故障诊断，表明MCKD 具有良好的降噪和信号增强效果；张洪梅等[8]将MCKD与CEEMDAN相结合提取滚动轴承微弱故障特征，取得了不错的效果。排列熵[9](Permutation Entropy，PE)由Bandt 等最早提出，是一种检测时间序列随机性和动力学突变的方法，具有计算简单、速度快以及抗干扰能力强的特点。周涛涛等[10]运用PE筛选经CEEMD分解后含噪声较多的分量，然后进行去噪，可以用定量的方式确定含噪较多的分量。王涛等[11]将EMD 与PE 相结合，有效提取了故障特征频率及其倍频，验证了PE的有效性。

以上研究充分说明了VMD、MCKD、PE 在滚动轴承故障诊断中的优势，但是这些方法大多应用于实验室内中高速轴承的故障诊断中，在低速轴承中应用较少。由于实际工况下低速重载轴承故障形式多样，干扰源众多，轴承振动信号十分复杂，为了验证以上方法在工程实际中的应用效果，本文提出了基于VMD 和MCKD 并且结合PE 的低速轴承故障诊断方法，通过仿真信号和工程实测信号对其进行了验证。

1 方法原理

1.1 最大相关峭度解卷积(MCKD)

MCKD方法的思想是通过寻找一个最佳FIR滤波器f使周期信号滤波后的峭度值最大，从而突出信号的连续冲击成分，达到增强信号的目的。

相关峭度的定义为：

式中：M为移位数，T为解卷积周期。MCKD算法的目标函数为：

式中：f为滤波器向量，f=[f1,f2,…fL]T；L为滤波器长度。

为了得到使CKM(T)取得最大值的最优滤波器，求解方程，令：

最终的滤波器系数以矩阵的形式表述为：

其中：

MCKD的算法流程如下：

（1）确定冲击信号的周期T、移位数M和滤波器长度L；

（2）计算输入信号x的XT、XT0、(X0XT0)-1；

（3）求滤波器输出信号y；

（4）根据y计算γm与β；

（5）更新滤波器系数f；

如果滤波前后信号ΔCKM(T)＜ε，则停止迭代，跳回步骤（3）。ε为控制迭代终止的较小正数。

1.2 变分模态分解(VMD)

变分模态分解是一种完全非递归的信号处理方法，可以将非平稳信号x(t)分解为一系列围在中心频率周围的模态分量信号，其约束变分模型为：

式中：{uk}:={u1,…,uK}和{ωk}:={ω1,…,ωK}分别为所有分量及中心频率的集合，δ(t)为脉冲函数，K为分量个数。

通过引入Lagrange 乘子和二次惩罚项将式（5）转化为非约束性变分问题。增强的拉格朗日算子L可以描述为：

式中：α为惩罚因子，λ为Lagrange乘子。

利用交替乘子算法对式（6）求最优解：

VMD的算法流程为：

（3）k=k+1，直到k=K，根据式（7）和式（8）更新

（4）更新λ；

（5）根据式（9），ε为一个大于0 的正数，代表精度，判断是否满足收敛条件，是则分解结束；否则返回过程（2）。

1.3 排列熵（PE）

PE具有计算简单、速度快以及抗干扰能力强的特点，对于非线性数据有很好的适用性。其计算方法如下：

对于一离散时间序列{x(n),n=1,2,…,N}，对其进行相空间重构，得到矩阵：

式中：i=1,2,…,K,K为重构分量的数目；d为嵌入维数；τ为延迟时间，对矩阵中的每一行序列进行升序排列，得到一组新序列S(l)={j1,j2,…,jm}，l=1,2,…,k,k≤m!

对于经升序排列后的新序列，d维相空间映射有d!种排列的可能性，S(l)是其中的一种排列形式。

计算每一种序列出现的概率P1、P2、…、PK，则排列熵可表示为：

当Pl=1/d时，Hp(d)达到最大值ln(d!)，对排列熵进行归一化处理可得：

排列熵的大小表示时间序列的复杂性和随机性程度，PE 值越大则时间序列就越随机，反之则时间序列越规则。嵌入维数d和延迟时间τ的选取对排列熵的计算结果会产生影响，根据Bandt的研究，嵌入维数d一般取3～7。参考文献[12]的研究结果，本文嵌入维数d取为6，延迟时间τ取为1。

1.4 VMD-PE-MCKD故障诊断

VMD 方法在其被提出后就在滚动轴承故障诊断中得到相关应用，但大多应用于实验室里中高速轴承的故障诊断，对于低速轴承的应用研究较少，且在工程实际应用效果有待验证。实际工程中，设备运行时周围会伴随强背景噪声，MCKD 可以有效抑制信号中的噪声干扰，为采用VMD方法处理信号创造更佳的条件。同时依据排列熵值越大，信号越随机，故分量包含噪声更多的原理，可更准确筛选分量。经过多次试验后发现，当VMD分解层数K取为4、惩罚因子α取为2 000时，信号分解效果最佳。采用MCKD方法时滤波器长度L取为400，解卷积周期T取为100。

所提方法的流程图如图1所示。

图1 方法流程图

2 仿真信号分析

为了验证本文方法诊断低速滚动轴承故障的有效性，利用仿真信号进行分析，仿真信号x(t)由正弦信号x1(t)、滚动轴承外圈模拟信号x2(t)以及噪声信号x3(t)组成。

仿真信号中f1=2 500 Hz，设定低速轴承外圈故障频率为f0=8 Hz，x(t)的时域波形如图2所示。由于存在噪声干扰，仿真波形十分杂乱，有用信息被淹没，仅通过时域图难以识别故障特征。信号的包络谱如图3所示，同样因为噪声，图中存在较多干扰的频率成分，直接通过包络谱也难以提取故障特征频率，无法判断故障。

图2 仿真信号时域波形

图3 仿真信号的包络谱

采用本文的方法对信号先进行VMD 分解，VMD 分解模态数K=4，惩罚因子α=2 000。分解后所得分量如图4所示。然后计算各分量的排列熵值，如表1所示。根据按照排列熵筛选IMF 分量原则，选取排列熵值较大的IMF2、IMF3、IMF4进 行MCKD去噪。

表1 VMD分解各分量排列熵

图4 经VMD分解后所得的IMF分量

经MCKD去噪后的时域波形如图5所示。对比MCKD去噪前的时域波形可以看出，经过MCKD处理后的信号周期性冲击成分突出，从去噪后的信号图中可以看出与故障特征频率8 Hz 相近的频率成分被提取出来，相对应的还有其2～7 倍频、9 倍频、11 倍频处存在明显峰值，因此可以判定低速滚动轴承外圈出现故障。对去噪后的信号进行重构，然后对重构后的信号进行包络谱处理，如图6所示。

图5 MCKD去噪后的时域波形

图6 经VMD-PE-MCKD处理后的包络谱

为了进一步验证本文方法的优势，对仿真信号采用EEMD 分解，各分量的排列熵值如表2所示。选取排列熵值较大的IMF1、IMF2、IMF3进行MCKD去噪。对去噪后的分量进行重构后的包络谱如图7所示。从包络谱中仅可提取出与故障特征频率对应的1、5、8 倍频，这些频率幅值不明显且易受其他频率的干扰，诊断效果不太理想。

图7 EEMD-PE-MCKD包络谱

表2 EEMD分解各分量排列熵

对仿真信号采用CEEMDAN 方法处理，分解出14 个分量，各分量的排列熵值如表3所示。选取排列熵值较大的IMF1、IMF2、IMF3进行MCKD 去噪。对去噪后的分量进行重构后的包络谱如图8所示。

表3 经CEEEMDAN分解后各分量排列熵

图8 CEEMDAN-PE-MCKD包络谱

从包络谱中仅提取出与故障特征频率对应的4、5、10倍频，诊断效果不及本文所提方法。

3 工程实际信号分析

以某钢厂混合机托轮轴承为研究对象进行分析，轴承型号为23272 CA/W33，规格为Φ5.1m×24.5 m。将振动加速度传感器安装在轴承座上，混合机托轮轴承及传感器安装图见文献[13]。轴承在连续运行一段时间后，检测到异常情况，现对采集到的数据进行分析。采样频率为625 Hz，采样时间为20 s，共采集12 500个样本点。混合机运行时轴承转速为23.5 r/min，计算得到的轴承故障特征频率如表4所示。

表4 故障特征频率/Hz

首先采用EEMD方法分解信号得到13个分量，各分量的排列熵值如表5所示。选取排列熵值较大的IMF1、IMF2、IMF3进行MCKD处理后再进行重构，重构后信号的包络谱如图9所示。从包络谱中找到了与保持架故障频率对应的2倍频及其4倍频，但其4倍频幅值较低，因此故障特征提取效果一般，不能准确判定故障。

图9 EEMD-PE-MCKD包络谱

表5 经EEMD分解各分量排列熵

采用本文方法将信号进行VMD 分解得到4 个分量，分解后的分量如图10所示。其中，各分量的排列熵值如表6所示。取排列熵值较大的IMF3、IMF4重构后进行MCKD去噪，经MCKD去噪后的信号时域波形如图11所示。可见重构信号经MCKD降噪后冲击成分明显得到增强。然后再进行包络谱处理，如图12所示。从包络谱中找到了与保持架故障频率对应的1～4 倍频，特征提取效果明显。同时，与滚动体故障特征频率对应的1倍频和2倍频也被提取出来，因此可以判定保持架与滚动体发生故障。

图10 VMD分量图

图11 分量重构后经MCKD去噪的信号时域波形

图12 经VMD-PE-MCKD后的包络谱

表6 经VMD分解各分量排列熵

经CEEMDAN分解后各个分量的排列熵值如表7所示。取排列熵值较大的前5 个分量进行MCKD去噪，然后再进行包络谱处理，如图13所示。

图13 CEEMDAN-PE-MCKD包络谱

表7 经CEEEMDAN分解后各分量排列熵

从包络谱中仅找到保持架故障频率对应的2倍频及其4 倍频。这些频率幅值低，混杂在其他频率之中，由于受到非常严重的干扰，故障特征提取效果差，因此不能准确判断故障。

4 结语

针对低速重载轴承信号的非平稳、非线性特性以及故障频率低、故障信号微弱的特点，提出了一种基于变分模态分解以及最大相关峭度解卷积的故障诊断方法。运用VMD将信号分解为一系列分量，结合排列熵方法筛选出有用分量，对筛选出的分量进行MCKD去噪，很好解决了信号中存在噪声干扰的问题。通过对低速外圈故障信号进行仿真分析，验证了本方法的有效性。对混合机这类典型低速重载轴承信号的分析表明该方法可以有效提取出故障特征频率。