相空间稀疏化的信号压缩感知与重构方法*

温广瑞, 栾日维, 任延晖, 马再超

(西安交通大学机械工程学院 西安，710049)



相空间稀疏化的信号压缩感知与重构方法*

温广瑞, 栾日维, 任延晖, 马再超

(西安交通大学机械工程学院 西安，710049)

针对旋转机械振动信号受强噪声干扰导致传统FFT频域稀疏性差，难以进行正交匹配重构的问题，提出了相空间稀疏化结合正交匹配追踪(orthogonal matching pursuit, 简称OMP)的信号压缩感知(compressed sensing, 简称CS)方法。首先,对信号进行相空间重构(phase space reconstruction, 简称PSR)，并采用主分量分析(principal component analysis, 简称PCA)提取主要分量和重构信号，以提高信号的频域稀疏性；然后,采用随机高斯矩阵测量及压缩频域稀疏性得到优化的信号;最后,采用正交匹配追踪算法重构信号。仿真信号和转子典型不对中信号的分析结果表明，该方法可以提高受强噪声干扰的振动信号在频域内的稀疏性，实现转子振动信号的有效压缩和准确重构。

压缩感知； 相空间重构； 主分量分析； 正交匹配

引 言

旋转机械(如离心式压缩机、风力发电机组、大型汽轮发电机组等)是现代大型企业中的核心设备，在石油、化工等领域得到了广泛应用[1-3]。为保障转子、轴承和齿轮箱等核心部件长期安全稳定运行，实时监测其运行状态尤为重要。然而，当前旋转机械信息的采集通常表现出测点多、传感器多、数据量大的特点，对信号的传输带来了巨大压力[4-6]，迫切需求一种有效的数据压缩方法，以保证信号的高效传输。传统的压缩方法是将完整的信号样本投影到某个稀疏域(如频域)上，再对投影向量最重要的分量及其位置进行存储和传输，浪费了计算资源与存储空间。

文献[7-8]提出的压缩感知技术，突破了传统采样所必须遵循的乃奎斯特采样定理，为数据压缩提供了一种新的思路，并已广泛应用于图像编码、模式识别、无线通讯等方面。该理论指出，如果某个信号是稀疏的或在某个变换域内是稀疏的，那么就可以用一个与变换基不相关的测量矩阵将该信号投影到一个低维空间上，得到一组低维的测量样本，然后再通过求解一个关于最小范数的优化问题的方法，就可以从这些少量的测量样本中以高概率重构出原信号。但实际情况中，设备运行信息通常伴有强噪声，如果选择频域作为稀疏域，传统的FFT方法会导致信号频域内的稀疏性变差，因此需要一种有效的信号稀疏表示方法来满足压缩感知理论对信号稀疏性的要求。主分量分析[9]是一种经典的特征提取方法，其通过线性变换来实现特征压缩，达到用尽可能少的维数最大限度地表示原始特征信息的目的[10-11]，同时相空间重构[12-13]可以将信号从一维升到等效的多维信号空间，为主分量分析提供有效的数据支持，故将基于相空间表达模型的主分量分解引入到压缩感知理论的稀疏表示部分，提高信息在频域内的稀疏性，则有望达到压缩感知理论对信号稀疏性的要求。

基于以上分析，针对旋转机械转子振动信号信噪比低导致的频域内稀疏性差的问题，通过引入相空间表达模型的主分量分解方法，提高该类信号频域内的稀疏性，进而采用压缩感知理论对优化信号进行测量压缩与重构，并通过仿真信号和实验转子振动信号进行分析验证，以期为旋转机械转子振动信号的压缩与重构提供一种研究思路。

1 基本概念简介

1.1 相空间重构

给定一个时间序列x(n),n=1,2,…，N，根据Takens[13]提出的嵌入定理，选择某些固定的时间延迟上的数据点作为新的维数，可以重构出一个“等价”的相空间，从而将其从一维时间序列扩展成等价的多维时间序列，如式(1)所示

(1)

其中：t为延迟时间;m为嵌入维数;N为数据长度。

1.2 主分量分析

主分量分析是将含有噪声成分的混合信号中几个主要分量提取出来，而将不重要或者是噪声成分予以删除，从而用几个主要分量来表征混合信号中的绝大部分信息的一种统计分析方法。

实现该方法的关键步骤是求解混合信号中各分量的投影方向ν及各分量在该投影方向上所包含的信息量大小，计算式为

(2)

其中：C为由混合信号构造的协方差矩阵；λ为特征值，表示投影方向上的分量所含信息量的大小。

1.3 压缩感知

信号的稀疏表示是压缩感知的基础和前提。一般来说，一个时间序列x∈Rn×1，投影到Rn空间内的一组正交基Ψ上，得到一组表示系数s

(3)

如果s中非零元素个数为k，且k≪n，则该信号在变换域Ψ上是稀疏的。由于转子具有旋转运行的特点，其振动信号通常会表现出强烈的周期性，因此常选择三角函数基作为变换基，将信号由时域转换为频域，实现转子振动信号的稀疏表示。

随机高斯矩阵(矩阵元素满足N(0,1/m)的高斯分布)几乎与任何稀疏表示的基都不相关，且需要的测量值数目较少，是压缩感知理论中常用的随机测量矩阵之一。由此，信号的稀疏表示与测量压缩组成了压缩感知的线性测量过程

(4)

正交匹配追踪算法是一种贪婪迭代算法，即在每一次的迭代过程中，从过完备原子库里(即测量矩阵)选择与信号最匹配的原子来进行稀疏逼近，求出余量并保证余量与之前的每一个分量均正交。正交匹配追踪算法简单，运行速度快，精度高，是压缩感知理论中常用的信号重构方法之一。

2 基于相空间稀疏化方法的转子振动信号压缩感知

由于相空间重构结合主分量分析，可以提取信号中的主要成分，故将其代替传统FFT的频域稀疏化方法，从而达到使信号频域内更为稀疏的效果。此时，再利用合适的随机矩阵和重构算法，对频域内稀疏性得到优化的信号进行测量压缩和重构。

2.1 基于相空间稀疏化方法的特征提取

首先根据1.1节所述，对于长度为N的一维时间序列样本x= (x1,x2, …,xN)T进行相空间重构。由于Hanke1矩阵近似于方阵，滤波效果最佳[14]，信号在频域内的稀疏性也会达到最优，因此根据时间序列的数据长度N选择合适的参数t和m，可得

(5)

其中：m+n-1=N且n≥2，m≥2。

文献[10]指出，当m足够大时，对矩阵X进行零均值化处理，得到矩阵X及其协方差矩阵，即

(6)

解关于协方差矩阵C的特征方程，如式(2)所示。把特征值按由大到小的顺序排列。定义前p个特征值的累计贡献率，如式(7)所示，即

(7)

其中：η的大小用来衡量特征提取后信息的保留程度，越大则保留的信息越多。

从式(7)中选择对应特征值较大的若干投影方向组成变换矩阵D，将原始样本进行投影变换，得到提取特征后的p维样本，即

(8)

(9)

2.2 信号的测量压缩与重构

式(9)所示的恢复信号，其在频域内的稀疏性已经得到优化，可以直接进行信号的稀疏表示。对优化后的信号进行傅里叶变换，将信号转换至频域内，得到一组稀疏的傅里叶系数s，再进行信号的测量压缩。笔者选择M行N列的随机高斯矩阵作为测量矩阵Φ，根据式(4)对傅里叶系数s或者优化后的信号进行测量压缩，得到一组长度远小于原始数据长度的离散压缩样本y。最后，根据离散的压缩样本y对信号进行重构。本研究采用正交匹配追踪算法来重构信号。根据文献[15]，正交匹配追踪算法的步骤如下：

1) 初始化残差r0=y,指标集Λ0=Ø,迭代计数t=1;

2) 寻到满足下述最优化问题的指标λt

(10)

3) 扩充指标集和矩阵Λt=Λt-1∪{λt}及Φt=[Φt-1φλt]，Φ0为空矩阵；

4) 求解最小二乘问题

(11)

5) 计算新的信号估计和残差

(12)

6)t=t+1，若t

7) 重构信号x*的非零值指标为Λm中的元素，x*中第λj个元素的值等于xt的第j个元素。

2.3 振动信号的压缩感知流程设计

根据2.1，2.2节阐述的基于相空间稀疏化方法的振动信号压缩感知，可以归纳出其具体步骤如下：

1) 根据时间序列的数据长度，将转子振动信号按照式(5)构建Hankel矩阵；

2) 按照式(6)，将Hankel矩阵构建协方差矩阵，按照式(2)和(7)～(9)，采用PCA分析协方差矩阵，提取信号中的主要分量，滤掉噪声，得到频域内稀疏性得以优化的信号样本；

3) 对优化后的信号样本进行傅里叶变换，得到频域内的一组稀疏系数s；

4) 按照式(4)，选取随机高斯矩阵作为测量矩阵，对稀疏系数s进行压缩，得到压缩样本y；

5) 根据压缩样本y，采用正交匹配追踪算法对信号进行重构。

对上述步骤进行总结，得到基于相空间稀疏化方法的振动信号压缩感知流程框架，如图1所示。

图1 转子振动信号压缩感知流程
Fig．1 Flow chart of compressed sensing for rotor vibration signal

2.4 仿真信号验证

设置仿真信号的采样频率为1 024 Hz，采样长度为1 024，信号的频率组成主要有50，100 Hz，幅值分别为1，2，并在此基础上加入信噪比为4 dB的白噪声后，其信号形态及频谱如图2和图3所示。

图2 加噪仿真信号Fig．2 Simulation signal with noise

图3 加噪仿真信号频谱
Fig．3 The spectrum of simulation signal with noise

从图2可以看出，信号受到设置的白噪声干扰，各极值点变化较大，且无法察觉到周期分量的存在；从图3看出，信号的主要频率成分为50，100 Hz，与设置频率相对应，且白噪声表现出宽频特性，信号的稀疏性不佳。

首先根据2.3节归纳的转子振动信号的压缩感知流程步骤1，设置m=513，n=512，首先生成513×512的Hankel矩阵如下

(13)

对Hankel矩阵进行标准化处理并根据式(6)生成512×512的协方差矩阵；采用PCA分析协方差矩阵，依次按照式(2)和(7)～(9)，提取信号中的主要分量，滤掉噪声，得到频域内稀疏性得以优化的信号样本及其频谱，如图4和图5所示。

图4 稀疏性优化信号Fig．4 The signal with optimum sparsity in frequency domain

图5 稀疏性优化信号频谱
Fig．5 The spectrum of signal with optimum sparsity in frequency domain

对比图2和图4可以看出，优化后的信号可以消除仿真信号中的噪声，保留有用分量，突出仿真信号中的周期性规律；由图5可以看出，重构信号的频谱也集中在50和100Hz两个部分，没有了宽频白噪声干扰。因此，基于Hankel矩阵的PCA降噪方法能够很好地除仿真信号中混有的强噪声成分，在保留信号主要信息的基础上对信号频域内的稀疏性进行优化。

然后，再根据2.3节中归纳的步骤3对频域内稀疏性得到优化的样本进行傅里叶变换，得到一组傅里叶系数，采用随机高斯矩阵作为测量矩阵对这组傅里叶系数进行测量压缩，得到一组离散的复数压缩样本，对其取模后如图6所示。

图6 压缩样本(M=30)
Fig．6 Compressed sample(M=30)

最后，根据得到的压缩样本及随机高斯矩阵，同样依照2.2节中OMP算法的详细步骤对信号进行重构，得到重构后的信号及其频谱，如图7，图8所示。

图7 重构优化信号时域波形
Fig．7 The time domain waveform of reconstructed signal with optimum sparsity

图8 重构优化信号频谱
Fig．8 The spectrum of reconstructed signal with optimum sparsity

对比图4和图7可以看出，重构信号(图7)同样表现出了明显的周期性规律，且幅值与优化信号(图4)相吻合；对比图5和图8可以看出，重构信号的频谱集中在50Hz和100Hz处，幅值与优化信号基本吻合；而且，根据压缩样本长度M、原信号样本长度N和压缩率k计算公式

k=M/N*100%

(14)

由式(14)可知，该方法将原数据长度为1 024的信号压缩成长度为30的压缩样本，压缩率为2.93%。由此表明，基于相空间稀疏化的压缩感知方法是合理有效的。

3 应用情况分析

针对转子振动信号受到噪声干扰，频域内稀疏性不佳的问题，笔者在本特利转子试验台上模拟了不对中故障。试验台主要由转子系统和信号采集系统组成，如图9所示。转子系统主要有电机、转速调节器、转子、质量盘和支撑部件；信号采集系统包括三组位移传感器(两组用于采集转子振动信号，一组用于采集转速信息)、数据采集模块和计算机。

图9 本特利转子试验台Fig．9 Bently rotor test bench

故障信号的采样频率为1 024 Hz，数据长度为1 024，转子转速约3 000 r/min。该信号的时域波形及其频谱如图10、图11所示。

由图10可以看出，该信号的周期性受到强噪声干扰，极值点变化剧烈；从图11可以看出，频域的基频、2倍频占主要成分，噪声表现出宽频特性，信号的稀疏性不佳。

图10 转子不对中信号时域波形Fig．10 The time domain waveform of rotor signal with misalignment

图11 转子不对中信号频谱Fig．11 The spectrum of rotor signal with misalignment

如果直接对该信号进行傅里叶变换，得到一组傅里叶系数，然后根据2.3节中步骤4，采用随机高斯矩阵作为测量矩阵对这组傅里叶系数进行测量压缩，得到一组离散的复数压缩样本，对其取模后如图12所示。

图12 转子不对中信号压缩样本(M=50)Fig.12 Compressed sample of rotor signal with misalignment(M=50)

然后，根据得到的压缩样本及随机高斯矩阵，依照2.2节中OMP算法的详细步骤对信号进行重构，得到重构后的信号及其频谱，如图13、图14所示。

图13 重构信号时域波形Fig．13 The time domain waveform of reconstructed signal

图14 重构信号频谱Fig.14 The spectrum of reconstructed signal

对比图11与图14可以发现，重构信号的频率成分为50，100和 160 Hz，160Hz为无关频率成分，且50和100 Hz的频率成分在重构前后存在明显的误差。

如果采用笔者所提出的方法，则首先需要依据2.3节归纳的转子振动信号的压缩感知方法，设置m=513，n=512，生成513×512的Hankel矩阵如下

(15)

由对Hankel矩阵标准化处理并生成512×512的协方差矩阵；采用PCA分析协方差矩阵，依次按照式(3)和(7)～ (9)，提取信号中的主要分量，滤掉噪声，得到频域内稀疏性得以优化的信号样本及其频谱，如图15和图16所示。

图15 稀疏性优化后转子振动信号Fig．15 The rotor vibration signal with optimum sparsity

图16 稀疏性优化后转子振动信号频谱Fig．16 The spectrum of rotor vibration signal with optimum sparsity

由图15可以看出，优化后的信号表现出了明显的周期性规律；由图16可以看出，重构信号的频谱也集中在基频、2倍频位置，幅值准确，没有了宽频白噪声干扰。因此，基于Hankel矩阵的PCA降噪方法能够很好地去除转子振动信号中混有的白噪声成分，在保留信号主要信息的基础上对信号频域内的稀疏性进行优化。

然后，根据2.3节中归纳的步骤3对频域内稀疏性得到优化的样本进行傅里叶变换，得到一组傅里叶系数，进而根据步骤4，采用随机高斯矩阵作为测量矩阵对这组傅里叶系数进行测量压缩，得到一组离散的复数压缩样本，对其取模后如图17所示。

图17 稀疏性优化后信号压缩样本(M=50)Fig．17 Compressed sample of signal with optimum sparsity(M=50)

最后，根据得到的压缩样本及随机高斯矩阵，依照2.2节中OMP算法的详细步骤对信号进行重构，得到重构后的信号及其频谱，如图18、图19所示。

图18 稀疏性优化后重构转子振动信号时域波形Fig．18 The time domain waveform of reconstructed rotor vibration signal with optimum sparsity

图19 稀疏性优化后重构转子振动信号频谱Fig．19 The spectrum of reconstructed rotor vibration signal with optimum sparsity

对比图15和图18可以看出，重构信号基本再现了原信号，表现出了明显的周期性，与优化信号幅值相差不大；对比图16和图19可以看出，重构信号的频谱同样集中在基频、2倍频位置，与优化信号幅值基本一样；而且，该方法将原数据长度为1 024的信号压缩成长度为50的压缩样本，根据式(14)可得压缩率为4.88%。由此表明,基于相空间稀疏化的压缩感知方法能够对受到噪声干扰的转子振动信号进行测量压缩并准确重构。

4 结束语

针对旋转机械振动信号在采集过程中会受到强噪声的干扰，导致其频域内的稀疏性较差，直接应用压缩感知方法对其进行测量压缩并重构后，重构信号与原始信号会存在较大误差的问题，提出了相空间稀疏化的振动信号正交匹配重构方法。通过引入相空间表达模型的主分量分解方法，提高该类信号频域内的稀疏性，然后采用压缩感知理论对频域稀疏性得到优化的信号进行测量压缩，进而实现振动信号正交匹配重构。通过仿真信号和转子实际振动信号的分析与验证，表明该方法可为旋转机械各类受噪声干扰振动信号的压缩与重构提供一种新的研究思路。

[1] 冯志鹏，宋希庚，薛冬新. 旋转机械振动故障诊断理论与技术进展综述[J]. 振动与冲击，2001，20(4)：36-39．

Feng Zhipeng, Song Xigeng, Xue Dongxin. Survey of vibration fault diagnosis of rotational machinery [J]. Journal of Vibration and Shock, 2001, 20(4): 36-39．(in Chinese)

[2] 杜永祚，秦志英. 旋转机械动态信号全息谱分析[J]. 振动、测试与诊断，2002，22(2)：81-88．

Du Yongzuo, Qin Zhiying. Holo-spectrum analysis of rotating machinery dynamic signals [J]. Journal of Vibration, Measurement & Diagnosis, 2002, 22(2): 81-88．(in Chinese)

[3] 孟宗，闫晓丽，王亚超. 基于LMD和HMM的旋转机械故障诊断[J]. 中国机械工程，2014，25(21)：2942-2951．

Meng Zong , Yan Xiaoli, Wang Yachao. Rotating machinery fault diagnosis based on local mean decomposition and hidden markov mode [J]. China Mechanical Engineering, 2014, 25(21):2942-2951．(in Chinese)

[4] 杜金榜，王跃科，潘仲明，等. 旋转机械振动数据压缩及语音压缩技术的应用研究[J]. 计算机测量与控制，2006，14(12)：1594-1634．

Du Jinbang, Wang Yueke, Pan Zhongming, et al. Vibration data compression of rotary machines and application of speech compression technology [J]. Computer Measurement & Control, 2006, 14(12): 1594-1634．(in Chinese)

[5] 翁浩，高金吉. 旋转机械振动信号压缩小波基优化选取方法[J]. 振动、测试与诊断，2013，33(3)：437-444．

Weng Hao, Gao Jinji. Selection method of optimum wavelet base in vibration-signal compression of rotating machinery [J]. Journal of Vibration, Measurement & Diagnosis, 2013, 33(3): 437-444．(in Chinese)

[6] 徐敏强，张嘉钟，张国斌，等. 基于小波变换的旋转机械振动信号数据压缩方法的研究[J]. 振动工程学报，2000，13(4)：531-536.

Xu Minqiang, Zhang Jiazhong, Zhang Guobin, et al. Method of data compressing for rotating machinery vibration signal based on wavelet transform [J]. Journal of Vibration Engineering, 2000, 13(4): 531-536.(in Chinese)

[7] Donoho D L. Compressed sensing[J]. IEEE Transactions on Information Theory， 2006,52(4):1289-1306.

[8] Candè E J, Wakin M B. An introduction to compressive sampling[J]. IEEE Signal Processing Magazine, 2008,25(2):2l-30.

[9] Broomhead D S, King G P. Extracting qualitative dynamics from experimental data, Physica D[J]. Nonlinear Phenomena, 1986, 20(2-3)： 217-236 ．

[10] 李巍华，史铁林，杨叔子. 基于非线性判别分析的故障分类方法研究[J]. 振动工程学报，2005，18(2)：133-138．

Li Weihua, Shi Tielin, Yang Shuzi. Mechanical fault classification using nonlinear discriminant analysis [J]. Journal of Vibration Engineering, 2005, 18(2): 133-138. (in Chinese)

[11] 屈梁生，张西宁，沈玉娣. 机械故障诊断理论与方法[M]. 西安：西安交通大学出版社，2009:170-171．

[12] Pockard N H, Crutchfield J P, Farmer J D， et al. Geometry from a time series [J]. Phys Rev Lett, 1980, 45(3): 712-716.

[13] Takens F. Detecting strange attractors in turbulence [J]. Lecture Notes in Mathematics, 1981, 898(2): 361-381．

[14] 吴浩浩，罗志增. 基于构造Hankel矩阵的SVD陷波方法[J]. 计算机应用研究，2010，27(12)：4514-4516．

Wu Haohao, Luo Zhizeng. Signal notch method based on hankel matrix and SVD [J].Application Research of Computers, 2010, 27(12): 4514-4516．(in Chinese)

[15] Tropp J A, Gilbert A C. Signal recovery from random measurements via orthogonal matching pursuit [J]. IEEE Transactions on Information Theory, 2007, 53(12): 4655-4666．

*国家自然科学基金资助项目(51365051,51421004);教育部新世纪优秀人才支持计划资助项目(NCET-13-0461)；中央高校基本科研业务费专项资金资助项目

2015-04-06；

2015-08-14

10.16450/j.cnki.issn.1004-6801.2017.02.003

TH17； TP18

温广瑞，男，1976年7月生，博士、教授、博士生导师。主要研究方向为机械运行状态故障诊断及性能维护、现场动平衡理论及方法、远程及现场监测与诊断系统开发。曾发表《A new transient field balancing method of a rotor system based on empirical mode decomposition》(《Journal of Vibroengineering》2013,Vol.15,No.3)等论文。 E-mail：grwen@mail.xjtu.edu.cn