王宏超, 陈 进, 董广明
(上海交通大学机械系统与振动国家重点实验室, 上海 200240)
由于各种干扰源和噪声的影响,来自现场传感器的轴承故障信号非常复杂,直接对这些信号进行故障特征提取是非常困难的。盲源提取(Blind source extraction,BSE)技术的出现为这个问题的解决提供了新的技术和手段。BSE技术是近10年来新兴并引起广泛重视的一类信号处理方法,有时也叫盲信号抽取,它源于盲源分离(Blind source separation, BSS)技术。BSS技术在语音、通讯和医学工程等领域都有着成功的应用[1~3],但其在旋转机械故障诊断中的应用有着较大的局限性,其原因可归结为如下[4]:瞬时混合模型的幅值不确定性和顺序不确定性;适用于振动信号的多通道卷积混合盲分离算法迄今仍没有令人满意的解决方案;源物理相关性的不确定性;源数的估计存在偶然性;盲源分离要求系统具有可逆性,但实际工程中很多机械系统是不可逆的。BSE相对于BSS有以下优点[5]:充分利用被提取对象的先验信息;仅仅提取感兴趣的信号;无需估计源的数目;计算量小;特别适用于源数较多而感兴趣信号较少的情况(如机械系统)。综上所述,BSE在某种程度上弥补了BSS的缺点。在既有的BSE技术中,约束独立成分分析最具代表性,其最先由Wei Lu等人于2000年提出[5]。随后Wei Lu分别将CICA用于图像盲分离与医学信号盲分离[6,7],相对于传统独立成分分析(Independent component analysis, ICA)方法取得了不错的效果。Zhang Z L将CICA方法加以改进,提出了基于形态的CICA方法并将其用于微弱脑电信号的盲源提取[8],通过仿真合成信号和实际脑电信号验证了所提出方法的正确性及实际应用价值。在旋转机械故障诊断方面,CICA也取得了一定的应用[9],在文献[9]中Wang Z Y等将CICA用于滚动轴承的故障诊断中,并通过仿真信号和滚动轴承全寿命加速疲劳实验信号验证了CICA在旋转机械盲信号提取中的有效性。然而,上述CICA在盲源提取中的成功应用都需要以下两个前提条件:用于构建参考信号的目标源信号基本周期的精确估计;正确参考信号的构建。如文献[9]所验证,不正确参考信号的选取会造成误判或漏判的结果。此外,用于构建参考信号的目标源信号基本周期的估计值也会对提取结果造成很大的影响。然而在滚动轴承的实际工程应用中由于安装精度、制造误差及滚子相对于滚道的随机滑动等影响,都会造成目标源故障信号的实际基本周期与理论计算基本周期的偏差。这样就限制了CICA方法在旋转机械故障诊断中的应用。在文献[10]中,Barros等提出了一种基于最佳估计延时的能快速提取周期信号的盲提取方法。但所述方法本质上还是以目标信号的基本周期或基本周期的整数倍作为最佳延时。本文首先用基于二阶统计量的自相关函数算法估计目标源信号基本周期[10],再与基于高阶统计量的固定点算法相结合的盲提取方法成功用于复杂运行环境下滚动轴承故障信号的盲提取[11,12],此种方法在某种意义上来讲是对文献[10]所述方法的改进。此外,此种方法相对于CICA方法具有较强的鲁棒性:只需要大致估计目标信号的基本周期并将其作为所述方法相关步骤的基本周期即可。同文献[9],将所述方法用于滚动轴承全寿命加速疲劳实验故障信号的盲提取,并与文献[9]CICA在滚动轴承全寿命加速疲劳故障信号盲提取中的应用作以比较,得出前者相对于后者具有更强的目标源信号基本周期估计误差的容错性。
假设有一多维观测信号矢量x(t)可表示为
x(t)=As(t)+n(t)
(1)
式中s(t)=[s1(t),s2(t),…,sn(t)]T是n×1满足统计独立性的源信号矢量,x(t)=[x1(t),x2(t),…,xm(t)]T是m×1的混合信号矢量,矩阵A={aij}是由混合参数aij构成的m×n阶的混合矩阵。n(t)=[n1(t),n2(t),…nm(t)]T是与源信号矢量统计独立的加性噪声矢量。
(2)
式中 矩阵C=(ci,j)n×n被称为混合-分离复合矩阵,该矩阵的每一行与每一列有且仅有一个非零元素。瞬时混合模型的混合和分离过程如图1所示。
图1 瞬时混合和分离模型示意图
根据传感器数m和源数n的关系,瞬时ICA模型可划分为3种:
当m=n时,此时混合矩阵A是方阵,称为平方ICA模型;
当m 当m>n时,称为超定ICA模型。 本文是基于平方ICA模型展开讨论的。 此外,为提高计算速度及简化本文算法,首先对观测信号x(t)进行白化处理,即 (3) 式中V称为白化矩阵。以下所述所有观测信号x(t)均为白化信号。 对于要提取的目标源信号si,假设对一特定整数τ*有以下关系式成立 (4) 式中sj为其他非目标源信号,k代表时间,τ*即为最佳时间延时[10]。 J(W)=E{y(k)y(k-τ*)}=WTE{x(k)x(k-τ*)T}W (5) 在约束条件‖W‖=1下,最大化公式(5)所示的目标函数就可以提取出目标函数si。因为对于目标函数si,其关于延时为τ*的延时自相关函数是一个大于0的正值,而其他非目标源函数sj关于τ*的延时自相关函数都将为0。 忽略Rx(τ*)与Rx(τ*)T的微小差距,利用标准梯度算法,可得出文献[10]的算法 (6) 式中Rx(τ*)=E{x(k)x(k-τ*)T}。尽管文献[11]所述方法具有算法简洁、计算速度快等优点,但是其具有以下缺点: (1)即使目标源信号si关于τ*的时延自相关函数E{si(k)si(k-τ*)}>0,但不能保证其他所有非目标源信号sj关于τ*的时延自相关E{sj(k)sj(k-τ*)}=0(j≠i); (2)即使所有源信号严格独立不相关,但是在实际计算中由于用源信号有限点数的算术平均来代替源信号的数学期望,就会造成源信号的自相关计算值不为0[15,16],即:即使有下式成立 (7) 但也极有可能下式成立 (8) 综上所述,文献[10]所述基于二阶统计量的盲提取方法受上述两点限制,其盲提取性能的稳定性及精确性就受到了很大的影响,甚至会提取出非目标源信号。由于本文的研究是基于瞬时线性ICA基础上的,也就是假设源信号在物理意义上是相互独立的,所以可以用高于二阶统计量的更高阶的统计量去改善文献[10]所述的方法。以下是所述方法的介绍。 所用的方法大概分为3个步骤,如图2所示。 图2 所述方法流程图 各个步骤的具体计算过程如下: (1)目标源信号基本周期的估计:当滚动轴承发生故障时,其故障基本周期可以由理论计算频率的倒数得出; (2)初始目标源信号及分离矩阵计算: 式(5)所示的目标函数 (9) (10) (11) 由式(10)可得出式(11)的相应算法如下 (12) (3)初始分离矩阵W的进一步优化: (13) 为提高算法对野点和冲击噪声的鲁棒性,可以采用文献[11]的改进固定点算法 (14) 式(13)或(14)收敛时,即可得到最佳分离矩阵W。以下的仿真及实验结果均取式(13),(14)二者之间效果最好的。 图3 原始信号 图4 观测信号 图5 不同盲提取方法得到的结果 图3为4个信号,其中图3中的s1为准周期信号,s2为正弦信号,s3为余弦信号,s4为随机高斯白噪声。在Matlab中任意产生一随机4×4阶的矩阵A对4个信号进行混合。图4为混合后的观测信号。提取目标信号是s1。图5为用不同方法得到盲抽取结果。其中图5中的y1为文献[10]所述方法的盲抽取结果,得到了错误的抽取结果;图5中y2为文献[15]所述方法的提取结果;图5中的y3为本文所述方法的提取结果。虽然直观上看y2与y3要都为正确的目标源信号s1。可以用下式来对提取结果y2和y3的盲抽取精度进行比较 (15) 将s1和y2,s1和y3分别代入式(15)得到的PI值分别为22.4和47.6 dB,由此,所述方法相对于文献[15]所提出的方法具有更高的提取精度。 图6 基于不同参考信号的CICA抽取结果 图6是CICA方法用不同参考信号对观测信号的盲抽取结果。其中参考信号r1为正确的参考信号,将其与观测信号一起输入到CICA算法中,y11为抽取的结果。从y11可以看出CICA在正确输入参考信号的情况下可以很好地抽取出目标信号s1;不改变参考信号r1的方波周期,而改变方波的宽度即得到参考信号r2,将r2与观测信号一起输入到CICA算法中,y22为抽取的结果,可见抽取不出目标信号。同样不改变r1方波的宽度,而改变方波的周期得到参考信号r3,将r3与观测信号一起输入到CICA算法中,y33为抽取的结果,可见仍抽取不出目标信号。由此,CICA算法的能否成功应用在某种程度上取决于参考信号的能否正确设计,其使用受到了很大的限制。而分别将r1,r3的周期作为所述方法的基本估计周期(或基本周期的整数倍),均能提取出目标源信号s1。限于篇幅限制在此不再给出抽取结果。 同样采用文献[9]中的滚动轴承全寿命加速疲劳实验数据。实验具体实施步骤及数据采集参数均可见参考文献[9]。选用和文献[9]同一个故障轴承进行分析,文献[9]给出了故障轴承的全寿命周期的均方根值(Root mean square, RMS)值曲线图。本文给出故障轴承全寿命周期的峭度指标曲线图,如图7所示,因为峭度指标对滚动轴承发生故障时的冲击特征更为敏感。 图7 实验轴承全寿命周期的峭度值 对第2 304组数据进行分析。将3个加速度传感器采集到的信号作为观测信号,其时域图如图8所示(注:文献[9]的时域图纵坐标单位为加速度,而本文的纵坐标单位为电压)。 图9为图8所示信号相对应的包络解调谱,从图9中均无法很好地提取出滚动轴承内圈故障通过频率fi=246 Hz及其谐频:图9(a)虽然能大致提取出内圈故障通过频率,但其谐频及调制频率即转频提取效果并不好;图9(b),(c)的提取结果会造成误判或漏判。 将本文所述方法用于上述观测信号,提取出的信号的时域图及其包络解调谱图如图10所示,由包络解调谱可以看出其很好的提取出内圈故障通过频率fi=246 Hz及其谐频。在用所述方法对观测信号进行处理时,估计目标信号的基本周期设置为T=(1/fi)×fs=104(用点数表示,可参考文献[9])。 图11是基于不同参考信号的CICA盲提取结果,其中图11(a)的r1是正确参考信号:周期为T=104(采用点数表示,可参考文献[9]);方波宽度为20(点数),图11(b)为盲抽取的信号,图11(c)为11(b)的包络解调谱,由此看见,CICA很好地提取了内圈故障通过频率;改变图11(a)的周期,即T′=106,方波宽度不变,得到参考信号r2如图11(d);图11(e)为盲提取结果的时域图,图11(f)为图11(e)的包络解调谱,提取不出内圈故障频率;同样,不改变图11(a)的周期,而改变方波宽度(方波宽度由20设置为30)即得到参考信号r3如图11(g);图11(h),(i)分别为盲提取结果的时域图及包络解调谱,同样提取不出内圈故障特征频率。可见,CICA能否正确地提取故障特征频率,对精确参考信号的设计有着苛刻的要求(参考信号的周期及方波形状的微小改变都会对CICA提取结果造成很大的影响)。在工程实际应用中,由于设计精度、安装误差及滚动体的随机滑动都会造成滚动轴承理论计算故障频率与实际故障频率的误差。而本文所述的方法只需要大致估计目标源信号的基本周期即可。图12为将T=106作为目标源信号的估计周期,用本文所述方法得到盲提取信号的时域图及其相应的包络解调谱,仍可以提取出内圈故障特征频率。说明了本文所述方法相对于CICA方法有较强的目标源信号周期误差的容错性。 图8 实验观测信号 图9 实验观测信号的包络解调谱 图10 所述方法的盲提取结果 图13为文献[15]所述方法的提取效果,从其包络解调谱中无法提取中滚动轴承内圈故障通过频率。 图11 用不同参考信号CICA的提取结果 图12 用所述方法以T=106作为估计基本周期的提取结果 图13 用文献[15]所述方法的盲提取结果 本文将一种盲提取方法用于复杂运行环境下滚动轴承故障信号的盲提取中,取得了不错的结果。并通过仿真验证了所述方法相对于其他盲提取方法有较高的精度;通过仿真和实验验证了所述方法相对于CICA方法具有目标源信号基本周期误差容错性的优点;此外,所述方法只需要估计目标源信号的基本周期即可,而CICA方法不但需要目标源信号基本周期的精确估计,而且还需要构建准确的参考信号,其使用相对于所述方法具有较大的局限性。 参考文献: [1] Hyvarinen A, Oja E. Independent component analysis:algorithms and applications[J]. Neural Networks, 2000, 12:411—430. [2] Choi S J, Cichocki A, Park H M, et al. Blind source separation and independent component analysis:A review[J]. Neural Information Processing-letters and Reviews, 2005, 6(1):1—57. [3] James C J, Gibson O J. Temporally constrained ICA:An application to artifact rejection in electromagnetic brain signal analysis[J]. IEEE Transaction on Biomedical Engineering, 2003, 50(1):1 108—1 116. [4] Antoni J. Blind separation of vibration components:rinciples and demonstrations[J]. Mechanical Systems and Signal Processing 2005, 19:1 166—1 180. [5] Lu W, Rajapakse J C. Constrained independent component analysis, in:Advance in Neural Information Processing Systems, vol.13(NIPS 2000)[M]. MIT Press, Cambridge, MA, 2000:570—576. [6] Lu W, Rajapakse J C. Approach and applications of constrained ICA[J]. IEEE Transaction on Neural Networks, 2005, 16(1):203—212. [7] Lu W, Rajapakse J C. ICA with reference[J]. Neurocomputing, 2006, 69:2 244—2 257. [8] Zhang Z L. Morphologically constrained ICA for extracting weak temporally correlated signals[J]. Neurocomputing, 2008, 71:1 669—1 679. [9] Wang Z Y, Chen J, Dong G M, et al. Constrained independent component analysis and its application to machine fault diagnosis[J]. Mechanical Systems and Signal Processing, 2011, 25:2 501—2 512. [10] Barros A K, Cichocki A. Extraction of specific signals with temporal structure[J]. Neural Computation, 2001, 13(9):1 995—2 003. [11] Cichocki A, Amari S I. Adaptive Blind Signal and Image Processing. Learning Algorithms and Applications[M]. New York:John Wiley & Sons,2002. [12] Hyvarinen A, Oja E. A fast fixed-point algorithm for independent component analysis[J]. Neural Computation, 1997, 9(7):1 483—1 492. [13] Sabri K, Badaoui M E, Guillet F, et al. A frequency domain-based approach for blind MIMO system identification using second-order cyclic statistics[J]. Signal Processing, 2009, 89:77—86. [14] Belouchrani A, Meraim K A, Cardoso J F, et al. A blind source separation technique using second-order statistics[J]. IEEE Transactions on Signal Processing, 1997, 45(2):434—444. [15] Zhang Z L, Zhang Y. Robust extraction of specific signals with temporal structure[J]. Neurocomputing, 2006, 69:888—893. [16] Bermejo S. Finite sample effects in higher order statistics contrast functions for sequential blind source separation[J]. IEEE Signal Processing Letters,2005, 12(6):481—484.2 方法的提出
3 仿 真
4 实 验
5 结 论