SE-KSD 优化的FRFT-VMD 轴承故障诊断方法

2024-03-21韦明辉江丽霞涂凤秒姜蓬勃

机械科学与技术 2024年2期

韦明辉 , ，江丽霞 , ，涂凤秒 , ，姜蓬勃 ,

（1.西南石油大学机电工程学院，成都 610500；2.石油天然气装备技术四川省科技资源共享服务平台，成都 610500）

滚动轴承对机车等旋转器械是不可或缺的重要部件之一，但由于其工作环境致使其具有易损坏的特点，而滚动轴承的运行状况对于设备有着巨大影响，任何一个小故障都可能影响设备的正常运行从而给国家以及人民带来巨大损失。若能在滚动轴承失效之前便对其进行故障诊断，不仅能最大限度减少损失，还能降低设备的维修成本[1]。因此，对滚动轴承进行故障诊断具有必要性。然而，采集到的滚动轴承故障信号常混有大量干扰信号，使得所采集到的信号具有非平稳、非线性等特点。由于干扰信号的存在，导致在进行故障诊断时，无法从中快速辨别故障特征信号，使故障诊断的难度上升，时间增加。如何从这些信号中提取所需的故障特征信息成为许多学者研究的问题[2]。

Huang 等[3]提出了经验模态分解方法（Empirical mode decomposition，EMD），该方法根据信号在时间尺度上的局部特征结构，自适应地提取反应信号本质特征的固有模态分量，可该算法缺少理论支撑，会出现端点效应和模态混叠现象影响信号处理结果。吴振华等[4]将噪声辅助信号分析方法引入到EMD中，提出了集成经验模态分解（EEMD），虽然该方法在抑制模态混叠方面起到了一定作用，可效果不彻底，且会产生多个虚假分量，而且噪声的加入会导致重构误差增加。Dragomiretskiy 等[5]提出了变分模态分解（Variational mode decomposition，VMD），一种全新的非递归自适应信号处理方法，近年来，该方法成为了一个研究热点，但都只是在整数域上的改进。由于算法本身的局限性，过分解和欠分解的问题无可避免。

本文结合上述文献，提出了一种基于样本熵（Sample entropy，SE）和峭度均方差（Kurtosis standard deviation，KSD）优化的分数阶变分模态分解（Fractional variational mode decomposition，FRFT-VMD）和随机森林（Random forest，RF）分类器相结合的轴承故障诊断方法，首先，通过搜寻最小SE 值确定分数阶Fourier 变换的最佳阶数，并对信号进行最佳阶数分数傅里叶变换，其次，通过峭度均方差准则来确定变分模态分解的最优参数，将变换后的信号进行变分模态分解，之后，计算分解后信号的峭度和脉冲因子作为特征向量，最后，将特征向量输入至随机森林分类器进行故障识别分类。通过与不同参数的FRFT-VMD分解所得各模态分量对比，本文方法所得的模态分量具有更多故障特征信息，与VMD-RF 的故障诊断方法比，本文作者所提出的基于SE-KSD 优化的FRFT-VMD-RF 方法所得故障特征信息聚类性更强，具有更高的识别准确率。

1 分数阶Fourier 变换

1.1 分数阶Fourier 变换定义

分数阶傅里叶变换作为傅里叶变换在分数域空间上更为广义的存在方式，使得信号可以在以时间轴旋转到任意角度坐标轴上进行表示，在一定程度上对时域信息和频域信息进行融合，进而将其某部分特征进行突出显示[6]。

传统的傅里叶变换被定义为存在于信号空间中的连续线性算子ξ，其特征方程为

式（1）中，传统傅里叶变换所对应的特征值为λn=e-jnπ/2，其特征函数为Hermite-Gauss 函数Hn(t)e-t2/2，其中Hn(t)为n阶Hermite 多项式，表达式为

从式（1）看出，随机信号 ψn的傅里叶变换等价于其自身与复数 ψn之积。接下来将引入一个定义，该定义可由分数阶傅里叶变换基本定义进行直接推导得到。

定义1 设 ψn表示为Hermite-Gauss 函数，作为普通Fourier 变换中特征值为 λn的特征函数，并可用作构成有限信号空间的特征函数。然后，分数阶傅里叶变换可以定义为线性和满足

根据上述表达式，有限能量函数xn可以展开为Fourier 变换特征函数的线性叠加，其表达式为

式（5）被称为分数阶Fourier 变换核的频谱展开，Hermite-Gauss 函数 ψn须得满足：

式中： α=pπ/2 ，当p=1时，分数阶Fourier 变换变为传统Fourier 变换。同时用核函数的形式表示分数阶Fourier 变换，即：

式中：Kp（u,t）为分数阶傅里叶变换的核函数，且同时可以看出以p为变量的参数是以4 为周期的，因此F4n和F4n±2分别相当于恒等算子I和奇偶算子P。当p=0时，经分数阶傅里叶变换后得到的是原信号。其逆变换形式为

由图1 可知，在时频面上，将（t,ω）按逆时针方向作 α角度的旋转，变换到（u,v）域，便可以得到新的信号表示。因此，利用分数阶傅里叶变换，通过改变角度 α 或阶数p，可以观察到信号从时域到频域的所有特征信息。

图1 时频平面上分数阶傅里叶变换旋转示意图Fig.1 Schematic diagram of fractional Fourier transform rotation on the time-frequency plane

1.2 分数阶Fourier 变换阶次选择

在分数阶Fourier 域可更好解决时频域难以解决的数据可分性差等问题，但阶次的选择对其影响较大。

本文采用样本熵[7]为指标，度量数据在分数阶域的集中程度[8]。样本熵是用来评价混乱度的指标[9]，样本熵值越大，表明包含的有效信息越少，反之表明表明包含的有效信息越多[10]。

样本熵算法步骤如下：

1）令原始数据 {Xi}={x1,x2,···xn} ，长度为N，根据原始信号重构一个m维的向量，即

2）定义dij为xi与xj的距离，且为两者对应元素差值的绝对值的最大值，即

3）计算dij小于相似容限r的数目及此数目与距离dij的总数N-m-1 的比值，用Bmi(r)表示，即

4）计算Bmi(r)的平均值

5）对维数m+ 1 重复上述过程1）到5），得到(r) ，便可得到得到Bm+1（r）

6）样本熵的定义为

当N为有限数时，式（14）可表示为

2 变分模态分解

2.1 变分模态分解定义

VMD 算法核心思想是认为原始信号每个模态的绝大多数特征内容都是紧紧围绕在某一中心频率周围的，通过维纳滤波、频率混合和希尔伯特变换的方式将各子模态的带宽问题分解为多个约束问题，再利用一系列迭代求解的方法得到中心频率的最优解[11]。相较于传统EMD 方法相比，VMD 的数学推导基础更加坚实，具有更好的收敛表现和更出色的鲁棒性。

综上所述，VMD 最后是将给定信号f(t)分解成k个固有模态分量（IMF），并且保证分解后的IMF分量的总带宽和最小。该算法实现步骤如下；首先通过希尔伯特变换得到各个IMF 分量的解析信号，并且计算出相应的单边频谱[12]，即

之后预先估计各个IMF 分量的中心频率，通过频率混合使得各个IMF 的频谱调制到对应的基频带上，即

最后由高斯平滑指标计算得到调制信号梯度平方L2的L2范数，据此估计每个IMF 分量的带宽，通过约束方程来使得IMF 的带宽的和为最小，该方程表达式为

式中： dt为对函数求时间t的偏导数； {uk}为信号f(t) 分解得到的k个IMF 分量， {uk}={u1,u2,u3,···,uk}；{ωk} 为各IMF 的中心频率， {ωk}={ω1,ω2,ω3,···,ωk}。引入含有二次惩罚因子 α和拉格朗日乘数算子 λ的增广拉格朗日方程，使得上述变分约束方程转变为了非约束性变分方程，从而得到该模型的最优解：

因为惩罚因子和拉格朗日乘数算子的增加，从而使得模型在高斯白噪声存在的情况下具有更好的重构精度和更加严格的约束条件。通过交替方向乘子法对式（19）中的、及 λn+1进行交替更新，从而得到非约束变分方程的鞍点。

利用埃尔米特实信号的对称性，将式（21）转换为非负频率的半空间积分方程，同时进行二次优化可得的更新方程为

同理根据带宽之和最小的约束条件，再利用傅里叶等距变换将最小问题转换到频域中，得到与 λn+1交替更新方程：

VMD 采用递归循环的方式对原始信号进行层层剥离，有效避免了EMD 模态混叠的缺点，具有自适应的能力[13-14]。

2.2 变分模态分解算法参数选择

变分模态分解算法虽然克服了传统经验模式分解及其改进方法的缺点，但分解前需要设定分解层数K和惩罚因子α，参数的选择对分解结果的影响很大[15]。本文采用峭度均方差准则来对参数进行寻优[16]。峭度是反映随机变量分布的数值统计值，对滚动轴承振动信号类的冲击信号非常敏感[17]。均方差可以反映数据组中个体之间的离散程度，衡量分布程度的结果。通过峭度均方差准则则可反应VMD 分解后不同模态分量之间的差异性，峭度均方差值越大，表明不同模态分量间差异越大，反之则表明不同模态分量间差异越小。

3 轴承故障特征提取

本文提出一种基于样本熵和峭度均方差优化的分数阶变分模态分解的故障特征提取新方法。该方法的核心思想是通过样本熵搜寻分数阶傅里叶变换的最优阶次，把原始数据中可分性差的数据映射到合适的分数阶空间，再通过峭度均方差准则得到变分模态最优参数后在分数域对其进行变分模态分解。

实验数据由CWRU 电气实验室轴承数据中心提供，本文选取的轴承不同状态数据如表1 所示。

表1 轴承不同状态数据Tab.1 Different status data of bearings

在故障特征提取部分，选用外圈故障这一组数据，该数据故障特征频率为107.365 Hz，基于分数阶Fourier变换的变分模态分解算法流程图如图2 所示。

图2 基于分数阶Fourier 变换的变分模态分解算法流程Fig.2 Procedures of variational modal decomposition algorithm based on fractional fourier transform

步骤1 计算分数阶阶次P从0.01 依次增加0.01 至2 的样本熵值，样本熵最小值处的阶次即为分数阶傅里叶最优阶次，对数据进行最佳阶数分数阶Fourier 变换。

步骤2 计算α从1 000 依次增加100 到10 000，K从2 依次增加1 到11 的峭度均方差值，进行变分模态分解的分解层数K和惩罚因子α寻优，对分数阶Fourier 变换后的数据进行变分模态分解。

步骤3 对变分模态分解后的每个变分模态分量进行分数阶Fourier 逆变换。

步骤4 对每个变分模态分量进行包络谱分析。

计算数据在0 ～ 2 阶次分数阶Fourier 变换时的样本熵值，由图3 可知在X＝93 时，样本熵值最小，所以分数阶Fourier 最优阶次为0.93。

图3 分数阶Fourier 最优阶次Fig.3 Fractional Fourier optimal order

计算最优阶次分数阶傅里叶变换后惩罚因子α 由1 000 依次递增100 到10 000，分解层数K由2 依次递增1 到11 的峭度均方差值，由图4 可知，惩罚因子α增加了15 次，即为2 500，分解层数K增加了3 次，即为4，由此可知变分模态分解最优参数α和K分别为2 500 和4。图5 为基于SE-KSD 优化的FRFT-VMD 分解后的各模态分量。

图4 惩罚因子和分解层Fig.4 Penalty factor and decomposition layer

图5 分数阶VMD 分解各分量Fig.5 Fractional order VMD decomposition of each component

从图6 可看出，基于SE-KSD 优化的FRFT-VMD分解所得各模态分量中包含更多更明显的故障特征信息。对比不同参数的FRFT-VMD 所提取的特征频率的幅值，本文所提出的故障特征提取方法所得的特征频率及其相应的倍频的幅值明显大于其他参数FRFT-VMD 的幅值，表明本文所提方法可以为轴承故障诊断提供更加有效的故障特征信息。

图6 不同参数的FRFT-VMD 各分量包络谱图Fig.6 Envelope spectra of FRFT-VMD components with different parameters

4 实验验证

4.1 实验流程

为了进一步验证该方法在轴承故障诊断中的作用，将采用西储大学数据库数据以及自己搭建的实验平台所测数据进行验证，实验流程如图7 所示。

图7 实验流程图Fig.7 Experimental flow chart

4.2 数据库故障诊断

采用的实验数据来自西储大学轴承数据库，选用的轴承不同状态的数据如表1 所示。轴承各状态下原始时域图如图8 所示。

图8 轴承各状态下原始时域图Fig.8 Original time domains diagram of bearings in different states

为了说明基于SE-KSD 优化的FRFT-VMD 方法在轴承故障诊断中的作用，本文将峭度和脉冲因子作为故障特征向量，对基于SE-KSD 优化的FRFT-VMD 所提取的特征向量和基于VMD 所提取的特征向量进行K值聚类分析，从图9 聚类结果图可知，基于SE-KSD 优化的FRFT-VMD 所提取的特征向量具有较强的聚类性，不同故障之间的区分度更高，而基于VMD 所提取的特征向量聚类效果较差，故障类型之间区分不清且较分散。

为了进一步验证所提方法的有效性，本文采用峭度和脉冲因子作为故障特征向量，使用随机森林分类器[18-19]进行故障诊断。首先将不同类型的故障分为10 000 个数据为一组的10 组数据，总共分为40 组数据，对这40 组数据分别进行基于SE-KSD优化的FRFT-VMD 分解，再对所得变分模态分量提取峭度和脉冲因子作为故障特征向量，每个故障类型将会得到40 组峭度和脉冲因子特征向量，将各种故障类型中32 组作为随机森林分类器的训练集，8 组作为测试集数据对整体模型的分类准确率进行验证，并将各轴承状态进行编号，如表2 所示，设置随机森林分类器中的决策树数目为800，随后进行故障诊断。将相同数据使用变分模态分解后，提取相应故障特征向量，同样使用随机森林进行故障诊断。

表2 各状态轴承对应编号Tab.2 Corresponding numbers of bearings in each state

两种方法训练得到的随机森林结构图如图10所示。通过测试集的测试，结果显示两种方法的分类结果图如图11 所示，从图中可以看出，基于SE-KSD优化的FRFT-VMD 特征提取方法所提取的不同故障类型间的特征向量可以更好的进行识别并分类，分类结果图显示预测结果与期望值基本吻合，而基于VMD 特征提取方法的分类结果中，预测结果与期望值之间有较大出入。

图10 随机森林结构图Fig.10 Random forest structure

图11 分类结果图Fig.11 Classification results

两种方法分类结果如表3 和表4 所示。由表3和表4 可知，基于SE-KSD 优化的FRFT-VMD-RF故障诊断方法所得的分类总的准确率结果为96.875%，且对于编号1、2、3 的故障类型诊断准确率为100%。使用基于VMD-RF 故障诊断方法所得的分类总的准确率结果为78.1%。实验结果显示，基于SE-KSD 优化的FRFT-VMD 分解所得到的特征向量具有更好的区分度，且训练得到的分类器模型具有更好的泛化能力。

表3 基于SE-KSD 优化的FRFT-VMD-RF 分类结果汇总Tab.3 Summary of FRFT-VMD-RF classification results based on SE-KSD optimization

表4 基于VMD-RF 分类结果汇总Tab.4 Summary of VMD-RF classification results

4.3 实测数据故障诊断

本节由自行搭建的滚动轴承实验平台所测数据进行故障诊断实验，滚动轴承实验平台由旋转机械和数据采集系统两部分组成，如图12 所示。旋转机械部分由正常/缺陷滚动轴承和三相鼠笼式电动机组成，轴承安装在三相鼠笼式电动机的风扇端，转速通过光电传感器测得。数据采集系统由安装于电机风扇端的ICP 加速度传感器和INV 3060A 型24 位网络分布式采集仪（含采集内嵌服务软件Coinv DASP V10）组成。轴承损伤方式采用EDM。在内圈、外圈和滚动体上分别加工一条深度和宽度为0.2 mm 的直线，模拟内圈故障、外圈故障和滚动体故障。采样频率为12 kHz。实验平台如图13 所示。

图12 实验平台框图Fig.12 Experimental platform's block diagram

图13 实验平台图Fig.13 Experimental platform

通过实验平台获取了4 种状态下振动数据，长度均为119 000 点左右。将每种状态的振动数据进行截取分组，每组10 000 个采样点，每种状态20 组，共得到80 组数据，对这80 组数据分别进行基于SE-KSD 优化的FRFT-VMD 分解，再对所得变分模态分量提取峭度和脉冲因子作为故障特征向量，每个故障类型将会得到80 组峭度和脉冲因子特征向量，将各种故障类型中56 组作为随机森林分类器的训练集，24 组作为测试集数据对整体模型的分类准确率进行验证，设置随机森林分类器中的决策树数目为800，训练得到的随机森林结构图如图14 所示，通过测试集的测试，结果显示分类准确率为90.625%，分类结果如图15 所示。