0 引言

轴承在飞机动力系统、机体结构和操控系统中均有广泛应用^[1-3]。因其具有复杂的运动关系，在空中和地面工作中长期处于不同频率交变载荷作用下，易发生各类故障，在同时叠加恶劣的环境噪声后，导致故障信号具有显著的非线性和非平稳特征，进一步加大了故障诊断的难度^[4-5]。故轴承对飞机运行稳定和飞行安全具有至关重要的影响。

目前，国内外学者基于信号分解和特征提取等手段对轴承故障信号分析展开了广泛研究。文献^[6]针对风力机轴承信号采用集合经验模态分解（ensemble empirical mode decomposition，简称EEMD）和Teager能量算子解调，通过提取瞬时频率和包络谱特征，实现故障信号的诊断。文献^[7]基于灰狼优化算法（grey wolf optimizer，简称GWO）改进变分模态分解（variational mode decomposition，简称VMD），结合分形维数特征，对轴承非线性故障信号特征进行提取。文献^[8]对轴承故障信号首先采用互补集合经验模态分解（complementary ensemble empirical mode decomposition，简称CEEMD）方法进行处理，后进一步通过小波阈值方法实现对信号噪声成分剥离。LIU Wenyi ^[9]等利用改进后的局部均值分解（local mean decomposition，简称LMD）对风力机轴承非高斯和非线性特征进行分析。HAN Tian^[10]等采用Teager能量强化算子和小波方法，去除轴承信号中的噪声成分，而后利用完备集合经验模态分解（complementary ensemble empirical mode decomposition，简称CEEMD）获得信号中的有用信息，达到损伤的精准识别。

常见信号分解方法中，经验模式分解（empirical mode decomposition，简称EMD）和EEMD易出现模态混叠的问题^[11-13]，VMD分解结果受变量设置影响较大，不具备自适应性^[14-15]。Gilles以小波分析为基础，发展得到了经验小波变换（empirical wavelet transform，简称EWT），该方法在保留EMD频带划分优势的同时，极大地缓解了模态混叠，且计算量小和鲁棒性高^[16]。文献^[17]对进行EWT分解后的轴承故障振动信号提取多尺度熵，较EEMD，具有更优异故障特征提取效果。文献^[18]基于EWT和GG聚类算法，实现轴承故障的无监督诊断。JIANG Xingxing^[19]等根据峭度最大化原则筛选经EWT分解后分量，从而实现对轴承故障振动信号的降噪。

轴承故障振动信号通常含有大量的冲击成分^[19-21]。而峭度为冲击性重要的度量手段。因此，本文首先对振动加速度信号采用EWT分解，后计算各阶固有模态函数（intrinsic mode function，简称IMF）时频谱负熵和频域谱负熵，并提出一种基于比例系数的自适应平均谱负熵算法；而后基于GWO以峭度最大化为适应度函数，对自适应平均谱负熵的比例系数和IMF分量数进行择优选取，最终实现故障振动加速度信号重构，从而实现原始信号中无关成分和噪声的剔除。

1 经验小波变换

EWT基于小波理论框架，对时域信号傅里叶谱依据极大值分割频带，而后对各个频带构建小波滤波器^[16]。

1.1 信号分割

在信号分割中，首先将Fourier支撑[0，π]划分成若干个组成Λ_n=[ω_n-1，ω_n]，n=1，2，···，N（ω₀=0，ω_N=π），各组成边缘点为ω_n，而后确定信号频谱相邻极大值的中点，则∪^N_n=1Λ_n=[0，π]，详见图1。ω_n为Λ_n中点，过渡区域范围T_n=2τ_n。

在频带分割中采用了尺度空间法。f（x）在[0，π]范围内Fourier谱，则其尺度空间函数L（x，t）为

式中，$g（x，t）=\frac{1}{\sqrt{2\pi t}}{e}^{-x^2/2t}$为高斯核函数，t为尺度参数。

尺度空间函数L（x，t）极小值随尺度参数增加而减少。Gilles等在文献^[16]提出尺度空间直方图概念。各个局部极小值定义一个长度为L_i的尺度空间曲线，故寻找有意义的模态即要确定阈值T，如极小值曲线长度大于T，则该曲线模态即为有效模态。

尺度空间法即为利用高斯核函数对信号进行平滑处理。随平滑次数增加，虚假频带分界点会消失，真实分界点保留。故记录每次平滑后各尺度下极小值个数，形成尺度和频带分界点的二维图像。每个初始极小值点都将对应一条曲线，保留长度大于阈值T的曲线，其对应的分界点即为保留的频带分界点。

1.2 模式分解

Λ_n的经验尺度函数和经验小波函数为：

式中，0＜γ＜1，ω_n/τ_n=γω_n。

β（x）为：

信号细节系数${\varpi }_f^{\epsilon }（n，t）$为：

近似系数${\varpi }_f^{\epsilon }（0，t）$为：

基于式（5）和式（6）重构信号，得到：

式中，${\varpi }_f^{\epsilon }（n，\omega ）$和${\varpi }_f^{\epsilon }（0，\omega ）$分别为${\varpi }_f^{\epsilon }（n，t）$和${\varpi }_f^{\epsilon }（n，t）$Fourier变换结果。

信号分解结果如下：

2 谱负熵和峭度

2.1 谱负熵

时域谱负熵和频域谱负熵^[22]用于量化信号的冲击性和周期性。

x（n）谱熵可表示为：

式中，H_ε为谱熵；ε_x为信号平方包络；f_C为中心频率；Δf为带宽。

时域谱负熵为：

式中，I_ε为时域谱负熵。

频域谱负熵为：

式中，I_E为频域谱负熵；E_x为信号平方包络谱。

当信号冲击性强时，时域谱负熵升高。当周期性强时，频域谱负熵增大。为综合反映信号特性，分析时通常采用平均谱负熵：

式中，I_1/2为平均谱负熵。

然而实际信号中冲击成分和周期成分多数并非对称分布。故本文提出自适应平均谱负熵，通过比例系数调整两类谱负熵权重：

式中：I_ρ为自适应平均谱负熵；ρ为比例系数。

2.2 峭度

因多数轴承出现故障后，振动信号冲击成分增多，峭度值会随之增加。该指标常作为轴承故障诊断和信号提取重要依据，且通常故障程度和峭度大小存在正相关^[20，23]。

峭度为归一化的4阶中心距，可量化信号的非高斯性：

式中，K为峭度；μ为信号均值；σ为信号标准差；x为信号。

3 灰狼优化算法

GWO为典型群智能算法，通过模拟狼群高效捕猎行为，对目标函数进行优化^[24]。

狼群在捕猎中通常按等级划分为4个等级。GWO设置α、β、δ和ω，代表不同等级的狼，对应最优解、次优解、第三优解和其他结果。狼群捕猎行为如图2所示。

假设狼群数量为L，需对p个变量进行优化，则狼个体表示为${\mathit{W}}_i=\left(w_{i1}，w_{i2}，\cdots ，w_{ip}\right)$。

狼群对猎物，即优化目标的围猎过程可用下式表示：

式中，t为迭代数；${\mathit{w}}_q（t）$为最优解；$\mathit{A}\cdot \left|\mathit{C}\cdot {\mathit{w}}_q（t）-{\mathit{w}}_i（t）\right|$为包围步长；A为收敛系数; C为摆动系数。

A和C可表示为：

式中，r₁和r₂为介于0~1间随机向量；a为控制系数，随迭代次数由2衰减到0。

根据确定α狼、β狼和δ狼的位置W_α、W_β和W_δ，确定其他灰狼ω：

GWO算法伪代码如图3所示。

4 故障信号分析流程

针对轴承振动加速时域量由EWT分解后，采用GWO对自适应平均谱负熵比例系数和重构IMF分量个数进行优化，实现信号重构。优化以峭度最大化为适应度函数。

如图4所示，故障信号分析过程如下：

1）对轴承振动加速度时域量进行EWT分解；

2）对各IMF分量计算时域谱负熵和频域谱负熵；

3）初始化GWO优化中狼群位置，个体变量为自适应平均谱负熵比例系数和重构IMF分量个数；

4）针对每个个体，根据比例系数计算各IMF分量自适应平均谱负熵，以此降序排列分量；

5）针对每个个体，根据分量个数重构原始信号；

6）针对每个个体，计算重构信号峭度值作为适应度函数；

7）依据适应度函数确定α狼、 β狼和δ狼；

8）更新狼群位置；

9）确定是否达到迭代次数或精度目标，如达到则输出最优结果，反之重新计算个体适应度函数；

10）根据最优结果，重构轴承振动加速度信号。

5 信号分析

5.1 分析信号

利用西储大学公开的轴承数据库对本文提出的型号分析方法展开验证。

轴承振动加速度信号采集的试验台如图5所示。故障类型包括外圈故障、内圈故障和滚动体故障。轴承类型为SKF6205。采集信号为振动加速度，频率为12 kHz。加速度传感器安装于驱动轴末端。其中外圈传感器在6点钟方位上。

为分析本文信号分析方法，对原始故障信号增加信噪比为3的白噪声，以模拟实际工作状态下数据受环境噪声污染情况。

故障信号时域变化如图6所示。由图可知外圈故障和内圈故障存在较明显的冲击特征。滚动体故障信号冲击特征较弱，从时域波动中难以观察到显著的周期性。添加噪声后，信号更为无序混乱，部分冲击成分被淹没在噪声中。

故障信号频域分布如图7所示。由图7可知，外圈故障在374 Hz、1 436 Hz、2 794 Hz和3 440 Hz处出现明显峰值，在各阶特征频率位置存在明显边频带。内圈故障主要以360 Hz、2 589 Hz和3 372 Hz为峰值，并在两侧同样存在边频带。滚动体故障在360 Hz和3 728 Hz处峰值较大，在2 500 Hz和3 500 Hz附近存在边频带。各类故障信号在各个频段内都均匀加入白噪声成分，导致部分峰值被淹没，进一步加大特征频率提取和故障诊断的难度。

正常轴承振动加速度时域变化和频率分布如图8所示。信号时域变化中无明显的冲击特征。频域分布以358 Hz、1 036 Hz和2 102 Hz为特征频率。以上频率附近无边频带。

正常和故障信号峭度如表1所示。轴承各类故障皆会导致振动加速度信号峭度值显著提高。噪声添加前，外圈故障、内圈故障和滚动故障信号峭度较正常信号分别提高约7.5倍、8.9倍和2.79倍。正常信号无大量冲击成分，在加入噪声后，其峭度近似不变。因含噪的故障信号中冲击成分被噪声淹没，故峭度值有所下降。

5.2 EWT分解

正常轴承振动加速度信号EWT分解后频域划分结果如图9所示。无噪信号和含噪信号分别得到38个和43个IMF分量。

故障轴承振动加速度信号EWT分解后频域划分结果如图10所示。由图可知，EWT分解结果对特征频率可以准确识别，对该处频率和附近边频带进行较准确的划分。含噪信号IMF 分解层数多于无噪信号。外圈故障无噪和含噪信号分别分解成28个和46个IMF分量。内圈故障分解为40个和46个IMF分量。滚动故障分解层数分别为33个和45个。

正常和故障轴承振动加速度信号EWT分解后，各IMF分量峭度、时域谱负熵和频域谱负熵结果如图11所示。峭度和时域谱负熵对信号冲击成分较为敏感。频域谱负熵用于量化周期成分。轴承信号中同时包含冲击性和周期性。而当其出现故障时，冲击成分将显著提高，具体表现为故障的轴承信号各IMF分量的峭度和时域谱负熵显著提高，且两者具有相同的变化趋势，而与频域谱负熵呈负相关。正常轴承信号较为平稳，各阶IMF分量的峭度和两类谱负熵波动平缓。峭度受噪声影响较大。对比各故障情况，含噪信号的峭度变化明显大于时域谱负熵。由此可见，时域谱负熵具有更优异的抗噪性能。因此基于时域谱负熵和频域谱负熵作为信号筛选依据，可更好地避免噪声影响。针对故障信号，加噪前后，冲击性最大的IMF分量不变。峭度和时域谱负熵虽具有相同的波动趋势，但两者变换幅度不尽相同，由此可见，其对冲击成分的灵敏程度存在差异。

5.3 GWO优化

在对轴承振动加速度信号进行EWT分解后，对自适应平均谱负熵比例系数和重构IMF分量个数进行优化，以达到重构信号峭度最大为目标。GWO优化最大迭代次数为20次，狼群个体数量为5只，各故障信号适应度函数变化和优化参数迭代过程如图12所示。随迭代次数增加，重构信号峭度逐渐提高，并高于原始信号峭度，故冲击成分被成功提取。由图12可知，比例系数和IMF分量个数随优化过程不断调整。在迭代后期适应度（峭度）趋于平稳，故优化过程逐渐收敛至最优结果。

重构信号峭度如图13所示，由图可知，重构后信号峭度总体有明显提高，冲击成分被有效提取。其中外圈故障含噪信号峭度提升71.2%，内圈故障含噪信号峭度提升60.5%。含噪信号峭度增大幅度大于无噪信号。其中内圈故障和外圈故障冲击性更强。滚动体故障冲击成分较少，重构后峭度提升不明显。

信号重构前后时域变化如图14所示。重构后信号中冲击特征更为明显。针对含噪信号，原始信号中白噪声成分在重构后被大量剔除。重构信号信噪比如表2所示，外圈故障、内圈故障和滚动体故障重构后，信噪比分别提高106.7%、30.0%和96.7%。

信号重构前后频域分布如图15所示。正常信号和故障信号在360 Hz附近存在周期成分。通过信号重构，该频率处信号被较好过滤。而各类故障导致的冲击成分主要集中在2 000~3 000 Hz，其特征频率和两侧边频带被成功提取。

6 结论

轴承为飞机重要的机械零件，其故障具有复杂的信号特征，具有明显的冲击性。为实现故障准确诊断，本文利用EWT故障信号进行分解后，基于GWO算法，以峭度最大化为适应度函数对重构IMF分量数量和自适应平均谱负熵比例系数进行优化，实现对信号重构，以期提取故障特征冲击成分，获得主要研究结论如下：

1）EWT对基于轴承信号频率峰值，实现IMF分量自适应划分；

2）时域谱负熵和峭度皆可量化冲击性，但前者数值对噪声敏感性较低；

3）以重构IMF分量和自适应平均谱负熵比例系数为优化目标，通过GWO优化后对轴承信号重构可有效提取故障冲击成分；

4）轴承故障信号重构后，峭度大幅提高，信噪比明显改善。

参考文献