改进PCA算法及其在转子特征提取中的应用*

李伟光, 郭明军, 杨期江, 赵学智, 李国臣

(1.华南理工大学机械与汽车工程学院 广州，510640) (2.广州航海学院轮机工程学院 广州，510725) (3.东莞职业技术学院实训中心 东莞，523808)

引 言

开展大型滑动轴承试验台的设计和特性研究，对各种类型、结构和参数的大型滑动轴承的研究具有重要意义[1]。 笔者所在团队自主研发了一种大型变支点滑动轴承试验台，采用一种基于特征值差分谱理论的PCA信号处理方法来提纯其轴心轨迹，进而识别其工作状态。

滑动轴承的失效形式多种多样，常用的诊断方法包括油样分析、声发射检测和振动分析等，其中通过监测转子位移信号的振动分析方法最为简单有效[2]。采集转子同一轴截面相互垂直布置的两个位移传感器的振动位移信号，将其合成轴心轨迹，轴心轨迹的辨识主要研究轴心轨迹图像的特征提取和识别问题，其前提也需要有清晰的轴心轨迹图作为依据[3]。然而实际采集的振动信号通常会受到干扰因素的影响，所以，需要对原始信号合成的轴心轨迹进行提纯，其本质是对原始振动信号进行降噪处理后合成轴心轨迹[4]。常用的提纯方法包括：数字或模拟低通滤波法[3]、小波变换和小波包变换[5]、粒子群算法、EMD降噪[6]及形态滤波[7]等。上述降噪方法一般会存在相位偏移、降噪畸变大及频带选择不明确等问题[4]，因此，改进已有方法或寻求新的方法都是一种可行的方案。

近年来，主成分分析技术被广泛应用于消噪、故障诊断、特征提取及模式识别等领域[8-11]。刘永斌等[8]根据累积贡献率选取主元个数，对双层圆柱壳体机械噪声数据进行降维。尚前明等[9]将PCA直接应用到船舶柴油机的故障监测中，通过PCA实现柴油机热工参数的降维，从而准确识别出柴油机的异常状态。Li等[10]采用PCA方法对风力机滚动轴承故障信号的特征矩阵进行降维，通过滤除特征中的冗余信息使得支持向量机(support vector machine,简称SVM)的分类取到了更加有效的信息，精度得以提高并极大减少了计算量。Seghouane等[11]针对主载荷矢量的密集结构使得主载荷矢量的降维难以解释的问题，提出了一种自适应块稀疏PCA方法。

主成分分析的关键是有效主元个数的确定，目前多数研究都是根据累计贡献率[8-11]取定某个阈值的方法选择有效主元的个数，阈值越大主成分的个数就越多，保留的信息量就越大，从而可能保留的噪声成分也越多[12]。赵学智等[13]提出用奇异值差分谱来描述有用信号和噪声的奇异值性质差异性，对于去除直流分量的信号，根据其差分谱首个峰值的位置可以自动选择有效奇异值的个数，并能有效抑制噪声成分的影响。将差分谱方法与 PCA技术融合, 用于描述主成分和次要成分的协方差矩阵特征值的差异性，以及探讨奇异值与特征值之间的关联等两个问题目前还鲜有文献报道。笔者结合装配了新型结构滑动轴承的试验对这两个问题开展研究, 并将具体的 PCA 算法应用到大型滑动轴承试验台的转子轴心轨迹的提纯上，提纯效果优于传统PCA算法。

1 Hankel矩阵方式下PCA信号处理原理

1.1 PCA信号分解原理

主成分分析[14]是指用k个n维的新变量y1,y2,…,yk来线性表示n维初始变量x1,x2,…,xm(m≥k),使得新变量的方差最大或降维损失最小,即有

(1)

其中：系数αi=(αi1,αi2,…,αim)T(i=1,2,…,k)为协方差矩阵C中降序排列的第i特征值λi对应的特征向量，且αi满足

(2)

将式(1)改写为分量形式有

(3)

其中：yi∈R1×n,i=1,2,…,k。

协方差矩阵C的特征方程为

Cαi=λiαi

(4)

协方差矩阵C的计算式为

C=

(5)

其中：cov(xi,xj)=E[(xi-E(xi))(xj-E(xj))T]。

1.2 Hankel矩阵构造

在式(5)中矩阵X采用Hankel矩阵，可利用一维零均值的离散信号a=[a(1),a(2),…,a(N)] ，按照以下方法构造Hankel矩阵

(6)

矩阵X称为重构吸引子轨迹矩阵，也称为Hankel矩阵。

1.3 奇异值与特征值的关系推导

对任意实矩阵X,其奇异值分解表示[12]为

X=UDVT

(7)

其中：U=(u1,u2,…,um)∈Rm×m，V=(v1,v2,…,vn)∈Rn×n分别为左奇异矩阵和右奇异矩阵且都为正交矩阵,其中ui∈Rm×1，vi∈Rn×1分别称为左奇异向量和右奇异向量；D=diag(σ1,σ2,…,σr)为对角矩阵,其元素为按降序排列的奇异值,即σ1≥σ2≥…≥σr≥0,r=min(m,n)为矩阵X的秩。

(8)

根据式(7)可得

(9)

根据U，V的正交性有

VTV=E

(10)

其中：E为单位矩阵。

又因为

(11)

把式(9～11)代入式(8)得

(12)

对比式(8)与式(12)可得

(13)

式(13)表明, 奇异值和特征值之间是一种平方关系。

1.4 有效主元分量选取

(14)

其中:bi形成的序列B=(b1,b2,…,br-1)为特征值差分谱。

差分谱反映了相邻两个特征值之间的变化趋势，如果特征值在某个位置变化较大，在差分谱中将出现一个峰值，而在所有的峰值当中必然存在一个最大值bk。此时，k处发生一个最大突变，这种突变反映了有效主成分和次要成分的相关性的差异，代表主成分和次要成分的分界，据此确定的k值即为主成分的个数。

为了研究奇异值差分谱与特征值差分谱的联系，将式(13)代入式(14)可得

(15)

从式(15)可以看出，特征值差分谱可以通过奇异值的平方差计算出来 ，这种平方运算可以将奇异值差分谱放大，从而使主成分峰值特征更加突出。从这一点而言，特征值差分谱是奇异值差分谱理论的升华。

1.5 分量信号重构

求解式(4)可得协方差矩阵的特征向量αi和λi，将所有特征值降序排列有λ1≥λ2≥…≥λm。用特征向量αi同时左乘式(3)等号两边，并只选择前k个主成分相加可得

(16)

(17)

2 基于差分谱理论的PCA算法

2.1 算法流程

根据特征值差分谱的性质，提出一种基于特征值差分谱理论的PCA算法(简称差分PCA算法)，具体步骤如下：

1) 对经零均值化处理的一维离散信号a，按式(6)构造Hankel矩阵X；

2.2 信号分析实例

为了验证文中算法，任意构造一个信号

x(t)=1.5sin(10πt+0.5)+2.0sin(20πt+0.8)

(18)

以1 024 Hz的采样频率对信号x(t)进行离散化,并叠加信噪比为0.15 dB的高斯白噪声,结果如图1所示。

图1 原始含噪信号Fig.1 Original noisy signal

图2 使用差分谱提取不同分量信号的过程Fig.2 Process to extract signals using difference spectrum

为了弄清楚差分谱不同峰值位置代表的含义，笔者逐个进行分析。从图2(a)可知，差分谱的第1个峰值出现在第2个序号，因此利用前2个分量进行重构，结果如图2(c)所示，图中虚线为理想信号中的2.0sin(20πt+0.8)成分，该提纯信号代表原始信号x(t)中的2.0sin(20πt+0.8)信号分量。另外已知，最大峰值序号对应的前4个分量提纯的信号包含了x(t)中的所有频率成分，那么利用最大峰值与其前面的第1峰值之间的第3，4个分量重构的信号必然为x(t)中的1.5sin(10πt+0.5)；利用第3，4个分量重构的信号如图2(d)所示，其中虚线为理想的1.5sin(10πt+0.5)，说明该提纯信号为原信号x(t)中的1.5sin(10πt+0.5)，提纯结果与前述分析一致。

综上所述，根据特征值差分谱最大峰值位置能够可靠地提取纯信号中的主成分，通过不同峰值之间的分量信号的组合可以提取出不同的频率成分，且频率对应的幅值越大其对应的峰值位置越靠前。

3 试验分析

3.1 试验装置简介

试验采用自主研发的滑动轴承试验台装置[1]，如图 3 所示，主要包括润滑系统、驱动系统、机械结

图3 试验台实物图Fig.3 Rotor test-bed

构等部分组成。

3.2 位移传感器测点布置

在转子两端垂线斜45°位置各布置一个电涡流传感器，分别标记为di(i=1,2,3,4)，对应的信号为Di(i=1,2,3,4)，其中靠驱动装置端安装的传感器如图4所示。试验所用的便携式LMS SCADASD多功能数据采集系统如图5所示。

图4 驱动装置端位移传感器测点布置Fig.4 Location of displacement sensors near the driver

图5 LMS数据采集系统Fig.5 LMS data acquisition system

3.3 振动信号采集

笔者对某次试验中(采样频率为 1 024 Hz)试验台主轴转速分别为 450 r/min 和 1 080 r/min 的驱动端主轴振动位移信号(D1和D2)进行分析，其时域波形及频谱图分别如图 6、图7 所示。在两种转速工况下，信号的频谱成分在整个频带内存在随机噪声和系统工作时激发出的高频谐波分量。

图6 信号去除直流分量的结果(450 r/min)Fig.6 Signals with D C components removed at 450 r/min

图7 信号去除直流分量的结果(1 080 r/min)Fig.7 Signals with D C components removed at 1 080 r/min

4 轴心轨迹提纯

4.1 采用差分PCA算法提纯轴心轨迹

利用处理过的D1，D2两组信号合成的转子轴心轨迹如图8所示，但不能清晰表明轴心轨迹。

图8 原始轴心轨迹图Fig.8 Original orbits of the rotor

采用笔者提出的差分PCA算法,分别利用处理过的振动位移信号D1和D2构造513×512维的Hankel矩阵Xi(i=1,2,3,4)；然后，按式(5)分别计算协方差矩阵Ci(i=1,2,3,4)，并对其进行特征值分解，经换算得到原始信号的特征值及其差分谱曲线如图9所示(E1，E2分别代表信号D1，D2的特征值)，图中只给出了前30个值的结果。由图9可知：前4个特征值较大，第5个及其以后的特征值相对较小; 差分谱中前两峰值的幅值较大。这与转子旋转时的1倍频与2倍频幅值最为突出、且2倍频与1倍频幅值接近[16]的结论吻合。

图9 特征值及其差分谱曲线Fig.9 Eigenvalue curves and their difference spectrum

若选取差分谱中第2个峰值对应的前4个分量，按照上述步骤对信号进行同样的处理，结果如图11所示。由图11可知，此时提取的是原始信号的1X(7.5 Hz)及3X(22.5 Hz)成分，而与图10的结果比较可知，前2个峰值之间的第3，4分量对应的是3X(22.5 Hz)，从而验证了2.2节中的“通过差分谱不同峰值之间的分量信号组合可以提取出不同的频率成分”的结论。

图11 前4个分量信号的重构结果Fig.11 Reconstruction results of the first four components

利用图10，11中的信号合成的轴心轨迹分别如图12(a，b)所示。由图12(a)可知轴心轨迹为长短轴相差不大的椭圆形，说明转子存在轻微不平衡现象。由图12(b)可知轴心轨迹为花瓣形，其特征是各小瓣到瓣心的距离r较大且互相分离。

图12 工况1的提纯轴心轨迹Fig.12 The purified axis trajectories at condition 1

图13 前4个分量信号的重构结果Fig.13 Reconstruction results of the leading four components

由图13中的信号合成的轴心轨迹如图14(a)所示，该图为各花瓣相交的花瓣形，但瓣心到花瓣的距离r是减小的，由此表明相对于工况1，工况2转子的轴心轨迹波动量更大，但是轴心振动幅值是减小的。本试验台装配的是一种不同于常规圆轴承的特殊结构滑动轴承，由于本研究主要侧重于该轴承结构转子的振动轴心轨迹特征的提取，对低转速下产生花瓣型轴心轨迹的机理需要进一步深入研究。若选取差分谱中前2个峰值之间对应的第3，4个分量，按照上述步骤对信号进行同样的处理，然后合成的轴心轨迹如图14(b)所示，轴心轨迹的椭圆形长短轴与工况1相近，表明与工况1的动不平衡量基本一致。两组试验未调整过轴承与转子系统的任何参数，其动不平衡量不会发生改变，但不同组试验测试采集的数据存在误差，会导致微小的差异。

图14 工况2的提纯轴心轨迹Fig.14 The purified axis trajectories at condition 2

4.2 与传统PCA算法提纯效果比较

传统PCA算法具有依赖于研究者的信号分析经验的缺陷，其根据累积贡献率[8-11]来确定主成分的个数，累积贡献率的表达式为

(19)

其中：Li为累计贡献率；k为主特征值个数；m为有效特征值个数；δ为某个选定的阈值。

现以工况1的信号为例，利用传统PCA方法进行处理。得到的特征值分布图与图9(a，b)中的一致。由图可知，当i≤14时有λi>1，故取m=14。为了便于分析，将k=1,2,…，6时，相应的特征值与累积贡献率列于表1中，此时分别利用提纯得到前k个分量合成的轴心轨迹如图15所示。

表1 前6个特征值及其贡献率

图15 使用传统PCA算法提纯的轴心轨迹Fig.15 Purified axis orbits with traditional PCA algorithm

由图15可知，贡献率不同，提纯得到的轴心轨迹形状也不同。对比图12与15可知，图15(b，d)分别与图12(a，b)基本一致，说明只有当k=2，4时，根据贡献率确定的分量才是原始信号的主成分，而当k取其他值时，提纯的结果并不准确。

在实际应用当中，累积贡献率的阈值并不容易确定，严重依赖于研究者的信号分析经验，这使得其应用具有一定的局限性。而笔者提出的特征值差分谱方法，一方面可以根据差分谱最大峰值位置自动筛选主成分的个数；另一方面，通过不同谱峰之间的分量信号的组合可以提取出不同的频率成分。这些是传统PCA方法所不具备的。

5 结束语

笔者提出利用协方差矩阵特征值差分谱的概念来描述有效主成分与次要成分的特征值差异性，根据差分谱的最大峰值位置可自动选择有效主成分的个数，从而解决了主成分分析的关键问题。研究了Hankel矩阵方式下PCA方法处理信号的原理，并提出一种基于特征值差分谱理论的PCA算法，通过仿真信号验证了该算法的有效性。研究表明：通过差分谱不同谱峰之间的分量信号的组合可以提取出不同的频率成分。将笔者提出的差分PCA算法用于装配了特殊结构滑动轴承的转子轴心轨迹提纯，效果优于传统PCA算法，可直观表明转子的动不平衡特征状态。