面向滚动轴承故障诊断的VBIVA算法

Variational Bayesian independent vector analysis algorithm for rolling bearing fault diagnosis

摘要

在实际工程中，采集到的滚动轴承故障信号往往来自多个源。多个故障信号在传播路径中相互耦合形成了复合故障信号，使故障诊断问题变得更加复杂。如果直接对复合信号进行分析，那么提取到的故障特征中往往存在多源的故障频率，导致无法正确判断故障出现的位置。针对这一问题提出了变分贝叶斯独立向量分析（variational Bayesian independent vector analysis，VBIVA）算法，并将该算法应用于故障诊断。通过与独立向量分析（independent vector analysis，IVA）算法以及变分贝叶斯独立分量分析（variational Bayesian independent component analysis，VBICA）算法的仿真对比，证明VBIVA算法有效地解决了复合故障信号的盲源分离及故障诊断问题。

关键词

滚动轴承 / 故障诊断 / 盲源分离 / 独立向量分析 / 变分贝叶斯

Key words

rolling bearing / fault diagnosis / blind source separation / independent vector analysis / variational Bayesian

滚动轴承^［1］是旋转机械的重要部件，常安装在驱动电机^［2］、齿轮箱^［3］等部位，且在汽车、航空航天等领域有着广泛的应用。滚动轴承起着减少摩擦和支撑载荷的作用，它的健康程度决定了整个设备能否平稳安全运行。因此，对滚动轴承进行故障诊断至关重要^［4］，这是保证旋转机械安全运行、降低维护成本的关键^［5］。

在过去的二十年里，学者们提出过许多基于物理变量的轴承故障诊断技术，这些物理变量包括磁通^［6］、电机电流^［7］、振动信号^［8］、声信号^［9］等。其中，从机械振动信号中提取故障特征是故障诊断的一种主要方法^［10-11］。在实际工程中，一个轴承故障往往会导致轴承出现新的故障，最终采集到的信号是含有噪声的复合故障信号^［12-13］，这使故障诊断的难度加大。因此，对复合故障信号进行盲源分离是解决问题的新思路。盲源分离算法是一种信号分离技术，仅凭借已知的混合信号就可以估计出未知的源信号。早在20世纪80年代，Herault等^［14］首次提出了BSS理论。近年来，盲源分离算法在故障诊断上得到了广泛的应用。Lu等^［15］提出了一种适用于轴承故障诊断的谐波向量分析（harmonic vector analysis，HVA）方法。Cherrad等^［16］使用了盲源分离方法分析轴承振动信号。文献［17-18］使用独立分量分析（independent component analysis，ICA）方法对振动信号进行轴承故障诊断。王立炜等^［19］将FastICA用于复杂路径下的外圈滚动轴承故障诊断。但是，滚动轴承故障信号作为一种典型的非平稳信号，很容易被周围噪声淹没，对盲源分离和故障诊断有很大影响。此外，多源的故障信号与传感器之间复杂的传输路径也使问题更加复杂化。

因此，在这样复杂的情况下，传统的盲源分离算法和特征提取方法的诊断效果会受影响。为了解决复合故障信号难以诊断的问题，本文提出了一种变分贝叶斯独立向量分析（VBIVA）算法，并将其应用于滚动轴承复合故障诊断的研究。首先，利用最小化辅助函数的IVA方法对复合故障信号进行初步分离；然后，利用变分贝叶斯理论推断出合适的参数，构建源信号模型、噪声模型及混合矩阵模型，对分离信号进一步更新；最后，对分离故障信号进行希尔伯特变换，提取相应的故障特征并完成故障诊断。在实验中对模拟复合滚动轴承振动信号进行盲源分离及故障诊断，证明了VBIVA方法的有效性。

1 理论基础

1.1 盲源分离模型

盲源分离是指根据已知的观测信号来估计未知的源信号。在传感器与滚动轴承振动源信号之间存在着复杂的路径，可以使用线性混合模型作为振动信号的混合模型，如式（1）所示。

\mathit{x} ( t ) = \mathit{A} ( t ) \mathit{s} ( t ) + \mathit{n} ( t ) （1）

式中： \mathit{x} ( t ) = [ x_{\mathrm{1}} ( t ) , x_{\mathrm{2}} ( t ) , \dots , x_M {(t)]}^{\mathrm{T}}为 M维的观测复合滚动轴承故障信号； \mathit{s} ( t ) = [ s_{\mathrm{1}} ( t ) , s_{\mathrm{2}} ( t ) , \dots , s_N {(t)]}^{\mathrm{T}}为 N维的源信号； \mathit{n} ( t ) = [ n_{\mathrm{1}} ( t ) , n_{\mathrm{2}} ( t ) , \dots , n_M {(t)]}^{\mathrm{T}}为 M维的噪声信号； \mathit{A} ( t )为 N \times M维的混合矩阵。为了方便书写，将 \mathit{x} \left(t\right) \mathrm{、} \mathit{s} \left(t\right) \mathrm{、} \mathit{y} \left(t\right) \mathrm{、} \mathit{A} \left(t\right)简化 \mathit{x} \mathrm{、} \mathit{s} \mathrm{、} \mathit{y} \mathrm{和} \mathit{A}。

1.2 变分贝叶斯学习

变分贝叶斯学习通过确定统计模型的后验分布来对未知样本做出预测。在混合故障信号 \mathit{x}的指导下，可以建立模型并通过计算后验分布来恢复源信号。将负自由能作为变分贝叶斯学习的目标函数，可以表示为

F ( \mathit{x} ) = {〈\mathrm{l}\mathrm{o}{\mathrm{g}}_{\mathrm{2}}p(\mathit{x},\Psi )〉}_{p^{\text{'}}(\Psi )} + H [ p^{\text{'}} ( \Psi ) ] （2）

式中： \Psi为模型 \Omega中参数的集合； p^{\text{'}} ( \Psi )近似于 p ( \mathit{x} , \Psi )的后验分布； {〈q〉}_p符号是指在概率密度函数 p前提下 q的数学期望值。KL散度用于计算两个概率密度间差异，它是非负的，由此可知 \mathrm{l} \mathrm{o} {\mathrm{g}}_{\mathrm{2}} p ( \mathit{x} ) \ge F ( \mathit{x} )。为了使KL散度最小化，就需要得到最大的 F ( \mathit{x} )。

2 变分贝叶斯独立向量分析

在故障诊断问题中，由于轴承的信号存在随机性^［20］并包含噪声，给复合故障信号的盲源分离及故障诊断造成一定困难。为了解决这些问题，本文提出了变分贝叶斯独立向量分析（VBIVA）算法，主要包含IVA算法和变分贝叶斯学习两部分。

2.1 IVA算法

传统的IVA算法假设故障信号间相互独立，其优化目标就是最大化分离信号之间的统计独立性。在使用共轭梯度法对目标函数进行更新的过程中需要设定步长，其大小决定了收敛速度和稳定性。较大的步长会有更快的收敛速度，但也可能导致收敛结果发散。与传统的IVA算法不同，本文设计了一种基于最大期望（expectation maximization，EM）^［21］算法思想的辅助函数，避免在迭代过程中选择到不合适的更新步长。通过最小化辅助函数的方法，得到具有最强独立性的解混矩阵。设计目标函数 J ( \tilde{\mathit{W}} )满足式（3）要求。

J ( \tilde{\mathit{W}} ) = \mathrm{m} \mathrm{i} {\mathrm{n}}_{\mathit{Q}} \Phi ( \mathit{W} , Q ) （3）

式中： \mathit{Q}为辅助变量； \mathit{W}为分离矩阵； \Phi ( \mathit{W} , \mathit{Q} )为目标函数 J ( \tilde{\mathit{W}} )的辅助函数，其表达式为

\Phi ( \mathit{W} , \mathit{Q} ) = \frac{\mathrm{1}}{\mathrm{2}} {\displaystyle \sum _{n=\mathrm{1}}^N} {\mathit{w}}_n^H{\mathit{Q}}_n{\mathit{w}}_n - \mathrm{l} \mathrm{o} {\mathrm{g}}_{\mathrm{2}} \left[\mathrm{d}\mathrm{e}\mathrm{t}(\mathit{W})\right] + C （4）

{\mathit{Q}}_n = \mathit{E} \left[{\phi }_{n}x{x}^H\right] （5）

式中： \mathrm{d} \mathrm{e} \mathrm{t} \left[\bullet \right]为矩阵行列式运算符号； H为共轭转置； Q_n为第 n行的辅助变量； {\mathit{w}}_n为解混矩阵 \mathit{W}的第 n行向量； C为常数；权重因子 {\phi }_n用于最小化辅助函数，从而间接最小化目标函数，它的表达式为

{\phi }_n = \mathrm{l} \mathrm{o} {\mathrm{g}}_{\mathrm{2}} ( {‖y_n‖}_{\mathrm{2}} ) = \frac{\mathrm{1}}{\mathrm{2}} \mathrm{l} \mathrm{o} {\mathrm{g}}_{\mathrm{2}} ( {\displaystyle \sum _{t=\mathrm{1}}^L} y_n ( t ) + \mathrm{1} + \epsilon ) （6）

式中： {‖\bullet ‖}_{\mathrm{2}}为 L \mathrm{2}范数；为了保证权重因子 {\phi }_n不为零，加入了一个数值极小的常数 \epsilon。解混矩阵 \mathit{W}的更新规则如式（7）、（8）所示

{\mathit{w}}_n^{{}^{(i)}} = ( {\mathit{W}}^{(i)} {\mathit{Q}}_n^{{}^{(i)}} {)}^{-\mathrm{1}} I_n （7）

{\mathit{w}}_n^{{}^{(i+\mathrm{1})}} = \frac{{\mathit{w}}_n^{{}^{(i)}}}{\sqrt[]{({\mathit{w}}_n^{{}^{(i)}}{)}^H{Q}_n^{{}^{(i)}}{\mathit{w}}_n^{{}^{(i)}}}} （8）

式中： \mathit{I}为单位矩阵；符号 {\bullet }^{(i)}为第 i次迭代结果。迭代终止条件是 。

2.2 变分贝叶斯学习

混合故障信号经过基于辅助函数的IVA算法初步分离后，得到分离信号 y、解混矩阵 \mathit{W}及其相对应的混合矩阵 \mathit{A}。本文引入了变分贝叶斯学习，在重建信号过程中对分离信号进一步更新。

源信号模型中的每个源信号 s_i有相对应的高斯成分，所有参数都包含 t个高斯成分。那么源信号的概率密度分布可以写为

p ( s | q , \theta ) = {\displaystyle \prod _{i=\mathrm{1}}^N} {\displaystyle \sum _{t=\mathrm{1}}^{T_i}}{\pi }_{i,t}N\left(s_i|{\mu }_{i,t},{\beta }_{i,t}\right) （9）

式中：源信号 s_i个数为 N，那么就有 T_i = N个高斯成分； \theta = [ {\theta }_{\mathrm{1}} , {\theta }_{\mathrm{2}} , \dots , {\theta }_N ]为高斯模型的参数， {\theta }_i = [ {\pi }_i , {\mu }_i , {\beta }_i ]； {\pi }_i为混合权重，其先验概率分布取狄利克雷分布 p ( \pi ) = {\displaystyle \prod _{i=\mathrm{1}}^N} D({\pi }_i,{\lambda }_i)； {\mu }_i为 s_i的期望，其先验概率分布 p ( \mu ) = {\displaystyle \prod _{i=\mathrm{1}}^N} {\displaystyle \prod _{t=\mathrm{1}}^{T_i}}N({\pi }_{i,t}|m_i,{\tau }_i)； {\beta }_i为精度值，定义呈Gamma分布的先验概率分布 p ( \beta ) = {\displaystyle \prod _{i=\mathrm{1}}^N} {\displaystyle \prod _{t=\mathrm{1}}^{T_i}}\Gamma ({\beta }_{i,t}|b_{{\beta }_i},c_{{\beta }_i})； m \mathrm{、} \tau \mathrm{、} b \mathrm{、} c为设定的超参数。

假设噪声 n的均值为零并且服从高斯分布，那么其概率分布密度可表示为

p ( n ) = N ( n \left|\mathrm{0},\Lambda ;M\right. ) （10）

式中： \Lambda为噪声精度，定义其先验分布 p ( \Lambda )为Gamma分布， p ( \Lambda ) = {\displaystyle \prod _{j=\mathrm{1}}^M} \Gamma ({\Lambda }_j|{b_{\Lambda }}_j,{c_{\Lambda }}_j)。

定义混合矩阵 \mathit{A}的先验分布 p ( \mathit{A} )服从高斯分布，其表达式表示为

p ( \mathit{A} ) = {\displaystyle \prod _{i=\mathrm{1}}^N} {\displaystyle \prod _{j=\mathrm{1}}^M}N({\mathit{A}}_{ji}\left|\mathrm{0},{\alpha }_{ji}\right.) （11）

式中： {\alpha }_{ji}为混合矩阵 \mathit{A}中混合系数 a_{ji}的方差，根据超参数 b_a和 c_a更新。

2.3 VBIVA算法推理过程

在VBIVA模型中，需要先初始化模型，设定 {\lambda }_i = \mathrm{5} , m = \mathrm{0} , \tau = \mathrm{1} , b_{\beta } = \mathrm{1} \mathrm{000} , c_{\beta } = \mathrm{0.001} , b_a = \mathrm{1} \mathrm{000} , c_a = \mathrm{0.001} , b_{\Lambda } = \mathrm{1} \mathrm{000} , c_{\Lambda } = \mathrm{0.001}。通过IVA算法得出初步的混合矩阵和分离信号，对初步分离的信号进行k-means聚类，得到模型的权重、均值和精度。

初始化模型后，以最大化目标函数 F ( x )为目标开始变分贝叶斯学习。最大化 F ( x )也就是最大化 p^{\text{'}} ( \Psi )，通过贝叶斯推理可以确定后验概率 p ( \Psi | x , \Omega )。使用变分近似法求得目标函数 F ( \mathit{x} )，如式（12）所示。

F ( \mathit{x} ) = {〈\mathrm{l}\mathrm{o}{\mathrm{g}}_{\mathrm{2}}p(\mathit{x},\Psi |\Omega )〉}_{p\mathrm{\text{'}}(\Psi \theta )} + H \left[p\mathrm{\text{'}}(\Psi )\right] （12）

设定更新阈值 \epsilon = \mathrm{0.000} \mathrm{1}，当更新后的目标函数 F ( x )小于这个阈值时，变分贝叶斯学习停止。在学习结束后，根据模型中获得的混合比例权重 {\pi }_i和均值 {\mu }_i恢复源信号。算法结构如图1所示。

3 复合故障信号诊断实验分析

滚动轴承通常由内圈、外圈、滚动体和保持架4部分组成，本文主要研究了内圈、外圈及滚动体的复合故障。使用仿真信号进行实验，仿真轴承的轴承节径为29.05 \mathrm{m} \mathrm{m}，滚动体直径为6.75 \mathrm{m} \mathrm{m}，滚动体个数为8，接触角设置为0 °，转速为900 \mathrm{r} / \mathrm{m} \mathrm{i} \mathrm{n}。通过计算得到模拟的故障脉冲信号的内圈故障频率 f_I = \mathrm{73.907} \mathrm{H} \mathrm{z}，外圈故障频率 f_O = \mathrm{46.037} \mathrm{H} \mathrm{z}，滚动子故障频率 f_B = \mathrm{30.52} \mathrm{H} \mathrm{z}。

观测到的复合故障信号经过希尔伯特变换（Hilbert transform，HT）后可提取出故障特征，以倍频出现的故障频率作为故障诊断的依据。如图2所示，在没有噪声的情况下，直接提取出内圈故障频率和外圈故障频率，在混合故障信号2和3中，虽然可以找到滚动子故障和外圈故障的频率，但它们的幅值明显低于内圈故障频率，被淹没在其中。这意味着传统的包络谱诊断技术不能准确地提取出复合故障特征频率，而本文提出的方法可以简单快速地对复合故障信号进行诊断。

为了证明本文方法的有效性，将本文提出的方法分别和变分独立分量分析（VBICA）^［22］、独立向量分析（IVA） ^［23］进行比较。图3是IVA算法^［23］的分离结果，提取到外圈的故障特征中包含了内圈的故障特征，同理内圈故障及滚动子信号也掺杂了其他故障特征，即复合故障信号没有被完全分离。而本文提出的VBIVA算法分离出的故障信号可以清楚地辨别故障特征，如图4所示。

含有噪声的复合故障信号更接近真实的滚动轴承振动信号。使用VBICA算法对其进行复分析，结果如图5~6所示，可以观察到故障频率周期性的变化，显然复合信号彻底分离。然而，分离故障信号中仍然有噪声的存在。如图7~8、表1所示，本文提出的VBIVA算法不仅可以准确地提取故障特征，还可以分离出干净的故障信号。对于复合故障信号的分离和诊断，本文提出的VBIVA算法效果更好。本文使用了均方误差（mean square error，MSE）作为评价指标，对分离效果进行比较。 \mathrm{M} \mathrm{S} \mathrm{E} = E ( {(y-x)}^{\mathrm{2}} )，其值越小表示结果越好。

4 结论

本文提出了一种基于辅助函数的IVA和变分贝叶斯学习的VBIVA算法，首先计算解混矩阵的辅助函数，省去了选择步长以及更新步长的步骤，得到初步分离的混合故障信号及解混矩阵；然后，为了进一步获得接近源信号的分离信号，引入了变分贝叶斯学习，通过更新源信号模型、噪声模型及分离矩阵模型，解决了信号受噪声干扰问题；最后，本文对分离的故障信号进行希尔伯特变换，提取故障特征，对故障进行诊断。仿真实验结果表明，本文提出的VBIVA算法有效地解决了轴承复合故障的盲源分离问题，特别是克服了工程环境中不同振动信号混合、传播路径复杂和噪声信号淹没等问题。

参考文献

[1] Wang Y， Tang B P， Qin Y，et al.Rolling bearing fault detection of civil aircraft engine based on adaptive estimation of instantaneous angular speed［J］.IEEE Transactions on Industrial Informatics，2020，16（7）：4938-4948.

[2] Ren B Y， Yang M， Chai N，et al.Fault diagnosis of motor bearing based on speed signal kurtosis spectrum analysis［C］//2019 22nd International Conference on Electrical Machines and Systems （ICEMS）.Harbin：IEEE，2019：1-6.

[3] Han Q K， Wang T Y， Ding Z，et al.Magnetic equivalent modeling of stator currents for localized fault detection of planetary gearboxes coupled to electric motors［J］.IEEE Transactions on Industrial Electronics，2021，68（3）：2575-2586.

[4] Han P P， He C B， Lu S L，et al.Rolling bearing fault diagnosis based on variational mode decomposition and cyclostationary analysis［C］//2021 International Conference on Sensing，Measurement & Data Analytics in the era of Artificial Intelligence （ICSMD）.Nanjing：IEEE，2021：1-6.

[5] 冯葛豪，胥永刚，孙国栋， 谱聚类傅里叶分解及其在轴承故障诊断中应用［J/OL］.轴承，2023，1-10［2023-09-25］.

[6] Zamudio R I， Osornio R R A， Antonino D J A，et al.Magnetic flux analysis for the condition monitoring of electric machines：a review［J］.IEEE Transactions on Industrial Informatics，2022，18（5）：2895-2908.

[7] Hoang D T， Kang H J.A motor current signal-based bearing fault diagnosis using deep learning and information fusion［J］.IEEE Transactions on Instrumentation and Measurement，2020，69（6）：3325-3333.

[8] Gangsar P， Tiwari R.Signal based condition monitoring techniques for fault detection and diagnosis of induction motors：a state-of-the-art review［J］.Mechanical Systems and Signal Processing，2020，144（3）：106908.

[9] Glowacz A， Tadeusiewicz R， Legutko S，et al.Fault diagnosis of angle grinders and electric impact drills using acoustic signals［J］.Applied Acoustics，2021，179（2）：108070.

[10] Wang W B， Lei Y G， Chen Z X，et al.A diagnosis method of cleanliness defect in rolling element bearing manufacturing［C］//2022 Global Reliabi-lity and Prognostics and Health Management （PHM-Yantai）.Yantai：IEEE，2022：1-6.

[11] Lv M R， Liu X， Li X.Application of pattern filter method on fault diagnosis of rolling bearing vibration signals［C］//2019 3rd International Confe-rence on Electronic Information Technology and Computer Engineering （EITCE）.Xiamen：IEEE，2019：523-527.

[12] Jia S Y， Qi Y S， Li Y T.A composite fault diagnosis method of gearbox based on an enhanced deconvolution algorithm［C］//2021 CAA Symposium on Fault Detection，Supervision and Safety for Technical Processes.Chengdu：IEEE，2021：1-6.

[13] Qin A S， Mao H L， Hu Q，et al.Bearing fault diagnosis method based on ensemble composite multi-scale dispersion entropy and density peaks clustering［J］.IEEE Access，2021，12（9）：24373-24389.

[14] Herault J， Jutten C.Space or time adaptive signal processing by neural network models［C］//AIP Conference Proceedings 151 on Neural Networks for Computing.Snowbird：ACM，1987：206-211.

[15] Lu J T， Yin Q T， Li S M，et al.Rolling bearing composite fault diagnosis method based on harmonic vector analysis［C］//2022 Global Reliability and Prognostics and Health Management （PHM-Yantai）.Yantai：IEEE，2022：1-6.

[16] Cherrad M L， Bendjama H， Fortaki T.Vibration analysis for defective bearings by blind source separation［C］//2021 International Congress of Advanced Technology and Engineering （ICOTEN）.Taiz：IEEE，2021：1-4.

[17] 吴金钟，艾延廷，陈英涛，基于盲源分离技术的航空发动机轴承故障诊断［J］.滨州学院学报，2022，38（2）：27-35.

[18] 陈剑，刘圆圆，黄凯旋，基于奇异值分解和独立分量分析的滚动轴承故障诊断方法［J］.计量学报，2022，43（6）：777-785.

[19] 王立炜，栾孝驰，沙云东，基于快速ICA的复杂路径下滚棒轴承外圈划伤故障特征提取［J］.机械设计与制造，2021（12）：77-81，87.

[20] Wu C B.Application of adaptive length morphological filtering in gearbox bearing fault recognition in quantum domain［C］//2020 International Conference on Computer Information and Big Data Applications （CIBDA）.Guiyang：IEEE，2020：459-464.

[21] Przyborowski M， Ślęzak D.Approximation of the expectation-maximization algorithm for Gaussian mixture models on big data［C］//2022 IEEE International Conference on Big Data （Big Data）.Osaka：IEEE，2022：6256-6260.

[22] Choudrey R A， Roberts S J.Variational mixture of Bayesian independent component analyzers［J］.Neural Computation，2003，15（1）：213-252.

[23] Kim T， Eltoft T， Lee T W.Independent vector analysis：an extension of ICA to multivariate components［C］//Independent Component Analysis and Blind Signal Separation.Berlin：Springer，2006：165-172.

基金资助

辽宁省自然科学基金(2022-MS-299)

航空科学基金(201933054002)

辽宁省教育厅项目(LJKMZ20220529)