改進(jìn)停止閾值準(zhǔn)則的振動(dòng)信號(hào)降噪方法研究

楊坤，王桂寶，王楠，王磊，李一飛

（1.陜西理工大學(xué)機(jī)械工程學(xué)院，陜西漢中 723001；2.陜西省工業(yè)自動(dòng)化重點(diǎn)實(shí)驗(yàn)室，陜西漢中 723001）

0 前言

機(jī)械振動(dòng)信號(hào)蘊(yùn)含著設(shè)備運(yùn)行過程中的重要信息，在監(jiān)測旋轉(zhuǎn)機(jī)械的運(yùn)行狀態(tài)時(shí)，通常需要采集設(shè)備機(jī)組主要部件的振動(dòng)信號(hào)數(shù)據(jù)，如滾動(dòng)軸承、齒輪箱、加速器等。并且根據(jù)這些數(shù)據(jù)，運(yùn)用信號(hào)變換等理論把信號(hào)轉(zhuǎn)換到頻域內(nèi)來處理，繼而進(jìn)行設(shè)備故障的監(jiān)測與診斷［1］。在旋轉(zhuǎn)機(jī)械設(shè)備的監(jiān)測環(huán)境中，設(shè)備內(nèi)部及外界都會(huì)產(chǎn)生較強(qiáng)的噪聲，這些噪聲信號(hào)不可避免地會(huì)摻雜到純凈的振動(dòng)信號(hào)中，使得振動(dòng)信號(hào)特征提取難度加大，繼而會(huì)對(duì)信號(hào)分析處理結(jié)果產(chǎn)生負(fù)面影響［2］。

近年來，壓縮感知理論為信號(hào)降噪這一方向提供了新思路。呂麒鵬等［3］在提取滾動(dòng)軸承原始振動(dòng)信號(hào)的故障頻率時(shí)，由于該信號(hào)存在著強(qiáng)背景噪聲，嚴(yán)重影響提取故障頻率的效率，甚至還會(huì)提取到錯(cuò)誤的故障頻率，在使用壓縮感知理論對(duì)原始振動(dòng)信號(hào)進(jìn)行降噪處理后，故障提取效率、準(zhǔn)確性以及計(jì)算效率得到大幅提高。郭亮等人［4］利用Lasso 算法對(duì)壓縮信號(hào)進(jìn)行稀疏重構(gòu)，恢復(fù)出振動(dòng)信號(hào)，對(duì)于時(shí)變信號(hào)具有很好的壓縮能力，并且能夠有效地去除信號(hào)中的噪聲成分。

在信號(hào)重構(gòu)階段使用最廣泛的是正交匹配追蹤算法，但是該算法的迭代誤差會(huì)在迭代循環(huán)中繼續(xù)傳導(dǎo)，導(dǎo)致重構(gòu)階段所得信號(hào)的精度比較低［5］。此外，在信號(hào)降噪的迭代過程中，觀測向量內(nèi)積序列的主要成分為純凈信號(hào)的信息。經(jīng)過多次迭代后，純凈信號(hào)成分逐漸被提取出來。迭代閾值選取不當(dāng)?shù)脑?，重?gòu)后的信號(hào)將包含大量的噪聲成分［6］。

為了解決上述問題，本文作者提出一種改進(jìn)停止閾值準(zhǔn)則的正交匹配追蹤（Improved Orthogonal Matching Pursuit，IOMP）算法，將其應(yīng)用于振動(dòng)信號(hào)的降噪階段，并且根據(jù)不同類型信號(hào)的稀疏系數(shù)分布情況，計(jì)算出不同的迭代終止閾值ε，在經(jīng)典重構(gòu)算法中加入一種篩選判斷條件。此判斷機(jī)制可以避免誤差的持續(xù)迭代，提高了重構(gòu)信號(hào)的精確度。

目前在使用小波降噪理論處理非平穩(wěn)振動(dòng)信號(hào)時(shí)，需要設(shè)定閾值對(duì)分解系數(shù)進(jìn)行處理，而閾值的設(shè)置恰當(dāng)與否，很大程度上影響著降噪效果［7］。文中方法還可避免小波閾值設(shè)定的問題，通過對(duì)振動(dòng)信號(hào)進(jìn)行壓縮投影和重構(gòu)實(shí)現(xiàn)信號(hào)的降噪。因此，使用優(yōu)化的降噪算法對(duì)噪聲信號(hào)進(jìn)行去噪處理，會(huì)提高后續(xù)信號(hào)處理結(jié)果的準(zhǔn)確性。

1 基于壓縮感知理論的信號(hào)降噪數(shù)學(xué)模型

1.1 噪聲信號(hào)與觀測矩陣的構(gòu)造

首先定義一個(gè)不含噪聲的振動(dòng)信號(hào)x∈RN，定義另一個(gè)噪聲信號(hào)為n∈RN，將這兩個(gè)信號(hào)進(jìn)行和運(yùn)算得到含噪信號(hào)c＝（x＋n）∈RN。

選取一個(gè)觀測矩陣，對(duì)原始信號(hào)先進(jìn)行觀測降維。文中采用的觀測矩陣Φ∈RM×N（M ＜N）是隨機(jī)高斯測量矩陣［8］，該矩陣的每個(gè)元素獨(dú)立服從均值為0、方差為1／M的高斯分布。即

含噪信號(hào)c在觀測矩陣Φ下進(jìn)行觀測投影，得到的觀測量y為

1.2 稀疏系數(shù)的分布

目前，用來進(jìn)行信號(hào)稀疏表示的方法很多，對(duì)不同類型的信號(hào)使用不同的稀疏字典其效果是千差萬別的［9］。而文中的研究對(duì)象主要是滾動(dòng)軸承的振動(dòng)信號(hào)，這種類型的信號(hào)在離散余弦變換［10］下具有不錯(cuò)的稀疏性，但是噪聲信號(hào)不屬于稀疏信號(hào)，故DCT（Discrete Cosine Transform）稀疏字典對(duì)噪聲信號(hào)沒有稀疏化的能力。在后續(xù)的觀測階段，振動(dòng)信號(hào)會(huì)保存下來，然而噪聲信號(hào)是不能恢復(fù)重構(gòu)的。

DCT 字典矩陣Ψ∈RN×N的元素可通過式（3）進(jìn)行計(jì)算；

由壓縮感知理論知，若一維信號(hào)x2∈RN在DCT字典矩陣下是稀疏的，則信號(hào)x2在該字典矩陣下的稀疏表示為

式中：Ψ-1∈RN×N為稀疏矩陣；s∈RN×1為系數(shù)向量。

所以振動(dòng)信號(hào)x和噪聲信號(hào)n在DCT 字典矩陣下稀疏表示分別為式（5）（6）所示：

式中：s1和s2分別表示純凈振動(dòng)信號(hào)和噪聲信號(hào)。

在Ψ對(duì)應(yīng)空間上的稀疏系數(shù)向量：

將式（5）和式（6）代入式（2）中可得

令A(yù)＝Φ·Ψ-1，則

式中：A∈RM×N稱為感知矩陣。

得到s1后，繪制該稀疏系數(shù)的分布直方圖，再用正態(tài)分布檢驗(yàn)的方法檢驗(yàn)其是否服從正態(tài)分布；若服從，用式（11）（12）計(jì)算出s1的數(shù)學(xué)期望μ和標(biāo)準(zhǔn)差σ。

式中：n為信號(hào)的采樣數(shù)目；pk是每個(gè)稀疏系數(shù)s1出現(xiàn)的概率。

再以3σ準(zhǔn)則計(jì)算出的迭代終止閾值ε：

當(dāng)更新殘差的二范數(shù)小于迭代終止閾值ε時(shí)，即，則需要停止迭代循環(huán)，否則重構(gòu)信號(hào)會(huì)包含大量的噪聲，降噪效果會(huì)大幅降低。不同類型信號(hào)的稀疏系數(shù)分布情況是不同的，以3σ準(zhǔn)則計(jì)算出的迭代終止閾值ε能將99.74%的稀疏系數(shù)按照大小關(guān)系區(qū)分開［11］。

在通常情況下，判斷感知矩陣A是否滿足約束等距條件是一個(gè)NP-hard 問題，因此文獻(xiàn)［12］指出RIP 的等價(jià)條件為：若觀測矩陣Φ與稀疏表示字典矩陣ψ不相關(guān)，即要求Φ的行向量不能由ψ的列向量表示，且ψ的列向量不能由Φ的行向量表示。如果矩陣Φ和矩陣ψ不相關(guān)，就可以應(yīng)用相關(guān)重構(gòu)算法得到稀疏向量s1的估計(jì)值，進(jìn)而得到降噪后的重構(gòu)信號(hào)。

1.3 OMP 重構(gòu)算法及其改進(jìn)方法

前面兩個(gè)小節(jié)構(gòu)造了觀測矩陣和稀疏字典矩陣，滿足振動(dòng)信號(hào)可以使用壓縮感知理論降噪的前提條件之后，接下來需要用重構(gòu)算法求解如下的最優(yōu)化問題：

為了求解此問題，可以將式（14）轉(zhuǎn)化為式（15），這可以將零范數(shù)轉(zhuǎn)化為l1范數(shù)問題

令A(yù)＝ΦΨ，則式（15）可寫為

與全樣本一致，管理層能力對(duì)研發(fā)投入影響在國有組和民營組均表現(xiàn)為1%顯著性水平的抑制作用。二者相關(guān)系數(shù)分別為-0.0755和-0.0518，這說明，雖然我們不能忽視高技術(shù)企業(yè)中國有經(jīng)濟(jì)占據(jù)較高比重的合理性，但國營企業(yè)繁冗的決策機(jī)制與治理架構(gòu)，加之經(jīng)營目標(biāo)的多元性，可能使高能力管理者進(jìn)行高風(fēng)險(xiǎn)創(chuàng)新活動(dòng)的意愿更低，驗(yàn)證了假設(shè)1b。

對(duì)于上述問題，目前較多使用OMP 算法來求解，得到稀疏向量的估計(jì)值后，代入式（4）便可恢復(fù)出原始純凈的振動(dòng)信號(hào)。

但是在OMP 算法的迭代運(yùn)算中，迭代誤差將會(huì)在下次循環(huán)中繼續(xù)傳導(dǎo)，導(dǎo)致該算法的穩(wěn)定性較差［13］，因此本文作者將一種改進(jìn)停止閾值準(zhǔn)則的OMP 運(yùn)用到壓縮感知的解碼端。為了改進(jìn)該算法的重構(gòu)效果，在求解最小二乘解的步驟之后加入了篩選判斷條件［14］，并加入迭代終止閾值ε。此判斷機(jī)制可以避免誤差的持續(xù)迭代，提高了重構(gòu)信號(hào)的精確度；迭代終止閾值ε能夠確保只將振動(dòng)信號(hào)的主要信息恢復(fù)出來，把噪聲成分在觀測投影的過程中丟棄，并且在重構(gòu)的過程中不提取出來，使得重構(gòu)信號(hào)既保留原始純凈信號(hào)的主要特征，又具有較好的去噪效果，且不降低重構(gòu)精度；因此確定迭代終止閾值ε非常重要。

首先給出改進(jìn)終止閾值準(zhǔn)則的OMP 算法的步驟如下：

由上述步驟可以得到如圖1 所示的改進(jìn)終止閾值準(zhǔn)則的OMP 算法流程。由此基于壓縮感知理論的振動(dòng)信號(hào)降噪模型建立完成，下文將采用仿真信號(hào)與實(shí)測信號(hào)相結(jié)合的方法探究所提改進(jìn)方法的可行性和有效性。

圖1 改進(jìn)終止閾值準(zhǔn)則的OMP 算法流程Fig.1 OMP algorithm flowchart to improve termination threshold guidelines

2 試驗(yàn)驗(yàn)證

2.1 仿真信號(hào)驗(yàn)證

為了模擬滾動(dòng)軸承振動(dòng)信號(hào)的真實(shí)工況，可以定義仿真信號(hào)［15］進(jìn)行數(shù)值模擬的仿真分析，如式（20）所示：

式中：x（t）＝2sin（2π×ft1nT）×5.5sin（2π×ft2nT）＋4sin（2π×ft3nT）；n（t）為SNR 為2 dB 的高斯白噪聲；采樣間隔T＝1／4 000 s；采樣數(shù)n＝1 000；組成x（t）的頻率分量分別為ft1＝50 Hz，ft2＝1 600 Hz，ft3＝20 Hz。x（t）與n（t）的信號(hào)波形如圖2、3所示。

圖3 高斯白噪聲信號(hào)的時(shí)域波形Fig.3 Time domain waveform of a white Gaussian noise signal

然后將這兩種信號(hào)進(jìn)行和運(yùn)算得到帶噪信號(hào)y（t）的波形，如圖4 所示。3 種信號(hào)的定義完成后，則需對(duì)其稀疏性進(jìn)行分析。首先利用DCT 字典矩陣將信號(hào)進(jìn)行稀疏表示，得到純凈信號(hào)與高斯噪聲信號(hào)稀疏系數(shù)分別如圖5、6 所示，可以看出：這兩種信號(hào)經(jīng)過字典稀疏表示后，原始信號(hào)x（t）大量的稀疏系數(shù)都是較小值，接近于零值上下浮動(dòng)，只有極個(gè)別的數(shù)值較大；反觀噪聲信號(hào)的稀疏系數(shù)并沒有在零值上下浮動(dòng)，數(shù)值較大，接近零的值很少。這表明DCT字典矩陣對(duì)x（t）的稀疏表示效果非常好，對(duì)n（t）并沒有很好的稀疏性。

圖4 帶噪信號(hào)時(shí)域波形Fig.4 Time domain waveform with noisy signal

圖5 純凈仿真信號(hào)的稀疏表示系數(shù)向量Fig.5 Sparse represented coefficient vector for pure simulated signal

圖6 高斯噪聲信號(hào)的稀疏表示系數(shù)向量Fig.6 Sparse represented coefficient vector for Gaussian noise signal

為了詳細(xì)掌握稀疏系數(shù)s1的分布范圍，則需要畫出該組系數(shù)的頻數(shù)分布直方圖如圖7 所示。可清晰地發(fā)現(xiàn)，大量數(shù)值分布在零值附近，并且分布情況也近似滿足正態(tài)分布曲線。

圖7 純凈仿真信號(hào)的稀疏表示系數(shù)分布Fig.7 Sparse represented coefficient distribution of pure simulated signal

圖8 純凈仿真信號(hào)的正態(tài)分布檢驗(yàn)圖Fig.8 Normal distribution test plot of pure simulated signal

然后分別應(yīng)用正交匹配追蹤（Orthogonal Matching Pursuit，OMP）算法和改進(jìn)停止閾值準(zhǔn)則的正交匹配追蹤算法對(duì)帶噪信號(hào)進(jìn)行降噪，得到兩種方法的降噪信號(hào)波形如圖9、10 所示?？芍焊倪M(jìn)停止閾值準(zhǔn)則的正交匹配追蹤算法所重構(gòu)的信號(hào)明顯比正交匹配追蹤算法所重構(gòu)信號(hào)的時(shí)域波形更接近于原始純凈的信號(hào)。

圖9 OMP 算法降噪信號(hào)時(shí)域波形Fig.9 OMP algorithm noise reduction signal time domain waveform

圖10 改進(jìn)停止閾值準(zhǔn)則OMP 算法降噪信號(hào)波形Fig.10 Improved stop threshold criterion OMP algorithm noise reduction signal waveform

為了進(jìn)一步評(píng)判降噪效果，可以定義信噪比（Signal to Noise Ratio，SNR）來衡量降噪性能的高低：

計(jì)算結(jié)果如表1 所示，可知：降噪前信號(hào)的信噪比為-1.025 7 dB，使用OMP 重構(gòu)算法降噪后信號(hào)的信噪比為1.292 5 dB，使用改進(jìn)停止閾值準(zhǔn)則的OMP重構(gòu)算法降噪后信號(hào)的信噪比為2.054 9 dB。通過對(duì)比可以看出，使用文中降噪方法可有效提升信號(hào)的信噪比，并且重構(gòu)信號(hào)所需的觀測量遠(yuǎn)小于其他降噪方法，從而驗(yàn)證了所提降噪方法的有效性。

表1 不同重構(gòu)算法降噪信號(hào)的信噪比Tab.1 Signal to noise ratio of denoised signals by different reconstruction algorithms

2.2 實(shí)測振動(dòng)信號(hào)驗(yàn)證

旋轉(zhuǎn)機(jī)械設(shè)備的運(yùn)行狀態(tài)復(fù)雜多變，這些設(shè)備的振動(dòng)信號(hào)大多屬于非平穩(wěn)信號(hào)。這類信號(hào)的統(tǒng)計(jì)特性都隨著時(shí)間變化而變化，因此降噪難度比較高。

下面使用實(shí)測振動(dòng)信號(hào)來檢驗(yàn)方法的效果，試驗(yàn)中采用的實(shí)測信號(hào)來源于美國凱斯西儲(chǔ)大學(xué)軸承實(shí)驗(yàn)中心［16］，軸承型號(hào)為6205-2RS JEK SKF 深溝球軸承，采樣頻率為12 kHz，轉(zhuǎn)速為1 750 r／min，負(fù)載746 W（2 HP），該振動(dòng)信號(hào)的時(shí)域波形如圖11 所示。在該信號(hào)中加入高斯白噪聲形成強(qiáng)噪聲背景下的帶噪信號(hào)，如圖12 所示。

圖11 實(shí)測振動(dòng)信號(hào)波形Fig.11 Measured vibration signal waveform

圖12 帶噪信號(hào)時(shí)域波形Fig.12 Time domain waveform with noise signal

實(shí)測振動(dòng)信號(hào)在DCT 字典下的稀疏表示如圖13所示，可以發(fā)現(xiàn)：0-400 之間的稀疏表示系數(shù)值不為零的較多，400-1 000 之間的系數(shù)接近于零值上下浮動(dòng)。這是因?yàn)閷?shí)測信號(hào)屬于隨機(jī)信號(hào)，這類信號(hào)的統(tǒng)計(jì)特性都隨著時(shí)間變化而變化，因此稀疏性不比仿真定義的信號(hào)好。這就為重構(gòu)算法提出了苛刻的條件。

圖13 實(shí)測振動(dòng)信號(hào)的稀疏表示系數(shù)向量Fig.13 Sparse represented coefficient vector for measured vibration signal

為了進(jìn)一步檢驗(yàn)該組數(shù)據(jù)是否滿足正態(tài)分布，需要對(duì)其進(jìn)行正態(tài)分布檢驗(yàn)，得到圖14，其中橫軸為變量區(qū)間，縱軸為累計(jì)概率分布，紅色點(diǎn)劃線是以稀疏系數(shù)s1的數(shù)學(xué)期望μ＝0.074 0 和標(biāo)準(zhǔn)差σ＝0.001 2構(gòu)建的標(biāo)準(zhǔn)正態(tài)分布線，如果藍(lán)色數(shù)據(jù)樣本點(diǎn)越接近紅色點(diǎn)劃線，則說明樣本點(diǎn)的正態(tài)特性越好。由圖14 可知：70%～80%數(shù)據(jù)樣本點(diǎn)都分布在紅色點(diǎn)劃線上，20%～30%的數(shù)據(jù)樣本點(diǎn)遠(yuǎn)離紅色點(diǎn)劃線。因此該組數(shù)據(jù)是服從μ＝0.074 0、σ＝0.001 2 的正態(tài)分布，從而得到迭代終止閾值ε＝0.434。

圖14 實(shí)測振動(dòng)信號(hào)的正態(tài)分布檢驗(yàn)圖Fig.14 Normal distribution test plot for measured vibration signal

然后分別應(yīng)用正交匹配追蹤算法和改進(jìn)停止閾值準(zhǔn)則的正交匹配追蹤算法對(duì)帶噪信號(hào)進(jìn)行降噪，得到兩種方法的降噪信號(hào)波形如圖15、16 所示?？芍焊倪M(jìn)停止閾值準(zhǔn)則的正交匹配追蹤算法所重構(gòu)的信號(hào)明顯比正交匹配追蹤算法所重構(gòu)的信號(hào)在縱坐標(biāo)幅值的大小和波形方面更接近原始的實(shí)測振動(dòng)信號(hào)。

圖15 OMP 算法降噪信號(hào)時(shí)域波形Fig.15 OMP algorithm noise reduction signal time domain waveform

圖16 改進(jìn)停止閾值準(zhǔn)則的OMP 算法降噪信號(hào)Fig.16 Improved OMP algorithm for noise reduction based on stop threshold criterion

定量地評(píng)判降噪效果，則可以定義信噪比（SNR）和均方誤差（Mean Square Error，MSE）兩個(gè)指標(biāo)來衡量降噪性能。信噪比的定義為式（21），均方誤差見式（22）：

經(jīng)過計(jì)算，降噪前信號(hào)的信噪比為-1.268 dB，降噪后信號(hào)的信噪比是2.541 8 dB、均方誤差MSE為0.030。由此可知，文中方法可以有效提升實(shí)測振動(dòng)信號(hào)的信噪比，減小均方誤差，從而得證文中所提的降噪方法同樣適用于實(shí)測滾動(dòng)軸承振動(dòng)信號(hào)的降噪。

為了比較經(jīng)典重構(gòu)算法與文中所提方法降噪效果的區(qū)別，下面使用子空間追蹤法（Subspace Pursuit，SP）、正交匹配追蹤法（Orthogonal Matching Pursuit，OMP）、分段弱正交匹配追蹤法（Stage wise Weak Orthogonal Matching Pursuit，SWOMP）3 種重構(gòu)算法在相同條件下進(jìn)行信號(hào)降噪。

結(jié)果如表2 所示。降噪后的信號(hào)信噪比最大，則認(rèn)為該方法的降噪效果越好；均方誤差越小，則認(rèn)為降噪后的信號(hào)越接近原始純凈信號(hào)。SP、OMP、SWOMP 三種方法降噪信號(hào)的信噪比分別為0.330 8、1.664 2、1.990 7 dB，均方誤差MSE 分別為0.074 2、0.055 7、0.045。顯而易見，改進(jìn)停止閾值準(zhǔn)則的正交匹配追蹤法優(yōu)于SP、OMP、SWOMP 三種降噪方法，可以獲得較好的降噪效果。

表2 不同重構(gòu)算法降噪信號(hào)的SNR 和MSETab.2 SNR and MSE of noise reduction signal with different reconstruction algorithms

3 結(jié)論

本文作者針對(duì)OMP 算法迭代誤差會(huì)在迭代循環(huán)中繼續(xù)傳導(dǎo)、終止閾值選取不當(dāng)導(dǎo)致信號(hào)成分和噪聲成分發(fā)生混疊等問題，提出一種改進(jìn)停止閾值準(zhǔn)則的正交匹配追蹤算法，建立基于壓縮感知理論的振動(dòng)信號(hào)降噪模型。

通過對(duì)仿真信號(hào)的降噪分析，提升了帶噪信號(hào)的信噪比，驗(yàn)證了該方法的有效性；然后對(duì)深溝球軸承的實(shí)測振動(dòng)信號(hào)進(jìn)行降噪處理，在強(qiáng)噪聲背景下信號(hào)的信噪比提升到2.541 8 dB、均方誤差減小到0.030。并且和其他重構(gòu)算法降噪后重構(gòu)信號(hào)的結(jié)果作比較，表明文中所提方法的有效性與優(yōu)越性，并且信號(hào)經(jīng)過稀疏表示后，所需要處理的數(shù)據(jù)量遠(yuǎn)小于其他降噪方法。