相空間稀疏化的信號壓縮感知與重構方法*

溫廣瑞, 欒日維, 任延暉, 馬再超

(西安交通大學機械工程學院 西安，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)；中央高?；究蒲袠I務費專項資金資助項目

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