CN110245595A - 基于时频谱图和二维主成分分析的机械轴承故障诊断方法 - Google Patents

Abstract

本发明涉及机械故障诊断、模式识别技术等领域，具体为一种基于时频谱图和二维主成分分析的机械轴承故障诊断方法，包括接收原始信号数据并进行预处理，获得原始信号数据的初始向量；对初始向量进行短时傅里叶变换时频分析，得到原始信号所对应的短时傅里叶变换三维时频谱图；对短时傅里叶变换三维时频谱图进行灰度图转化，得到二维的时频谱图；对二维时频谱图进行维数约减，提取得到二维时频谱图所对应的低维特征矩阵；利用t‑sne算法对二维主成分分析的输出矩阵进行二次维数约简处理；采用t‑sne算法的输出矩阵训练支持向量机，完成故障诊断模型的构建，输入测试数据即可输出诊断结果；本发明减少了计算时间，并提高机械故障诊断的精度。

Description

基于时频谱图和二维主成分分析的机械轴承故障诊断方法

技术领域

本发明涉及机械故障诊断、模式识别技术等领域，具体为一种基于时频谱图和二维主成分分析的机械轴承故障诊断方法。

背景技术

现代机械设备是多领域融合、多尺度耦合的复杂系统，导致机械设备容易发生故障。而轴承作为机械设备中的关键零件，也是易发生故障的位置。通常情况下机械设备的故障以轴承振动信号的形式表现出来，轴承振动信号中包含了大量的设备运行状态信息，并且轴承振动信号的测试方便、快捷、低廉，因此基于轴承振动信号的故障诊断方法是目前国内外应用最广泛的故障诊断方法。

目前，应用于机械设备轴承故障诊断的各种方法大部分都是利用信号处理技术从轴承原始振动信号中提取的统计特征进行故障诊断，这些特征不能完整地表征轴承振动信号中包含的关于机械设备运行状态的信息，导致故障诊断的效果并不理想。

近年来，图像识别领域已经受到越来越多的关注，并且取得了一些应用成果，如：人脸识别、车牌识别等，但是如何将图像识别技术应用到机械设备轴承故障诊断仍然是值得我们深入探讨的崭新课题。

通常来说，轴承振动信号是一维数据，而图像是二维或者三维数据，要将图像识别技术应用于轴承机械故障诊断，首先必须将一维振动数据转化为二维图像。时频分析方法，如短时傅里叶变换等，能够将一维振动信号转化为三维的时频谱图，将三维的时频谱图进行灰度图转换就能得到二维的时频谱图。时频谱图能够更好地表征轴承振动信号，可以反映出信号成分随时间的变化情况，因此从时频谱图中提取的轴承机械故障特征比从轴承原始振动信号中提取的机械特征包含更多的故障信息。

然而，二维时频谱图的维数过高、存在“维数灾难”的问题，且现有机器学习算法，如支持向量机等，无法直接采用二维图像数据来进行模型训练，无法识别出机械故障。

发明内容

为了更好地提取轴承振动信号的中包含的故障信息，获得了更好的模式识别效果，本发明提出一种基于时频谱图和二维主成分分析的机械轴承故障诊断方法，包括如下步骤：

S1：接收机械轴承的振动信号数据并进行预处理，并获得机械轴承的振动信号数据的初始向量；

S2：对初始向量进行短时傅里叶变换时频分析，得到机械轴承的振动信号数据所对应的短时傅里叶变换三维时频谱图；

S3：对短时傅里叶变换三维时频谱图进行灰度图转化，得到二维的时频谱图；

S4：利用二维主成分分析对二维时频谱图进行维数约减，提取得到二维时频谱图所对应的低维特征矩阵；

S5：利用t-分布领域嵌入(t-distributed stochastic neighbor embedding，t-sne)算法对二维主成分分析的输出矩阵进行二次维数约简处理；

S6：输入机械轴承的振动信号数据，利用t-分布领域嵌入算法的输出矩阵来训练支持向量机，完成故障诊断模型的构建；将待测机械轴承的振动信号数据输入故障诊断模型，即可得到诊断分类结果。

进一步的，接收机械轴承的振动信号数据并进行预处理包括：接收到的机械轴承的振动信号数据为机械轴承的振动信号向量，对振动信号向量进行零-均值规范化，向量中的每个幅度值减去整个向量幅度值的均值之后再除以整个向量幅度值的标准差，得到了经过处理之后的振动信号初始向量x₀，初始向量x₀包括1024个点。

进一步的，所述对初始向量进行短时傅里叶变换时频分析包括：

S_x(ω,t)＝∫x(τ)m(τ-t)e^-jωtdt；

其中，t为时间，ω为角频率，τ为变量，S_x(ω,t)为原始信号对应的短时傅里叶变换三维时频谱图。

进一步的，所述对短时傅里叶变换三维时频谱图进行灰度图转化包括将短时傅里叶变换三维时频谱图中的所有元素归一化到灰度区间，得到灰度图像A，灰度图像A中第i行第j列的元素y_ij表示为：

y_ij＝(|x_ij|-x_min)/(x_max-x_min)×255；

其中，x_ij为S_x(ω,t)矩阵中的第i行第j列的元素；x_max为S_x(ω,t)矩阵中模值最大的元素，x_min为S_x(ω,t)矩阵中模值最小的元素，|·|为求模运算。

进一步的，所述二维主成分分析包括直接利用灰度图图像矩阵构建图像协方差矩阵，并导出特征向量进行图像特征提取，得到二维时频谱图所对应的特征矩阵，具体表示为：

Y＝AX；

S_x＝E[(Y-E(Y))(Y-E(Y))^T]；

J(X)＝tr(S_x)；

G_t＝E[(A-E(A))^T(A-E(A))]；

J(X)＝X^TG_tX；

Y_k＝AX_k,k＝1,2,...,d.；

其中，Y为灰度图像矩阵A的投影特征向量，X为投影酋向量，S_x为训练样本投影特征向量的协方差矩阵，E[·]为求期望，G_t为图像协方差矩阵，{X₁,...,X_d}为最优投影轴对应J(X)最大d个的特征值，Y_k为最优投影轴对应的投影特征向量。

进一步的，二次维数约简处理包括：

其中，qi_j表示y_i和y_j之间的配对相似性，y_k表示第k个低维数据点。

进一步的，采用t-分布领域嵌入算法的输出矩阵训练支持向量机包括；

其中，||·||²表示二范数，w为超平面的法向量，C为惩罚参数，ε_i为松弛变量，y_i为对第i个超平面分类的标签；b为位移量，即到超平面的距离；φ(x_i)为经过核函数处理后的样本点x_i，N表示支持向量机的维度。

本发明采用二维主成分分析来对二维时频谱图进行维数约简，将二维时频谱图压缩为一维特征向量，从而实现故障特征的提取，相较于主成分分析，二维主成分分析直接由原始图像矩阵计算图像协方差矩阵，而无需将原始图像矩阵构建出一维向量再计算图像协方差矩阵，大大减少了计算时间；同时，采用二维主成分分析直接从二维时频谱图中提取的故障特征包含更加丰富的故障信息，可以提高机械故障诊断的精度。

附图说明

图1为本发明基于时频谱图和二维主成分分析的故障诊断方法的流程示意图；

图2为本发明实施例中的滚动体轻微故障的初始向量示意图；

图3为本发明实施例中的滚滚动体中度故障的初始向量示意图；

图4为本发明实施例中的滚动体严重故障的初始向量示意图；

图5为本发明实施例中的内圈轻微故障的初始向量示意图；

图6为本发明实施例中的内圈中度故障的初始向量示意图；

图7为本发明实施例中的内圈严重故障的初始向量示意图；

图8为本发明实施例中的外圈轻微故障的初始向量示意图；

图9为本发明实施例中的外圈中度故障的初始向量示意图；

图10为本发明实施例中的外圈严重故障的初始向量示意图；

图11为本发明实施例中的滚动体轻微故障的三维时频谱图示意图；

图12为本发明实施例中的滚动体中度故障的三维时频谱图示意图；

图13为本发明实施例中的滚动体严重故障的三维时频谱图示意图；

图14为本发明实施例中的内圈轻微故障的三维时频谱图示意图；

图15为本发明实施例中的内圈中度故障的三维时频谱图示意图；

图16为本发明实施例中的内圈严重故障的三维时频谱图示意图；

图17为本发明实施例中的外圈轻微故障的三维时频谱图示意图；

图18为本发明实施例中的外圈中度故障的三维时频谱图示意图；

图19为本发明实施例中的外圈严重故障的三维时频谱图示意图；

图20为本发明实施例中的滚动体轻微故障的灰度图示意图；

图21为本发明实施例中的滚动体中度故障的灰度图示意图；

图22为本发明实施例中的滚动体严重故障的灰度图示意图；

图23为本发明实施例中的内圈轻微故障的灰度图示意图；

图24为本发明实施例中的内圈中度故障的灰度图示意图；

图25为本发明实施例中的内圈严重故障的灰度图示意图；

图26为本发明实施例中的外圈轻微故障的灰度图示意图；

图27为本发明实施例中的外圈中度故障的灰度图示意图；

图28为本发明实施例中的外圈严重故障的灰度图示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明提出了一种基于时频谱图和二维主成分分析的故障诊断方法，包括如下步骤：

S1：接收机械轴承的振动信号数据并进行预处理，并获得机械轴承的振动信号数据的初始向量；

S2：对初始向量进行短时傅里叶变换时频分析，得到机械轴承的振动信号数据所对应的短时傅里叶变换三维时频谱图；

S3：对短时傅里叶变换三维时频谱图进行灰度图转化，得到二维的时频谱图；

S4：利用二维主成分分析对二维时频谱图进行维数约减，提取得到二维时频谱图所对应的低维特征矩阵；

S5：利用t-分布领域嵌入算法对二维主成分分析的输出矩阵进行二次维数约简处理；

S6：输入机械轴承的振动信号数据，利用t-分布领域嵌入算法的输出矩阵来训练支持向量机，完成故障诊断模型的构建；将待测机械轴承的振动信号数据输入故障诊断模型，即可得到诊断分类结果。

在本实施例中，主要对凯斯西储大学的轴承故障数据中的9种故障振动信号进行分析，包括：外圈轻微故障、外圈中度故障、外圈严重故障、内圈轻微故障、内圈中度故障、内圈严重故障以及滚动体轻微故障、滚动体中度故障、滚动体严重故障。

首先对机械轴承的振动信号数据进行预处理，得到初始向量x₀，x₀是关于时间t的函数；所述预处理，即对得到的振动信号向量进行零-均值规范化，是指振动向量中的每个幅度值减去整个向量幅度值的均值之后再除以整个向量幅度值的标准差，可以表示为：

其中，x表示振动信号向量中的第i个幅度值，表示初始向量中的第i个值，μ表示振动向量中所有幅度值的平均值，σ表示振动向量中所有幅度值的标准差；得到了经过处理之后的振动信号初始向量x₀，初始向量x₀包含1024个点；滚动体轻微故障的初始向量如图2、滚动体中度故障的初始向量如图3、滚动体严重故障的初始向量如图4、内圈轻微故障的初始向量如图5、内圈中度故障的初始向量如图6、内圈严重故障的初始向量如图7、外圈轻微故障的初始向量如图8、外圈中度故障的初始向量如图9、外圈严重故障的初始向量如图10。

对初始向量进行短时傅里叶时频分析，得到原始信号所对应的短时傅里叶变换三维时频谱图；所述短时傅里叶时频分析，是指以初始向量x₀为函数，对其进行如下的变换：

S_x(ω,t)＝∫x(τ)m(τ-t)e^-jωtdt

其中，t为时间，ω为角频率，τ为变量，S_x(ω,t)为得到的短时傅里叶变换三维时频谱图，S_x(ω,t)为一个复数矩阵。对初始序列x₀处理得到的数据体现为1024×1024矩阵的形式，将得到的短时傅里叶变换时频分解矩阵S_x(ω,t)中的元素记为x_ij；其中，滚动体轻微故障的三维时频谱图如图11、滚动体中度故障的三维时频谱图如图12、滚动体严重故障的三维时频谱图如图13、内圈轻微故障的三维时频谱图如图14、内圈中度故障的三维时频谱图如图15、内圈严重故障的三维时频谱图如图16、外圈轻微故障的三维时频谱图如图17、外圈中度故障的三维时频谱图如图18、外圈严重故障的三维时频谱图如图19。

将短时傅里叶变换三维时频谱图中的所有元素归一化到[0-255]的区间，将其转化灰度图像；所述灰度图像中包含了原始振动信号中的信息，归一化处理可以表示为：

y_ij＝(|x_ij|-x_min)/(x_max-x_min)×255；

其中，x_ij为S_x(ω,t)矩阵中的第i行第j列的元素，yi_j为A矩阵中第i行第j列的元素，x_max为S_x(ω,t)矩阵中模值最大的元素，x_min为S_x(ω,t)矩阵中模值最小的元素，|·|为求模运算，A为灰度图像矩阵。通过上述归一化处理，原短时傅里叶变换时频分解矩阵中的元素全部转化为[0-255]间的数，进而转化为灰度图；其中，滚动体轻微故障的灰度图如图20、滚动体中度故障的灰度图如图21、滚动体严重故障的灰度图如图22、内圈轻微故障的灰度图如图23、内圈中度故障的灰度图如图24、内圈严重故障的灰度图如图25、外圈轻微故障的灰度图如图26、外圈中度故障的灰度图如图27、外圈严重故障的灰度图如图28。

对灰度图进行二维主成分分析包括直接利用灰度图图像矩阵构建图像协方差矩阵，并导出特征向量进行图像特征提取，得到二维时频谱图所对应的特征矩阵，具体表示为：

Y＝AX；

S_x＝E[(Y-E(Y))(Y-E(Y))^T]；

J(X)＝tr(S_x)；

G_t＝E[(A-E(A))^T(A-E(A))]；

J(X)＝X^TG_tX；

Y_k＝AX_k,k＝1,2,...,d.；

其中，Y为灰度图A的投影特征向量，X为投影酋向量，S_x为训练样本投影特征向量的协方差矩阵，J(X)表示矩阵S_x的迹，E[·]为求期望，G_t为图像协方差矩阵，{X₁,...,X_d}为最优投影轴对应的J(X)中最大的d个特征值，E表示误差矩阵；Y_k为最优投影轴对应的投影特征向量。由实验可知，主成份个数的确定一般以其累计贡献率大于0.86而小于0.91时所对应的个数选取时最为适宜。在本实施例中，X₁所对应的特征值已达到0.86的累积贡献率，故最后得到的投影特征矩阵为Y₁。

得到的投影特征矩阵Y₁为向量，若直接用于训练支持向量机，维度依然过大，会导致模型训练时间过长和过学习的问题，故本发明使用t-分布领域嵌入算法对得到的投影特征向量再次降维，得到t-分布领域嵌入算法降维后的数据，训练支持向量机模型，并对测试集进行测试分类，得到最分类结果，采用t-分布领域嵌入算法的输出矩阵训练支持向量机包括；

其中，||·||²表示二范数，w为超平面的法向量，C为惩罚参数，ε_i为松弛变量，y_i为对第i个超平面分类的标签；b为位移量，即到超平面的距离；φ(x_i)为经过核函数处理后的样本点x_i，N表示支持向量机的维度。

机械设备故障诊断的目的是对故障类型进行分类，其本质是模式识别问题。机械设备故障诊断的主要步骤包括信号预处理、特征提取、故障识别等，而特征提取对故障诊断精度起着至关重要的作用。

为了验证本发明在机械故障诊断中的应用效果，通过凯斯西储大学的轴承故障数据来说明本方法的有效性。本发明可以用于包括但不限于轴承在内的机械故障诊断。

凯斯西储大学的轴承故障数据是利用电火花加工切槽来模拟轴承不同部位的单点故障，主要是通过设置不同的切槽宽度：0.007、0.014、0.021(单位：英寸，1英寸＝2.54厘米)，来模拟轴承轻度故障、中度故障以及严重故障。本次实验选用主轴转速1772rpm、载荷为1HP、采样频率为12000Hz的驱动端振动轴承的9种故障振动信号进行分析，包括：外圈轻微故障、外圈中度故障、外圈严重故障、内圈轻微故障、内圈中度故障、内圈严重故障以及滚动体轻微故障、滚动体中度故障、滚动体严重故障。分别在各种故障状态下测取100组样本，每组样本采样点数为1024点，并在原始数据中加入5db的白噪声。将其中50组数据作为训练样本，采用t-sne降维，其余50组数据作为测试样本(out of sample)。分别从每个样本的振动信号中提取本发明所述特征以及其余三种特征提取方法：AR模型(auto regressivemodel，简称AR模型)、时频特征和熵特征，都通过t-sne降维处理后，采用参数优化的支持向量机来进行故障诊断。表1展示了4种不同的特征提取方法在加入5db的白噪声最后的模式识别准确率。从表中可以看出，其他特征提取方法的模式识别准确率达不到本发明的分类精度。因而本发明在机械故障诊断中是具有很好的应用背景。

表1本发明与AR、时频、熵对应的模式识别准确率

尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同物限定。

Claims (7)

1.基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，包括以下步骤：

S1：接收机械轴承的振动信号数据并进行预处理，并获得机械轴承的振动信号数据的初始向量；

S2：对初始向量进行短时傅里叶变换时频分析，得到机械轴承的振动信号数据所对应的短时傅里叶变换三维时频谱图；

S3：对短时傅里叶变换三维时频谱图进行灰度图转化，得到二维的时频谱图；

S4：利用二维主成分分析对二维时频谱图进行维数约减，提取得到二维时频谱图所对应的低维特征矩阵；

S5：利用t-分布领域嵌入算法对二维主成分分析的输出矩阵进行二次维数约简处理；

S6：输入机械轴承的振动信号数据，利用t-分布领域嵌入算法的输出矩阵来训练支持向量机，完成故障诊断模型的构建；将待测机械轴承的振动信号数据输入故障诊断模型，即可得到诊断分类结果。

2.根据权利要求1所述的基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，接收机械轴承的振动信号数据并进行预处理包括：接收到的机械轴承的振动信号数据为振动信号向量，对振动信号向量进行零-均值规范化，向量中的每个幅度值减去整个向量幅度值的均值之后再除以整个向量幅度值的方差，得到了经过处理之后的振动信号初始向量x₀，初始向量x₀包括1024个点。

3.根据权利要求1所述的基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，所述对初始向量进行短时傅里叶变换时频分析包括：

S_x(ω,t)＝∫x(τ)m(τ-t)e^-jωtdt；

其中，t为时间，ω为角频率，τ为变量，S_x(ω,t)为原始信号对应的短时傅里叶变换三维时频谱图。

4.根据权利要求1所述的基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，所述对短时傅里叶变换三维时频谱图进行灰度图转化包括将短时傅里叶变换三维时频谱图中的所有元素归一化到灰度区间，得到灰度图像A，灰度图像A中第i行第j列的元素y_ij表示为：

y_ij＝(|x_ij|-x_min)/(x_max-x_min)×255；

其中，x_ij为S_x(ω,t)矩阵中的第i行第j列的元素；x_max为S_x(ω,t)矩阵中模值最大的元素，x_min为S_x(ω,t)矩阵中模值最小的元素，|·|为求模运算。

5.根据权利要求1所述的基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，所述二维主成分分析包括直接利用灰度图图像矩阵构建图像协方差矩阵，并导出特征向量进行图像特征提取，得到二维时频谱图所对应的特征矩阵，具体表示为：

Y＝AX；

S_x＝E[(Y-E(Y))(Y-E(Y))^T]；

J(X)＝tr(S_x)；

G_t＝E[(A-E(A))^T(A-E(A))]；

J(X)＝X^TG_tX；

Y_k＝AX_k,k＝1,2,...,d.；

其中，Y为灰度图像矩阵A的投影特征向量；X为投影酋向量；S_x为训练样本投影特征向量的协方差矩阵；E[·]表示求期望；J(X)表示矩阵S_x的迹；G_t为图像协方差矩阵；{X₁,...,X_d}为最优投影轴对应的J(X)中最大的d个特征值；Y_k为最优投影轴对应的投影特征向量，即二维时频谱图所对应的特征矩阵；上标T表示转置矩阵。

6.根据权利要求1所述的基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，所述二次维数约简采用t-分布领域嵌入算法对图像特征矩阵进行降维，二次维数约简处理包括：

其中，q_ij表示y_i和y_j之间的配对相似性，y_k表示第k个低维数据点。

7.根据权利要求1所述的基于时频谱图和二维主成分分析的机械轴承故障诊断方法，其特征在于，采用t-分布领域嵌入算法的输出矩阵训练支持向量机包括：

其中，||·||²表示二范数，w为超平面的法向量，C为惩罚参数，ε_i为松弛变量，y_i为对第i个超平面分类的标签；b为位移量，即到超平面的距离；φ(x_i)为经过核函数处理后的样本点x_i，N表示支持向量机的维度。