王國鋒,李玉波,秦旭達(dá),喻 秀,李啟銘
(天津大學(xué)機(jī)械工程學(xué)院,天津 300072)
機(jī)械故障診斷在現(xiàn)代工業(yè)生產(chǎn)中起著越來越重要的作用.對(duì)反映機(jī)械設(shè)備狀況的振動(dòng)信號(hào)進(jìn)行信號(hào)分析及特征提取,然后利用數(shù)學(xué)模型進(jìn)行故障識(shí)別是目前故障診斷的主要方法.
滾動(dòng)軸承是機(jī)械設(shè)備中故障常發(fā)部件,多變的工況條件以及設(shè)備自身的非線性,使反映軸承故障的動(dòng)態(tài)信號(hào)表現(xiàn)出較強(qiáng)的非平穩(wěn)特性.目前信號(hào)主要的時(shí)頻分析方法有短時(shí)傅里葉變換、Wigner-Ville 分布、小波尺度分析等.這些方法的時(shí)頻分辨率都要受到不確定性原理的限制,因此無法同時(shí)在時(shí)域和頻域均得到很好的分辨率[1].
故障識(shí)別模型以神經(jīng)網(wǎng)絡(luò)模型應(yīng)用最廣,但利用神經(jīng)網(wǎng)絡(luò)模型進(jìn)行故障識(shí)別需要大量的訓(xùn)練樣本,在有限的訓(xùn)練樣本情況下,泛化能力差,且優(yōu)化過程存在陷入局部極值的可能[2];另外,大部分神經(jīng)網(wǎng)絡(luò)模型處理的是靜態(tài)模式分類問題,而旋轉(zhuǎn)機(jī)械的故障行為是一個(gè)動(dòng)態(tài)變化的過程.
為此,筆者提出了應(yīng)用時(shí)變自回歸(time-varying autoregressive ,TVAR) 和 隱 馬 爾 科 夫 模 型 (hidden Markov model,HMM)的滾動(dòng)軸承故障診斷方法.時(shí)變自回歸模型是一種基于參數(shù)模型的時(shí)頻分析方法.它是在 AR(autoregressive)建模的基礎(chǔ)上發(fā)展起來的,TVAR 方法得出的時(shí)頻譜圖具有分辨率高、無交叉干擾項(xiàng)以及計(jì)算速度快等優(yōu)點(diǎn),非常適合處理非平穩(wěn)信號(hào),已在故障診斷領(lǐng)域有所應(yīng)用[3-5].隱馬爾科夫模型是信號(hào)動(dòng)態(tài)時(shí)間序列統(tǒng)計(jì)模型,能夠?qū)σ粋€(gè)時(shí)間跨度上的信息進(jìn)行統(tǒng)計(jì)建模和分類,具有嚴(yán)謹(jǐn)?shù)臄?shù)據(jù)結(jié)構(gòu)和可靠的計(jì)算性能.它能通過較少的樣本訓(xùn)練出可靠的模型,并按模式匹配原理,尋找與未知信號(hào)最相似的模式作為識(shí)別結(jié)果,現(xiàn)已成為語音識(shí)別的主流技術(shù)[6].傳統(tǒng)故障診斷方法一直停留在靜態(tài)分析基礎(chǔ)上,忽略了故障變化發(fā)展的動(dòng)態(tài)信息,而 HMM 用于故障診斷,可對(duì)故障變化的動(dòng)態(tài)信息進(jìn)行統(tǒng)計(jì)分析,能夠盡早地發(fā)現(xiàn)故障發(fā)展的跡象.與神經(jīng)網(wǎng)絡(luò)相比,保留了更多訓(xùn)練數(shù)據(jù)的統(tǒng)計(jì)信息,具有更高的識(shí)別率和魯棒性[7].
筆者結(jié)合了兩種模型的優(yōu)點(diǎn),研究了非平穩(wěn)信號(hào)的時(shí)變參數(shù)建模以及隱馬爾科夫建模故障識(shí)別過程,并通過實(shí)驗(yàn)驗(yàn)證該方法在滾動(dòng)軸承故障診斷中的可靠性.
{x(t),t=1,2,…,N }是一非平穩(wěn)時(shí)間序列,TVAR是利用先前值 x (t ? i)來估計(jì)現(xiàn)在值 x(t),可寫成
式中:x (t)為信號(hào)的采集值;ai(t)(i=1,2,…,p)為時(shí)變參數(shù);p 為模型階數(shù);e(t)為模型殘差,服從零均值正態(tài)分布.
將時(shí)變參數(shù){ai(t),i=1,2,…,p}作為一組基時(shí)間函數(shù) gj(t )的線性組合,即
式中:m 為擴(kuò)展維數(shù);aij為組合的時(shí)不變權(quán)值.
將式(2)代入式(1),得
式中
這樣模型時(shí)變參數(shù)的估計(jì)問題就變成常參數(shù) aij的估計(jì),從而把一個(gè)線性非平穩(wěn)問題轉(zhuǎn)化為一個(gè)線性時(shí)不變問題.
文獻(xiàn)[8]闡述了利用最小二乘法估計(jì)模型參數(shù) A的過程,在估計(jì)出參數(shù) A 之后就可以根據(jù)式(2)計(jì)算出時(shí)變參數(shù) ai(t)(i =1,2,…,p),繼而求解模型方差,即
瞬時(shí)頻率能夠通過 TVAR 參數(shù) ai(t )提取.用各個(gè)時(shí)刻的參數(shù)依次計(jì)算出每個(gè)時(shí)刻的參數(shù)譜圖,疊加在一起就可以構(gòu)造一個(gè)三維時(shí)頻譜圖,因此TVAR 譜函數(shù)定義為
HMM 是一個(gè)雙內(nèi)嵌式隨機(jī)過程,即 HMM 由兩個(gè)隨機(jī)過程組成:一是隱狀態(tài)轉(zhuǎn)移序列,它對(duì)應(yīng)一個(gè)單純的 Markov 過程;二是與隱狀態(tài)有關(guān)的觀測(cè)序列.其中隱狀態(tài)轉(zhuǎn)移序列是不可觀測(cè)的,只能通過另一個(gè)隨機(jī)過程的輸出觀測(cè)序列進(jìn)行推斷,所以稱之為隱馬爾科夫模型[9].
HMM 可以記為 λ={Nλ,Mλ,πλ,Aλ,Bλ}或簡(jiǎn)寫為λ={πλ,Aλ,Bλ}.其中 Nλ為模型中 Markov 鏈的狀態(tài)數(shù)目;Mλ為每個(gè)狀態(tài)對(duì)應(yīng)的可能的觀測(cè)值數(shù)目;πλ為初始概率分布矢量;Aλ為狀態(tài)轉(zhuǎn)移概率矩陣;Bλ為觀測(cè)值概率矩陣,對(duì)于連續(xù)的 HMM,Bλ是一組觀測(cè)值概率函數(shù).
圖1為1 個(gè)三狀態(tài)的離散HMM 的例子,由圖可以看出,1 個(gè)HMM 可以由1 組狀態(tài)轉(zhuǎn)移概率以及與每一個(gè)狀態(tài)有關(guān)的觀測(cè)分布概率來定義.
圖1 三狀態(tài)HMM示意Fig.1 Diagram of three-state HMM
一個(gè)標(biāo)準(zhǔn)的HMM 模型需要解決概率計(jì)算、模型訓(xùn)練和隱狀態(tài)估計(jì)3 個(gè)基本問題.
(1)概率計(jì)算.給定模型參數(shù) λ={πλ,Aλ,Bλ}、觀測(cè)序列 O={o1,o2,…,ot},計(jì)算觀測(cè)序列 O 在給定模型下的概率,即計(jì)算似然 P(O︱λ)或者對(duì)數(shù)似然ln,P(O︱λ),它表示參數(shù) λ 擬合數(shù)據(jù) O 的準(zhǔn)確程度.可通過前向-后向算法對(duì)HMM 進(jìn)行遍歷實(shí)現(xiàn).
(2)模型訓(xùn)練,即模型參數(shù)估計(jì)問題.也就是說,對(duì)于初始模型和給定用于訓(xùn)練的觀測(cè)序列 O={o1,o2,…,ot},如何調(diào)整模型 λ={πλ,Aλ,Bλ}的參數(shù),使其能夠最好地?cái)M合觀測(cè)數(shù)據(jù) ? argmaxλ=f (O?)λ,這可以通過EM 算法或Baum-Welch 算法實(shí)現(xiàn).
(3)隱狀態(tài)估計(jì).給定模型參數(shù) λ={πλ,Aλ,Bλ}和觀測(cè)序列 O={o1,o2,…,ot},基于某種最優(yōu)準(zhǔn)則估計(jì)最可能的隱狀態(tài)序列 Q={q1,q2,…,qt},即估計(jì)產(chǎn)生觀測(cè)序列的最可能經(jīng)過路徑,這可以通過 Viterbi 算法估計(jì)實(shí)現(xiàn).
圖2所示為利用 TVAR-HMM 進(jìn)行故障診斷的主要步驟.
圖2 TVAR-HMM故障診斷流程Fig.2 Flow chart of fault diagnosis based on TVAR-HMM
(1)對(duì)不同故障類型的振動(dòng)信號(hào)進(jìn)行 TVAR 建模,得到時(shí)頻譜圖.建模過程中,需要選定基函數(shù)gj(t)、模型階數(shù) p 和擴(kuò)展維數(shù)m ,文獻(xiàn)[10]詳細(xì)介紹了它們的選擇方法.
(2)基于能量法提取故障特征.對(duì)信號(hào) TVAR 建模后,故障信息完全體現(xiàn)在時(shí)頻譜圖上.用能量法提取故障特征可以更加有效地提取時(shí)頻圖上的故障信息.文獻(xiàn)[11]已研究了這種方法的有效性.其主要過程是:將時(shí)頻空間劃分為若干區(qū)域,分別計(jì)算其時(shí)頻局部能量,將時(shí)頻局部能量作為元素,組成特征向量.
(3)用特征向量訓(xùn)練 HMM 模型.特征向量為HMM 模型的觀測(cè)序列,故障類型為隱狀態(tài),1 種故障類型訓(xùn)練1 個(gè)HMM 模型.
(4)利用訓(xùn)練好的模型對(duì)測(cè)試樣本進(jìn)行識(shí)別,即對(duì)于未知的觀測(cè)序列,分別計(jì)算在各個(gè)模型下的輸出概率,然后進(jìn)行比較,選取輸出概率最大的模型,即代表識(shí)別出的故障類型.
為了驗(yàn)證本文方法在實(shí)際故障診斷中的應(yīng)用情況,筆者進(jìn)行了實(shí)驗(yàn)分析.實(shí)驗(yàn)在一個(gè)故障模擬平臺(tái)上進(jìn)行,包括驅(qū)動(dòng)電機(jī)、軸承座、扭矩加載器和采集設(shè)備,見圖 3.通過更換正常軸承和故障軸承進(jìn)行故障診斷實(shí)驗(yàn).實(shí)驗(yàn)測(cè)試的軸承是深溝球軸承,型號(hào)為6205-2RS,品牌為SKF,其內(nèi)圈直徑為 25,mm,外圈直徑為52,mm,滾子直徑為7.94,mm.軸承故障用電火花加工而成,故障直徑分別為 0.177,8,mm、0.355,6,mm、0.533,4,mm,3 個(gè)等級(jí),深度都為 0.279,4,mm.電機(jī)轉(zhuǎn)速為 1,797,r/min.振動(dòng)信號(hào)由加速度傳感器拾取,分別采集安裝無故障軸承、外圈故障軸承、內(nèi)圈故障軸承、滾子故障軸承時(shí)的振動(dòng)信號(hào),各測(cè)得數(shù)據(jù) 50 組,每組2,048 點(diǎn),采樣頻率為12,000,Hz.
圖3 實(shí)驗(yàn)裝置Fig.3 Experimental device
對(duì)無故障、內(nèi)圈故障、外圈故障、滾子故障4 種狀態(tài)信號(hào) TVAR 建模時(shí),參照文獻(xiàn)[10]中的論述,基函數(shù)選擇傅里葉函數(shù),采用AIC準(zhǔn)則確定模型階數(shù)p.在確定模型階數(shù)前必須先確定擴(kuò)展維數(shù) m,為了保證信號(hào)時(shí)頻譜的高分辨率和計(jì)算速度,一般模型的擴(kuò)展維數(shù)取為 8~12.本次實(shí)驗(yàn)選取擴(kuò)展維數(shù) m=10,依據(jù)AIC 準(zhǔn)則確定軸承無故障、內(nèi)圈故障、外圈故障、滾子故障狀態(tài)下模型階數(shù) p 分別為 11、15、15、12.
4 種類型的振動(dòng)信號(hào)TVAR 建模后,可以得到相應(yīng)的時(shí)頻譜,如圖 4(a)~(d)所示.由時(shí)頻譜中可以看出,不同故障類型軸承的振動(dòng)信號(hào)的頻率成分是不同的,并且軸承存在故障時(shí),時(shí)頻譜中存在時(shí)變性與沖擊性.圖 4(e)為軸承外圈故障信號(hào)時(shí)域波形,圖4(f)為軸承外圈故障信號(hào)經(jīng)傅里葉變換后得到的頻譜.通過圖4(d)~(f)的比較可以明顯看出圖4(d)與4(f)中信號(hào)的頻率成分保持一致,圖 4(d)與 4(e)中頻率成分發(fā)生沖擊的時(shí)間及時(shí)間間隔基本一致.這說明經(jīng) TVAR 建模得到的時(shí)頻譜具有較高的時(shí)頻分辨率.
基于能量法對(duì)時(shí)頻圖進(jìn)行特征提取時(shí),考慮到速度與精度雙重因素,經(jīng)過多次實(shí)驗(yàn)分析,選擇將時(shí)頻譜分為16 部分(過多會(huì)影響HMM 訓(xùn)練速度;過少會(huì)影響測(cè)試的準(zhǔn)確度),求每一部分能量的均值,將均值序列歸一化作為特征向量.通過特征提取,每種狀態(tài)得到 50 組 16 維的特征向量,4 種狀態(tài)共得到 4×50組特征向量.
圖4 4種狀態(tài)的時(shí)頻譜Fig.4 Time-frequency spectra of four states
建模時(shí),采用 4 個(gè)隱狀態(tài)代表無故障軸承、內(nèi)圈故障軸承、外圈故障軸承和滾子故障軸承,分別記作λ1、λ2、λ3、λ4.模型的觀測(cè)序列是 TVAR 建模后得到的 16 維特征向量.初始概率分布矢量、初始狀態(tài)轉(zhuǎn)移概率矩陣和初始觀測(cè)值概率矩陣均由 Matlab 中rand 函數(shù)隨機(jī)取得,然后歸一化.
選用120 組特征向量(每種狀態(tài)30 組)作為訓(xùn)練樣本,訓(xùn)練前需使用Lloyd 算法對(duì)特征向量進(jìn)行標(biāo)量量化.訓(xùn)練算法為Baum-Welch 算法. HMM 的訓(xùn)練過程中,隨著迭代次數(shù)的增加,最大對(duì)數(shù)似然估計(jì)值也在不斷增加,直到達(dá)到收斂誤差為止.訓(xùn)練結(jié)束后,得到4 個(gè)隱狀態(tài)對(duì)應(yīng)的4 個(gè)HMM 識(shí)別模型.圖5 給出了4 種狀態(tài)的HMM 訓(xùn)練曲線,所有狀態(tài)在迭代次數(shù)20 都達(dá)到了收斂,收斂速度快.
圖5 HMM訓(xùn)練曲線Fig.5 HMM training curves
在訓(xùn)練好各狀態(tài) HMM 以后,將剩余的 80 組特征向量(每種狀態(tài)20 組)作為測(cè)試樣本輸入到各狀態(tài)HMM 模型中(測(cè)試前需使用 Lloyd 算法對(duì)特征向量進(jìn)行標(biāo)量量化編碼),在各個(gè)模型下能夠各自輸出 1個(gè)對(duì)數(shù)似然概率值.對(duì)數(shù)似然概率值反映了特征向量與各 HMM 的相似程度,對(duì)數(shù)似然概率值越大,特征向量越接近該狀態(tài) HMM,特征向量屬于使輸出對(duì)數(shù)似然概率值最大的模型所對(duì)應(yīng)的故障類型.識(shí)別算法選用 Viterbi 算法.表 1 給出了部分測(cè)試樣本的測(cè)試結(jié)果,通過得到的對(duì)數(shù)似然概率值,可以準(zhǔn)確識(shí)別出軸承的故障類型;對(duì)于同一種故障類型,故障程度越高,得到的對(duì)數(shù)似然概率值越大.
為了比較基于HMM 和基于BP 網(wǎng)絡(luò)、RBF 網(wǎng)絡(luò)故障診斷模型的識(shí)別效果,將本實(shí)驗(yàn)的數(shù)據(jù)分別用BP 網(wǎng)絡(luò)、RBF 網(wǎng)絡(luò)故障診斷模型來識(shí)別.網(wǎng)絡(luò)的層數(shù)設(shè)置為3 層,輸入層神經(jīng)元個(gè)數(shù)為16 個(gè),輸出層神經(jīng)元個(gè)數(shù)為 4 個(gè).采用與基于 HMM 故障診斷模型同樣的訓(xùn)練數(shù)據(jù)和測(cè)試數(shù)據(jù)分別訓(xùn)練和測(cè)試 BP 網(wǎng)絡(luò)和RBF 網(wǎng)絡(luò),得出診斷結(jié)果.
圖6為本次實(shí)驗(yàn)中 3 種模型的故障識(shí)別正確率比較.由圖6 可見,基于 HMM 故障診斷模型的有效識(shí)別率明顯好于基于 BP 網(wǎng)絡(luò)或 RBF 網(wǎng)絡(luò)故障診斷模型.實(shí)驗(yàn)結(jié)果表明,基于 HMM 故障診斷模型可以成功地診斷旋轉(zhuǎn)機(jī)械的軸承故障,識(shí)別率較高,穩(wěn)定性好.
表1 部分測(cè)試結(jié)果Tab.1 Some testing results
圖6 HMM、BP、RBF識(shí)別有效性比較Fig.6 Diagnosis accuracies of HMM、BP、RBF
(1)TVAR 模型非常適合非平穩(wěn)信號(hào)的處理,時(shí)變參數(shù)體現(xiàn)了信號(hào)頻率的時(shí)變性和沖擊性.TVAR 建模分析得到的信號(hào)時(shí)頻譜的分辨率高,且無交叉項(xiàng)的干擾,計(jì)算速度快.
(2)HMM 模型可以成功地對(duì)旋轉(zhuǎn)機(jī)械的軸承故障進(jìn)行診斷,識(shí)別率高.由于其為動(dòng)態(tài)時(shí)間序列統(tǒng)計(jì)模型,非常適合識(shí)別信息量大、非平穩(wěn)性強(qiáng)、特征重復(fù)再現(xiàn)性不佳的故障信號(hào),在旋轉(zhuǎn)機(jī)械故障診斷中有很好的應(yīng)用前景.
[1]李 強(qiáng). 非平穩(wěn)信號(hào)特征提取新方法和實(shí)用診斷技術(shù)的研究[D]. 天津:天津大學(xué)機(jī)械工程學(xué)院,2003.Li Qiang. Research on New Method of Non-Stationary Signal Feature Extraction and Practical Diagnostic Techniques[D]. Tianjin:School of Mechanical Engineering,Tianjin University,2003(in Chinese).
[2]王奉濤. 非平穩(wěn)信號(hào)故障特征提取與智能診斷方法的研究與應(yīng)用[D]. 大連:大連理工大學(xué)機(jī)械工程學(xué)院,2003.Wang Fengtao. Research on Non-Stationary Signal Feature Extraction and Intelligent Diagnostic Method[D].Dalian:School of Mechanical Engineering,Dalian University of Technology,2003(in Chinese).
[3]陳 慧,張 龍,李 嶷. 基于 TVAR 的非平穩(wěn)工況轉(zhuǎn)子故障診斷技術(shù)研究[J]. 華東交通大學(xué)學(xué)報(bào),2006,23(5):116-118.Chen Hui,Zhang Long,Li Yi. TVAR based fault diagnosis for rotation machine under non-stationary condition[J]. Journal of East China Jiaotong University,2006,23(5):116-118(in Chinese).
[4]Wang Guofeng,Luo Zhigao. Fault identification and classification of rolling element bearing based on timevarying autoregressive spectrum [J]. Mechanical Systems and Signal Progressing,2008,22(4):934-947.
[5]張 龍,熊國良,柳和生,等. 基于時(shí)變自回歸模型與支持向量機(jī)的旋轉(zhuǎn)機(jī)械故障診斷方法[J]. 中國電機(jī)工程學(xué)報(bào),2007,27(9):99-103.Zhang Long,Xiong Guoliang,Liu Hesheng,et al. TVAR and SVM based fault diagnosis method for rotation machine [J]. Proceedings of the CSEE,2007,27(9):99-103(in Chinese).
[6]Abramovich Y I,Spencer N K,Turley M D E. Timevarying autoregressive(TVAR)models for multiple radar observations [J]. IEEE Transactions on Signal Processing,2007,55(4):1298-1311.
[7]Tai A H,Ching Wai-Ki,Chan L Y. Detection of machine failure:Hidden Markov model approach[J]. Computers& Industrial Engineering,2009,57(2):608-619.
[8]Hall M G,Oppenheim A V,Willsky A S. Time varying parametric modeling of speech [J]. Signal Processing,1983,5(3):267-285.
[9]Ocak H,Loparo K A. A new bearing fault detection and diagnosis scheme based on hidden Markov modeling of vibration signals [C]//IEEE International Conference on Acoustics,Speech,and Signal Processing. Salt Lake City,USA,2001,5:3141-3144.
[10]Zhan Y M,Jardine A K S. Adaptive autoregressive modeling of non-stationary vibration signals under distinct gear states(Part 1):Modeling [J]. Journal of Sound and Vibration,2005,286:429-450.
[11]劉 偉. 基于能量的時(shí)頻特征提取研究[J]. 紅外與激光工程,2004,33(3):296-299.Liu Wei. Research on the time-frequency feature extraction based on energy[J]. Infrared and Laser Project,2004,33(3):296-299(in Chinese).