黃文健,黃瑾珉,曹承昊,陳漢新
武漢工程大學機電工程學院,湖北 武漢 430205
齒輪箱是很多機械設備中重要的部件,其常被用來傳遞動力、改變轉速、分配動力等。在實際應用中齒輪箱的工作狀態與設備的運行狀態息息相關,而其中的齒輪,軸,軸承等零部件長期在高負荷狀態下運行,故障率較高,很容易引發機器故障。在齒輪箱部件失效中由齒輪故障引發的約占60%[1],因此對齒輪箱中齒輪故障診斷的研究就顯得很有必要了。
齒輪常見的故障有齒輪磨損、齒面膠合和擦傷、齒面接觸疲勞和斷齒等。在本次實驗中,選取齒輪裂紋來模擬齒輪故障,得到的實驗信號摻雜著許多噪聲,所以采用小波包變換對信號進行降噪[2],之后再提取特征。主元分析法(Principal Component Analysis,PCA)[3]是將多個線性相關變量壓縮為少數幾個不相關的變量的一種多元統計方法。葉永偉等[4]用PCA來對設備特征進行提取;杜振寧等[5]將小波包分解與PCA相結合來對軸承進行故障診斷。由于PCA具有概念簡單,計算方便等優點,因此常被用來對數據進行降維。
Wald[6]在1947年的著作中系統地介紹了序貫概 率 比 檢 驗(Sequential Probability Ratio Test,SPRT)。由于其具有簡單、高效等特性,因此近些年來在故障診斷方面使用非常普遍。Goodman等[7]將SPRT與自適應波形結合起來識別雷達信號;Ray等[8]將SPRT應用在核電站和飛行器的故障檢測與識別;CHEN等[9-10]將SPRT用于齒輪裂紋的故障診斷。SPRT與固定樣本檢驗不同之處在于這種算法需要更少的平均樣本量,但是卻具有更高的檢驗效率。
在齒輪的故障診斷中,為了避免直接從眾多特征參數中選取SPRT的參數,本文將PCA與SPRT相結合起來應用于齒輪裂紋的故障診斷。實驗中采集的信號在降噪之后,利用PCA來處理信號的特征數據并選擇合適的主元作為檢驗參數以提高檢驗的準確度與可靠性,通過分析最后得到的結果圖來進行識別與分類[11-12],結果表明該方法是有效的。
特征提取是從信號里提取有效信息的過程,適當的方法和有效的特征對后續過程很有幫助。在很多研究中對數據進行降維會使之后的計算分析更加的高效與簡便。在用PCA研究多變量問題時,實際上是通過線性變化將數據從高維空間映射到低維空間,低維空間中較少的新變量能集中地反映原始變量中包含的絕大部分信息。其基本方法是將數據矩陣的協方差矩陣的特征值按大小進行排序,以此順序計算各個主元及其貢獻率,之后選取滿足要求的主元從而達到降維的目的。目前PCA已廣泛用于特征提取、數據降噪以及故障重構等方面。
PCA的基本原理如下:設原始數據樣本包含m個n維特征向量 Xk=( x1,x2,… xn)T,則PCA具體步驟如下:
計算樣本的平均值:

求出樣本向量的協方差矩陣:

計算S的特征值 λ1,λ2,…,λn及其對應的特征向量 v1,v2,…,vn。
計算各主元的貢獻率和累計貢獻率:

選取前p個累計貢獻率滿足要求的主元作為映射矩陣,將原始信號變換到特征空間:

式(5)中,V=(v1,v2,…vp)。
SPRT的基本原理為:設 x1,x2,…,xn是一組獨立同分布的隨機變量序列。對于這樣的總體樣本值,提出兩個假設,零假設(H0∶θ=θ0)和備選假設(H1∶θ=θ1)。它們的聯合概率密度可以定義為:

其中j=0,1,檢驗的似然比計算如下:

式(1)中,f(x/θ)為條件概率分布;θ為分布參數。
在以前的檢驗中,當λn的值在閾值C附近時,對兩種假設的接受與否的判斷顯得不夠合理。而Wald改進的地方是當 λn遠小于C時接受 H0,當λn遠大于C時拒絕 H0,當 λn與C很接近時暫且不做出判斷,而是再繼續檢測λn+1,直到之后的似然比遠大于或遠小于C為止。
在SPRT中,事先給定兩類錯誤的概率α與 β的值,則閾值A、B可以按照下列公式來確定:

檢驗判斷準則為:代入觀測值x1,計算λ1(x),若 λ1(x)<B,則停止觀測并接受零假設 H0;若λ1(x)>A,則停止觀測并接受備選假設 H1;若A≤λ1(x)≤B,則繼續進行檢驗,直至似然比值滿足上述兩種情況為止。
選取正常齒輪F1和有裂紋的故障齒輪F2、F3進行實驗。齒輪裂紋全深度a=2.4 mm,裂紋全寬度b=25 mm,厚度為0.4 mm,齒輪裂紋的角度為45°。實驗中齒輪箱為空載狀態,轉速為800 r/min。3組齒輪的參數如表1所示。

表1 三種故障模式Tab.1 Three fault modes
圖1是齒輪箱的結構圖。在運行中齒輪3和4之間會產生振動,因此可以任選其中之一來模擬故障。本次實驗選擇用齒輪3來進行模擬實驗。

圖1 齒輪箱結構圖Fig.1 Structure of gearbox
齒輪箱信號采集系統如圖2所示,將加速度傳感器安裝在齒輪箱的水平垂直兩個方向上,通過動態模擬器采集振動信號,并通過信號分析器將信號輸入并存儲在計算機中,本文僅研究分析水平方向的信號。分別用S1、S2和S3表示正常情況與兩種故障情況下的信號。

圖2 齒輪箱信號采集系統Fig.2 Signal acquisition system of gearbox
小波包變換將同時分解信號的低頻和高頻部分,這有利于有效地區分信號中的噪聲。對實驗信號進行3層小波包分解,得到8個頻帶的分解信號后再重構信號。這樣重構后的信號會更準確。
振動信號的特征信息是通過一些特征參數來表達的,每個特征參數都能反映該信號的某種特性,平穩的振動信號表明齒輪箱運行正常,而異常的或是有較大波動的振動信號則表明齒輪箱出現了故障。在常用的時域特征參數中,均值能反映信號波動能量的大小;有效值是描述動態信號強度的指標;標準差是描述信號的波動分量;峭度指標能反應振動信號中的沖擊特征;波形指標是反映信號偏離高斯分布程度的指標;脈沖指標與峰值指標都是用于檢測信號中有無沖擊的指標;裕度指標常被用來檢測機械設備的磨損狀況[13-15]。本文選取這8個參數指標進行計算。假設xi=[x1,x2,…,xN]是一組離散待檢信號,N=8 192。每組取1 024個檢驗點,這樣可以得到7196組檢驗數據。
對實驗信號進行降噪之后再提取特征值,本文中n=1 024,特征參數的計算公式如下:


將經過預處理的振動信號代入公式中分別計算,得到齒輪箱在3種狀況下的3組信號的8種特征參數。
對上述3組數據進行PCA處理,根據得到的結果來選取待檢驗的序列。
經過PCA處理之后,根據貢獻率的高低來對處理后的主元進行降序排列,排序越靠前的主元保留的原有信息越多。在得到的處理結果中,通常是選取累計貢獻率高于85%的前幾個主元,這樣選取的主元才能包含足夠多的原有信息。
這些特征值經過PCA處理之后,前4個主元就能得到超過85%的信息。這里選取貢獻率最大的序列作為檢驗序列,記作 yi=[y1…yk]。其均值和標準差為:

經過以上處理之后的待檢序列大致符合正態分布。正常狀態下,該序列滿足零假設H0∶θ=θ0;故障狀態下,該序列滿足備擇假設H1∶θ=θ1。標準差保持為σ。當兩個假設都成立時,該序列的聯合概密度為:

式(20)和(21)中:P0i為零假設條件下的概率密度函數;P1i為備擇假設條件下的概率密度函數。
SPRT似然比為:

式(22)中:P0為零假設條件下的先驗概率;P1為備擇假設條件下的先驗概率。
對似然比公式進一步推導簡化可使得計算更加簡便:

簡化之后的閾值a=lnA,b=lnB。假設檢驗范兩類錯誤的概率分別取α=0.001,β=0.001,則計算可得a=6.91,b=-6.91。計算出似然比值后根據檢驗準則進行判斷。如果Δ<b,則停止檢驗并接受H0,即為正常狀態;如果Δ>a,停止檢驗并接受H1,即為故障狀態;如果b<Δ<a,則繼續比較下一個似然比值與a、b的關系,如此循環,直至得出結果。運用此方法對齒輪裂紋進行診斷的完整流程如圖3所示。

圖3 齒輪裂紋診斷流程圖Fig.3 Flowchart of gear crack fault diagnosis
利用均方根誤差(Root Mean Square Error,RMSE)來檢驗SPRT應用于此的可靠性。均方根誤差計算公式如下:

對3組不同的齒輪分別進行實驗,得到圖4(a)中的原始振動信號。其中S1、S2和S3分別為正常齒輪、25%和50%裂紋齒輪的振動信號。3組原始振動信號經過降噪之后,得到圖4(b)中的信號。
對比圖 4(a)和圖 4(b)可以看到,與原始振動信號相比,降噪后信號的幅值明顯有所減小,而且信號曲線也更為清晰。圖4(b)基本達到了預期效果,這也表明了小波包變換能很好地對原始振動信號進行降噪。
對3組實驗信號進行特征提取之后再用主元分析法處理,發現前4個主元包含超過85%的有效信息,信號S1、S2和S3的各主元貢獻率排序如表2所示。由PCA的處理結果可知這3組信號的第一主元的貢獻率分別達到了50.9%、50.2%、51.8%。由于排名第一的主元所占的比例最大,因此被選用來進行SPRT。

圖4 實驗信號:(a)原始振動信號,(b)去噪后的信號Fig.4 Experimental signals:(a)original vibration signals,(b)de-noised signals

表2 主元貢獻率及累計貢獻率Tab.2 Principal component contribution rate and accumulative contribution rate %
由似然比計算公式可知,平均值的變化對似然比影響較大。將正常信號S1的平均值記為μ0,故障信號S2和S3的平均值記為μ1。對信號進行檢測,結果如圖5所示。
當輸入的信號是 S1和S2時,檢驗結果如圖5(a)所示,從圖中可以判定,信號S1的似然比Δ1,2(S1)<b,即齒輪箱正常工作;信號 S2的似然比Δ1,2(S2)>a,即齒輪箱出現故障。當輸入的信號是S1和S3時,檢驗結果如圖5(a)所示,從圖中可以判定,信號S1的似然比Δ1,3(S1)<b,即齒輪箱正常工作;信號S3的似然比Δ1,3(S3)>a,即齒輪箱出現故障。

圖5 SPRT結果:(a)似然比(Δ1,2(S1),Δ1,2(S2),Δ1,3(S1),Δ1,3(S3))與序貫檢驗迭代次數的關系,(b)似然比(Δ2,3(S2),Δ2,3(S3))與序貫檢驗迭代次數的關系Fig.5 Results of SPRT:(a)relationship between SPRT likelihood ratios(Δ1,2(S1),Δ1,2(S2),Δ1,3(S1),Δ1,3(S3))and the number of sequential testing iteration,(b)relationship between SPRT likelihood ratios(Δ2,3(S2),Δ2,3(S3))and the number of sequential testing iteration
在實驗中發現不同程度裂紋的齒輪也可以通過這種方法加以區分,其方法基本相同。就本次實驗而言,將兩種故障信號S2和S3進行檢測,檢驗結果如圖5(b)所示,從圖中可以判定,信號S2的似然比Δ2,3(S2)<b,即齒輪具有輕度損傷;信號S3的似然比Δ2,3(S3)>a,即齒輪出現重度故障。所以根據檢驗結果可以對25%裂紋和50%裂紋的信號進行清晰地區分。
在實驗中,正常狀況下獲取兩組數據,分別進行降噪與PCA處理后共得到14 338個數據點。將這些數據分為16組,分別記為A1至A16。對16組數據分別進行序貫概率比檢驗,得到16個似然比值,記為B1至B16。計算這16個值之間的均方根誤差,即B1分別其他15個值比較,B2分別與其他14個值比較,…,B15與B16比較,得到15組比值,選取每組比值中的最大值形成誤差曲線。將故障狀況下的似然比記為a,計算故障狀況與正常狀況之間的誤差,即a與B1,B2…B16比較,將得到16個誤差值形成誤差曲線。通過比較誤差來驗證序貫概率比檢驗診斷齒輪裂紋故障的可靠性。
將正常信號S1的平均值記為 μ0,故障信號S2的平均值記為 μ1。將16組正常數據信號的似然比Δ1,2(S1)輸入程序中,得到了 16 條 normal曲線;將2組故障信號的似然比Δ1,2(S2),Δ1,2(S3)輸入程序中,得到了2條fault曲線,最后的結果如圖6(a)所示。圖6(a)表明16組正常曲線的趨勢基本相同,而故障信號的曲線與之相比則有較大的區別。計算16組正常信號的似然比之間的均方根誤差E1,2(S1),每組誤差中的最大值形成的誤差曲線如圖6(b)所示。分別計算故障信號S1和S2與16組正常信號似然比之間的誤差E1,2(S2)和E1,2(S3),得到圖6(b)中的另外兩條誤差曲線。在誤差方面,從圖6(b)中可以看出不同故障之間的誤差遠大于同種故障之間的,這表明序貫概率比檢驗可以很好地區分不同故障狀況。

圖6 SPRT及誤差結果:(a)似然比(Δ1,(2S1),Δ1,(2S2),Δ1,(2S3))與序貫檢驗迭代次數的關系,(b)均方根誤差值E1,(2S1),E1,(2S2)和E1,(2S3)Fig.6 Results of SPRT and RMSE:(a)relationship between SPRT likelihood ratios(Δ1,(2S1),Δ1,(2S2),Δ1,(2S3))and thenumber of sequential testing iteration,(b)RMSEs E1,(2S1),E1,(2S2)and E1,(2S3)
將正常信號S1的平均值記為 μ0,故障信號S3的平均值記為 μ1。將16組正常數據信號的似然比Δ1,3(S1)輸入程序中,得到了 16 條 normal曲線;將兩組故障信號的似然比Δ1,3(S2),Δ1,3(S3)輸入程序中,得到了2條fault曲線,最后的結果如圖7(a)所示。圖7(a)說明16組正常曲線的趨勢基本相同,而故障信號的曲線與之相比則有較大的區別。計算16組正常信號的似然比之間的均方根誤差E1,3(S1),每組誤差中的最大值形成的誤差曲線如圖7(b)所示。分別計算故障信號S1和S2與16組正常信號似然比之間的誤差E1,2(S2)和E1,3(S3),得到圖7(b)中的另外2條誤差曲線。在誤差方面,從圖7(b)中可以看出不同故障之間的誤差遠大于同種故障之間的,這表明SPRT可以很好地區分不同故障狀況。

圖7 SPRT及均方根誤差結果:(a)似然比(Δ1,(3S1),Δ1,(3S2),Δ1,(3S3))與序貫檢驗迭代次數的關系,(b)均方根誤差值E1,(3S1),E1,(3S2)和E1,(3S3)Fig.7 Results of SPRT and RMSE:(a)relationship between SPRT likelihood ratios(Δ1,(3S1),Δ1,(3S2),Δ1,(3S3))and thenumber of sequential testing iteration;(b)RMSEs E1,(3S1),E1,(3S2)and E1,(3S3)
由圖6和圖7可知,當以正常信號的平均值為參數μ0時,不論是以哪種故障信號的平均值為參數 μ1來建立模型,從最后的結果中都能對3種狀況進行清晰地區分。
本文提出了一種將PCA與SPRT相結合的檢測方法,并選取不同程度裂紋的故障齒輪與正常齒輪進行模擬試驗,來檢測此方法能否有效區分不同狀態下的齒輪。在預處理之后利用PCA對實驗信號所選取的特征參數進行降維處理,有效地減少了特征數量,使得之后的檢驗更加高效準確,最終有效地識別出了齒輪箱不同的狀態。此外,經過多次實驗后證明在其他的工作條件下,該方法也能對齒輪箱不同的狀態進行清晰地區分。因此將該方法應用于齒輪裂紋故障檢測與分類方面是有效而且可靠的。此方法也可供其他旋轉機械的診斷研究進行參考。