刘远红，蔡 煜，张彦生

（东北石油大学电气信息工程学院，黑龙江 大庆 163318）

1 引言

实际生产中对机器的运行状态进行监测尤为重要。在进行滚动轴承实时监测时会产生大量的冗余信息，给分析运算带来巨大的计算压力，导致无法及时发现问题、排除故障。

为了提取信号的特征频率并应用于故障诊断中，先后出现了：短时傅里叶变换[1]、小波变换[2]、S变换[3]、经验模态分解[4]等信号分解方法。但是单一的信号变换方法无法有效提取轴承信号的全面特征。与此同时，由于信号的数据量呈指数型增长，逐个处理庞大的轴承信号数据集也变得不切实际。于是许多降维算法相继被提出：PCA算法[5]、LDA算法[6]、ISOMAP算法[7]、LTSA算法[8]、LPP算法[9]、LLE算法[10]等。但这些基于向量处理的线性或非线性降维算法也会在不同程度上导致信号中特征信息的丢失。

近年来，基于张量的特征提取方法陆续得到了众多学者的关注。文献[11]中利用支持向量机回归方法机械制动传感器系统的故障诊断，而其中不包括故障类型的识别。张等人提出了多线性子空间与LLE算法相结合的MLSLLE算法，即在信号的多线性子空间中提取信号的初始特征，相当于对轴承信号首先实施一次降维处理，再利用LLE算法得到最终的降维结果[12]，但该方法需要新旧样本同时参与计算。文献[13]中利用张量的标准分解得到原始张量的因子矩阵和权重，将其重构回一维信号并提取信号特征。

设计了轴承振动信号张量模型的故障诊断方法。首先利用小波变换和相空间重构进行张量建模，通过HOSVD[14]挖掘出信号的重要特征，经过拟合增强不同故障类型曲线间的区分度，最后依据提出的区分度函数估计出各类数据间的最优分类点。在此基础上，基础提出了轴承信号张量模型故障分类方法。

2 轴承信号的张量建模与分解

信号在不同的变换方法下获得的结果是不同的，投影到不同的子空间时就可以获得不同的描述特征。由这些子空间构建的张量模型可以综合考虑信号的各方面特征，这时的张量模型可以视为一个多维矩阵，实现对信号的多角度观测，提高后续轴承故障诊断的精度。

为获取张量模型中轴承故障的特征信息，对张量模型进行分解，选取分解结果中具有代表性的部分作为轴承信号的初始特征。常用的张量分解有CP分解和Tucker分解。CP分解是指将张量分解为R组N个秩一向量外积的形式。但其分解的结果不唯一，求取R的最小值是一个NP问题。Tucker分解是指把N阶张量分解成一个N阶核心张量和N个因子矩阵的n模式乘积的形式，如式（1）所示。张量的n模式展开，如图1所示。其中因子矩阵反映了相应阶上的特征分布。

图1 张量不同维度展开形式Fig.1 Flattening a Tensor along Different Dimensions

高阶奇异值分解（High-Order Singular Value Decomposition，HOSVD）属于Tucker分解的特殊形式，其中得到的伴随矩阵是呈正交关系的。对轴承振动信号按照时间、频率、子空间三个不同维度建立三阶张量模型，利用HOSVD进行分解的具体过程为：将张量展开成一系列的矩阵，如图1所示。再对这些矩阵进行SVD分解，得到各维度上的对角奇异值矩阵、列空间展开项、行空间展开项。

式中：χ∈CΦ×Γ×Η—三阶轴承信号张量模型；Φ—模型的频率维度；Γ—模型的子空间维度；H—模型的时间维度。χ中的元素用xϕ，γ，η表示。由此，χ的等效模式展开矩阵为：

式中：XF、XS、XT矩阵—沿频率维度、子空间维度和时间维度切割得到的矩阵按行展开的形式。

由矩阵的奇异值分解知，存在正交矩阵U∈ℝm×m和V∈ℝn×n，使得矩阵X∈ℝm×n有X=UΛVT。其中，Λ=∈ℝm×n为伪对角矩阵，Λ=diag(σ1，σ2，σ3，…，σG)的对角元素是X的奇异值，G为矩阵X的秩数。并且有：

则（2）式在分解后得到：

为了选择尽可能少的元素来概括尽可能多的特征用于分析故障类别，在以上得到的分解结果中，选择ΛF，ΛT，ΛS矩阵中的对角元素作为初始特征。

3 轴承故障诊断系统

轴承故障诊断系统由特征提取子系统、特征优化子系统和故障诊断子系统三部分组成，虚线箭头表示已知故障类型的轴承振动信号的训练过程，实线箭头表示新采集的轴承振动信号的故障诊断流程，如图2所示。在特征提取子系统中，通过对轴承振动信号进行小波变换和相空间重构，建立其三阶张量模型，经HOSVD获得轴承信号的初始特征。在特征优化子系统中，对轴承信号的初始特征进行处理，得到初始特征曲线在最优分类点处的特征值。根据区分度函数，在已知故障类别的初始特征曲线上找到轴承故障的最优分类点，同时得到在最优分类点处各类轴承状态的特征值区间。将得到的特征值区间输入故障诊断子系统中，等待新的轴承振动信号经过特征提取子系统和特征优化子系统处理后得到其最优分类点处的特征值，在故障诊断子系统中进行比对，得出轴承的故障类型。

图2 轴承故障诊断系统结构Fig.2 Structure of Bearing Fault Diagnosis System

3.1 特征提取子系统

文献[12]利用小波变换和S变换对信号进行处理，构建轴承信号的张量模型。小波变换和S变换均为实现信号时频变换的信号分解方法，但S变换的计算复杂度高、效率低，因此只保留其中对信号的小波变换。

信号x(t)在小波基下进行展开，得到的连续小波变换为：

式中：ɑ—尺度参数；τ—平移参数。为了解决计算量问题，ɑ，τ均选为2j(j∈N+)的倍数，得到的离散小波变换为：

令轴承振动时间序列共有w个采样点，表示为{Wj}，其中，j=1…w，信号的小波变换结果表示为，分解层数为n。得到的小波变换空间为：

为了利用轴承振动一维时间序列估计出原系统的动态特性，引入信号的相空间重构[15]。使用延迟重构法，令轴承振动时间序列共有w个采样点，表示为{Wj}，其中j=1…w，相空间参数为嵌入维数m和时间延迟τ。得到重构的轴承信号相空间为：

式中：i—相空间中的第i个相点。

利用轴承信号的小波变换与相空间重构结果建立轴承信号的张量模型。通过小波变换获取轴承信号在时频上的特征。通过相空间重构估计出与原系统相同的动力学特征，避免实际工程中的多次采样时轴承初始校准位置不对齐的情况。

图3 张量表示Fig.3 Tensor Representation

在初始特征提取时，在HOSVD处理后的结果中，选用ΛS中的对角元素来代替庞大的张量模型，作为轴承振动信号的初始特征IF（theinitial features）。

3.2 特征优化子系统

为便于对初始特征进行观察和分析，寻找这些相似性作为故障诊断的依据，在二维空间中绘制初始特征曲线（IF曲线）。由式（3）可知IF向量中的元素均为非负元素，且按照递减的顺序排列。

将在特征提取子系统中得到的IF向量集描述为：

式中：C—轴承状态类别总数；pn—第n个轴承振动信号在特征提取子系统中得到的IF向量，每个IF向量含有（J+1）个元素。由于轴承信号在采集时伴有干扰，且轴承在工作中也存在持续地缓慢磨损，因此增加了IF曲线中无意义极值的分布。这些无意义极值点相比于其相邻极值点变化速度过快或过慢，使得C簇IF曲线出现少数的交叠现象。这些曲线是非平滑的，不利于通过函数运算进行处理。在保留IF向量中元素非负递减趋势的情况下，对式（10）进行抛物线拟合结果为：

拟合后的曲线簇Q去除了原IF曲线上具有干扰性质的极值点，同时保留了IF曲线的主体下降趋势。把Q投影在二维空间中，自变量区间为（0，J），则区间内一定存在一点j使得同类IF的Q（j）聚集度最高，非同类IF的Q（j）分散度最高。于是定义区分度函数L：

式中：L1—在j点处Q中所有IF曲线占用的总空间；L2—在j点处C簇曲线拟合结果分别占用的空间和。则最优分类点j可通过计算如下目标函数获得：

式中：j=argminL(j)—初始特征的最优分类点。

3.3 故障诊断子系统

在特征优化子系统中，对已知故障类型的IF曲线进行处理，通过区分度函数得到故障最优分类点，将其输入故障诊断子系统中。当机械设备采集到新的轴承振动信号时，在特征优化子系统中计算新的IF曲线在最优分类点处的特征值，将其与原有信号的故障分类信息比对，得出当前轴承的故障状态。

4 实验

4.1 不同IF提取方法的对比

使用江苏千鹏轴承实验平台采集轴承正常状态、外圈故障、滚珠故障和内圈故障的轴承振动信号，轴承转速1000r/min，采样频率5kHz，如图4所示。每类状态的轴承信号随机选取80个，同时建立两种不同的张量模型，对比初始特征效果。第一组采用MLSLLE算法中的张量建模方式，模型的子空间由轴承信号的小波变换和S变换构成。其中小波变换选用db3小波作为母小波，分解层数为7。S变换选取变换结果的前八层。第二组张量模型由特征提取子系统完成。

图4 轴承测试平台Fig.4 Bearing Test Platform

横坐标表示每个IF向量中的16个特征值（J=15），纵坐标表示IF值，如图5所示。图5（a）、图5（b）为第一组张量模型的IF曲线与其拟合结果。图5（a）中在从第四个IF值开始，即横坐标3的位置处，正常状态与内圈故障状态的IF值开始出现重叠，可用于区分类别的区间仅为（0，3）。经拟合后的IF曲线变得平滑（图5（b）），可区分空间扩大到（0，5）。图5（c）、图5（d）为第二组张量模型的IF曲线与其拟合结果。图5（c）中除仅有一条滚珠故障IF曲线在3处和内圈故障曲线相近，其余各类IF值在（0，11）区间内差距明显。图5（d）中各类IF曲线拟合后的可区分区间扩大到（0，13），同时忽略了个别IF曲线中的无用极值点，增大了各类IF曲线之间的可区分度。从上述实验可以看到，第二组张量模型中各类IF的区分度更高。同时在对IF曲线拟合后，可以增大不同类特征之间的差异，并能扩大可选取分类区间的范围。

图5 不同算法的初始特征Fig.5 Initial Features Extracted by Different Algorithms

4.2 轴承故障诊断系统效果验证

选用美国凯斯西储大学轴承故障诊断实验中心的轴承数据验证轴承故障诊断系统效果。数据集采集参数为：驱动端轴承型号为SKF6205，轴承损伤直径0.007英寸，采样频率为12kHz。400个轴承振动信号随机分成由320个信号构成的训练样本集和80个信号构成的测试样本集。训练样本用于确定IF最优分类点。测试样本用于验证故障诊断系统的准确性。

特征优化子系统中输出的训练样本和测试样本IF曲线，如图6所示。

图6 特征优化子系统的IF曲线Fig.6 IF Curve of Feature Optimization Subsystem

在（6，12）区间计算训练样本IF曲线的L值，如表1所示。在6处L值为0.3035，（6，8）区间L值逐渐下降，8处L值为0.2955，（8，12）区间内L值呈现递增趋势。L在8处取得最小值，确定此处为各轴承状态的最佳分类点并输入到故障诊断子系统中。测试样本集经特征提取子系统和特征优化子系统后，得到新的IF曲线，计算其在最优分类点处的IF值（图6（b））。在故障诊断子系统中，将其与训练样本的输入信息比较，实现故障类型的识别。各状态的训练样本和测试样本IF在最优分类点处的数值分布的范围，如表2、表3所示。可以清晰地看出各个轴承状态的IF值都有各自的分布区间且没有任何交集，符合用于轴承故障的分类的条件，如图7所示。测试样本在最优分类点处的IF值几乎全部落在训练样本的IF值分布区间中。个别在区间外的，也可以通过计算最相近IF，判断出故障类别，各类轴承状态的IF值分布区间差距大，不会对判断故障类别造成影响，证明了轴承故障诊断系统的实用性和准确性。

表1 IF曲线的区分度函数值Tab.1 Discrimination Function Value of IF Curve

表2 训练样本IF值统计Tab.2 IF Numerical Values of Training Samples

表3 测试样本IF值统计Tab.3 IF Numerical Values of Testing Samples

图7 故障诊断子系统对比结果Fig.7 Comparison Results of Fault Diagnosis Subsystems

同时，为提高轴承故障诊断的准确性，将历次实验中得到的已知故障类型的轴承初始特征并入训练样本中，丰富系统中训练样本的初始特征内容。用盒图统计了利用上述凯斯西储轴承信号数据进行的三次迭代实验结果，如图8所示。其中的样本集使用情况，如表4所示。样本集1作为训练样本，其余的作为迭代实验的测试样本。

图8 迭代实验结果Fig.8 Iterative Experimental Results

表4 迭代实验样本集Tab.4 Sample Set of Iterative Experiments

依实验的迭代顺序把最优分类点处IF统计量按照正常状态（1-4）、内圈故障（5-8）、滚珠故障（9-13）和外圈故障（14-16）的顺序进行排列，其中序号1、5、9、14对应的为四类状态的训练样本IF值统计量，其余为测试样本IF值统计量，如图8所示。可以看出各类轴承状态的IF区分度很大，且各次实验的结果一直保持在稳定的数值区间中，进一步证明了轴承故障诊断系统符合工业实用的需要。

5 结论

这里提出了应用于张量模型初始特征最优分类点的区分度函数，设计了一种利用HOSVD和区分度函数的轴承故障诊断系统。为验证轴承故障诊断系统的有效性，分别利用千鹏故障平台和凯斯西储大学的轴承数据进行实验。实验结果表明，故障诊断系统能够获得比MLSLLE算法中更优越的初始特征，且与最优分类点区分度函数结合后能够取得故障分类的高准确率。