基于VMD和排列熵的滚动轴承故障诊断研究*

杨 云，张昊宇，薛元贺,2，丁 磊

(1.华东交通大学电气与自动化工程学院，南昌 330033；2.中国铁路南昌局集团有限公司，南昌 330033)

0 引言

滚动轴承是旋转机械的重要组成部分，如果发生故障，在其工作时会造成安全隐患，因此判断出轴承的当前状态，并采取相应的措施处理十分必要[1]。

轴承的振动信号多为非线性、非平稳信号，经验模态分解(Empirical mode decomposition，EMD)作为早期的自适应分解算法与其他方法结合激起了人们的研究兴趣。王金东等[2]通过EMD分解轴承的振动信号最后结合SVM，可以较准确实现故障诊断；由于EMD算法本身的缺陷，马丽华等[3]提出基于集合经验模态分解(EEMD)和GG聚类的方法，最终的结果也证明了该方法的可行性。2005年由EMD理论更新出了一种新的自适应分解算法局部均值分解算法(Local mean decomposition，LMD)，侯高雁等[4]针对滚动轴承的振动信号，通过LMD形态学与EEMD形态学在故障提取中的对比研究，表明LMD对故障提取有着速度快，清晰度高的特点。人们在研究过程中，发现EMD、LMD两者算法本身存在局限性，无法解决端点效应和模态混叠现象。Dragomiretskiy Konstantin在2014年提出了变分模态分解(Variational modal decomposition，VMD)方法[5]，此方法可以避免端点效应、抑制模态变分模态分解和样本熵的特征提取方法,采用支持向量机进行故障识别。该方法能精确的实现故障诊断。熵作为一种构建特征向量的方法也广泛用于轴承的故障诊断领域[7-8]。

本文基于滚动轴承的故障运行机制提出一种变分模态分解和基于峭度准则排列熵构建特征向量的方法，通过SVM进行分类诊断。并通过实例信号进行分析，结果表明提出的方法可以实现滚动轴承的故障诊断。

1 变分模态分解原理及其参数选取

本节介绍变分模态分解原理和步骤，并且分析参数设置对分解结果的影响。

1.1 变分模态分解原理及步骤

1.1.1 变分模态分解原理

变分模态分解就是寻求K个估计带宽之和最小的模态函数，并且要求所有模态函数之和为原函数[6]，约束变分模型表达式：

(1)

式中，mk表示分解得到的K个IMF分量；ωk表示量的中心频率。

求解式(1)变分问题的最优解，引入增广lagrange函数：

(2)

式中，α为罚因子；λ为lagrange乘子。

利用交替方向乘子算法(ADMM)求取上述变分问题，最终结束迭代得到K个IMF分量。

1.1.2 变分模态分解步骤

变分模态分解算法的步骤如下：

(2)n=n+1，进入循环；

(3)根据更新公式进行更新mk,,ωk，直至分解个数达到K时停止循环；

(4)根据公式更新λ；

(5)给定精度ε，若满足停止条件，则停止循环，否则返回步骤(2)继续循环。

1.2 参数设置对分解结果的影响

变分模态分解算法包含的参数有分解尺度K、惩罚因子α、噪声容限和判别精度，研究发现，噪声容限和判别精度对变分模态分解的结果影响较小，本小节通过定一求二法分析确定变分模态分解参数，介绍K和α对分解的影响。

本小节采用西储凯斯大学轴承数据库数据进行分析，以采样频率为12 kHz下的驱动端轴承的滚动体故障数据做分析，图1为该故障信号的时域图，横坐标为时间，纵坐标为幅值。

图1 轴承滚动体故障时的振动信号

1.2.1 分解个数K对分解结果的影响

首先惩罚因子设定为1000，分别设置K值为2、3、4、5，最终分解个数和中心频率结果如表1所示。根据表1可知当K值小于3时，分解尺度明显不足；当K大于4时的中心频率2698和2840比较接近，很可能存在过分解现象，初步设定分解个数为4。

表1 不同分解个数K所对应的中心频率

1.2.2 惩罚因子α对分解尺度的影响

根据上面K值的确定，设定K值为4，再设置当α等于500、1000、1500、2000、2500时，各个分量的中心频率数值如表2所示。

表2 不同惩罚因子α对应的中心频率

根据表2所示可知α的值小于1000时对于分解结果出现欠分解现象，α的值大于1000时趋于稳定，因此再次设置α的值为600、700、800、900分析各个分量的中心频率发现α的值在700时不会出现欠分解现象，通过此方法分析滚动体故障信号大约K=4,α=700左右的时候通过变分模态分解算法分解结果最优。

通过上述定一求二的分析最终得到的4种状态参数组合如表3所示。

表3 定一求二设置得到的参数组合

2 基于峭度准则的排列熵特征向量构建方法

考虑到轴承故障时故障冲击随时间存在周期性，本文提出了基于峭度准则的排列熵特征向量构建方法，下面分别介绍其概念。

2.1 排列熵概念

排列熵是衡量以为时间序列复杂程度的熵，具体原理如下[9]：

(1)对X(i),i=1,2…,n的一个时间序列进行空间重构，得到矩阵如下所示：

(3)

式中，m为嵌入维数；τ为延迟时间；Z为重构相空间的向量个数。

(2)矩阵中每行可看为一个重构分量，将矩阵中的第j重构分量{x(j),x(j+τ),…,x(j+(x(j+(m-1)τ)}以升序的方法排列，得到：

x(j+(i1-1)τ)≤x(j+(i2-1)τ)≤…≤

x(j+(im-1)τ)

(4)

式中，i1,i2,…,im为重构分量中的每个元素所在列的索引。

如果重构分量存在相等值，那么：

x(j-(ip-1)τ)=x(j-(iq-1)τ)

(5)

则根据ip,iq的大小排序，如果ip

x(j-(ip-1)τ)≤x(j-(iq-1)τ)

(6)

因此，对于重构矩阵Y的任一重构分量Y(j)都将得到一组位置索引序列：

B(j)=(i1,i2,…,im)，j=1,2,3,…,k

(7)

式中，k≤m！，B(j)是符号序列的其中一种。

(3)算出每一位置索引序列出现的概率p1，p2…，pk，时间序列X={x(i)，i=1，2，3…，n}的不同位置索引序列的排列熵可以定义为：

(8)

当pi=1/m！时，Ep(m)就达到最大值lnm！

(4)为了方便各运行状态下的排列熵比较，通常用lnm!将Ep(m)进行归一化处理，即：

Ep=Ep(m)/lnm!

(9)

2.2 滚动轴承故障运行状态下模态分量的峭度分析

峭度为描绘波形尖峰度的参数，其数学描述公式为：

(10)

式中，a为信号的均值；σ为信号的标准差。

通过分析VMD分解后模态分量的峭度来分析滚动轴承故障状态下的运行机理，以西储凯斯大学数据库中的内圈故障为例进行分析。

选取内圈故障2048个采样点，通过上一节对不同状态下VMD参数组合[6,2000]，经过VMD分解滚动体故障信号得到的6个分解模态，根据峭度表述公式得到的不同模态峭度值如表4所示。

表4 不同模态峭度值

根据上表可以得知经过VMD分解后的模态分量当n=2、4、5、6时峭度值较高。

选取不同采样起始点，同样选取2048个采样点，进行VMD分解得到的模态分量峭度值如图2所示。

图2 不同采样起始点下的模态分量峭度值

根据上图分析在不同采样起始点下的峭度值同样遵循n=2、4、5、6时峭度值较高的规律，此为滚动轴承的故障运行机理。

2.3 基于峭度准则的排列熵特征向量构建方法

轴承无故障运行时，振动信号近似接近正态分布，此时的峭度值近似为3；当轴承发生故障时，其振动信号概率密度增大，信号幅值即偏离正态分布，峭度指标的绝对值越大，轴承的故障就越严重。因此，IMF中计算得到峭度绝对值越大，含有故障冲击的成分就越多[10]。

排列熵计算过程中，延迟时间τ和嵌入维度m的选取对于排列熵的计算结果有一定影响[11]，对于这两个参数的选取，本文跟经验选取τ=6，m=1。

由于第1节中分析得到的参数组合中滚动体故障时的K值最小，所以分别取经过VMD分解后内圈故障、外圈故障和正常状态下模态分量中峭度值较大的4个分量用以构建特征向量，具体做法如下：

(1)通过VMD算法分解振动信号，得到若干IMF分量；

(2)计算各个分量的峭度值，选取较大的4个分量；

(3)分别计算4个分量的排列熵值，以构建特征向量。

通过此方法得到的内圈故障用以构建特征向量的模态分量n为2、4、5、6；外圈用以构建特征向量的模态分量n为2、3、4、5；正常状态下用以构建特征向量的模态分量n为1、4、5、6。

3 滚动轴承状态识别流程

本文采用定一求二的方法确定VMD算法的参数组合[K,α]，经过分解得到各个IMF分量，针对故障信号的特点，采用峭度准则和排列熵方法构建特征向量，分析不同故障类型下的排列熵值，确定不同IMF分量下的特征向量，最后将构建的特征向量输入SVM中进行训练和模态分类，诊断流程如图3所示。

图3 故障诊断流程图

具体实现步骤为：

(1)获取实验数据，载入原始信号；

(2)通过定一求二的方法确定VMD的分解参数，并通过分解信号得到的各个IMF分量；

(3)计算不同状态类型下不同模态的峭度值，通过K值的局限，来确定用以构建特征向量的模态；

(4)通过计算不同模态的排列熵构建特征向量；

(5)建立SVM模型，将训练数据和待测数据输入其中得出结果。

4 实例分析

本文选用西储凯斯大学实验室的数据作为实例分析。以采样频率为12 kHz下的驱动端轴承故障数据做分析，分别采用损伤直径为0.177 8 mm的内圈、滚动体、外圈故障和正常状态下的振动信号，各选取4种状态振动信号的60组样本，其中训练样本40组，20组作为测试样本。

由于定一求二方法得到的参数组合(见表3)下4种状态中最小K值为4的限制，分别选取不同状态下的4个模态分量计算排列熵构建特征向量。

以不同采样点为例介绍不同特征向量的构建，分别通过VMD分解内圈故障、外圈故障、滚动体故障以及正常状态下的轴承振动信号，根据第3节的特征向量构建方法得到的4类特征如表5所示。

表5 特征向量构建值表

表5中标签1～4分别表示为正常状态、内圈故障、外圈故障以及滚动体故障滚动轴承状态类型；由于最小K值为4的限制，特征1～4分别表示为4个的模态分量对应的排列熵值。

由于特征向量数目过多在此不加以展示，将上述计算得到的40组样本输入到支持向量机中训练后，输入测试样本得到的识别结果如图4所示。

图4 支持向量机测试数据分类图

根据上图所示的结果可以看出，4种状态下测试数据分类的正确率分别为95%、100%、90%和95%，最终诊断的平均正确率为93.75%，因此本文提出的方法可以良好的实现滚动轴承的故障诊断。

5 结论

本文提出变分模态分解、排列熵以及支持向量机结合的滚动轴承故障诊断方法，通过实例信号实验分析得出结论如下：

(1)变分模态分解的参数设置对于分解结果尤为重要，良好的参数设置对于后期的故障诊断有着重要作用。

(2)通过峭度准则分析了滚动轴承的运行机理，并根据最终的分类结果可以得出利用峭度准则结合VMD和排列熵构建特征向量的方法，可以较好的实现故障诊断。