李宏坤，陳禹臻，張志新，周 帥

(1.大連理工大學(xué) 機械工程學(xué)院，大連 116023;2.大連理工大學(xué) 工業(yè)裝備結(jié)構(gòu)分析國家重點實驗室，大連 116023)

時頻分析方法是對非平穩(wěn)、非線性信號進行分析的有效方法。機械振動所構(gòu)成的時頻圖像含有描述機械工作狀態(tài)的特征信息［1］。為高效率地處理存放于像素矩陣中的數(shù)據(jù)信息，一般采取將矩陣進行分解的方法。分解后可以對原始矩陣中存放的大量數(shù)據(jù)進行壓縮和概括，進而對壓縮矩陣進行特征提取，得到機械設(shè)備的狀態(tài)信息。非負矩陣分解(Non-negative Matrix Factorization，NMF)［2］是目前國際上新的矩陣分解方法，并在理論和應(yīng)用方面已相對成熟。較完整的NMF理論由Lee等［2］提出，并應(yīng)用于人臉識別，這引起了科學(xué)界的廣泛關(guān)注。常用的矩陣分解方法有:主成分分析(PCA)、矢量量化(VQ)等。在所有這些方法中，原始的大矩陣V通過V=WH被近似分解為低秩的形式。這些方法的共同特點是:計算機的“認知”過程是基于數(shù)據(jù)整體的過程，即因子W和H在線性組合時可進行加運算，亦可進行減運算。這樣即使輸入的初始矩陣元素是全正的，分解后的結(jié)果也不能保證非負性。然而負值元素在實際問題中往往是沒有意義的［3］。例如圖像數(shù)據(jù)中不可能有負值的像素點，在文檔統(tǒng)計中負值也是無法解釋的。因此，探索矩陣的非負分解方法一直是很有意義的研究問題。正是如此，NMF方法才得到人們的廣泛關(guān)注。非負矩陣分解以計算機的“認知”過程是基于部分的過程為前提，通過添加“線性組合時只允許加運算不允許減運算”的限制條件，保證了分解結(jié)果的可解釋性。與PCA和VQ相比，NMF數(shù)學(xué)實現(xiàn)簡便、分解形式和分解結(jié)果上的可解釋、以及占用存儲空間少［4］，更重要的是非負限制與聯(lián)系部分形成整體的直覺概念是協(xié)調(diào)的，這與智慧生命對事物感知的生理過程是相符的［2］。目前它已廣泛應(yīng)用于諸多領(lǐng)域。

蔡蕾等［5］提出了基于圖像識別的設(shè)備運行狀態(tài)識別方法，并初步研究了稀疏性非負矩陣分解中各參數(shù)對識別率的影響。本文根據(jù)設(shè)備故障診斷中實際存在的問題，研究采用Hilbert時頻譜構(gòu)建圖像，應(yīng)用非負矩陣分解和主元分析對其進行特征向量提取，完成狀態(tài)的識別。并以滾動軸承在不同狀態(tài)下的識別為例驗證方法的有效性。

1 非負矩陣分解及PCA降維

1.1 NMF定義及算法實現(xiàn)

非負矩陣分解使計算機對整體的認知是基于其部分的認知。其定義為［2，6］:給定一個n×m圖像矩陣V，其每一列都是一組圖像數(shù)據(jù)，并且包含n個非負像素值。為了估計原始圖像矩陣，需要尋找兩個非負矩陣因子 W∈Rn×r和 H∈Rr×m，使得:

式中:W中的r列數(shù)據(jù)構(gòu)成基本圖像序列組，H中的每一列被認為是一個由構(gòu)成V中每一個圖像的基本圖像的線性組合系數(shù)組成的編碼且與在V中的每一個圖像序列一一對應(yīng)。因數(shù)分解中r值的選擇通常滿足:

WH的乘積即可被認為是V中數(shù)據(jù)的壓縮形式。線性組合只允許加法運算而不允許減法運算的約束將導(dǎo)致分解具有非負性。

NMF的求解問題實際上是一個最優(yōu)化問題，并利用乘性迭代的方法求解W和H。以W和H的非負性為初始條件，通過使目標函數(shù)F收斂至局部最大，得到非負矩陣V的近似分解［7］。本文使用的目標函數(shù)定義為 KL 散度(Kullback-Leibler divergence)［8］:

其中:分解矩陣W和H的更新法則如下:

當該目標函數(shù)達到局部最大值時，得到非負矩陣V的最優(yōu)近似分解。分解的精確性通過控制Viu與(WH)iu的商值來參與更新過程。更新法則不但維持著W和H的非負性，而且能有效杜絕W列上數(shù)據(jù)的“退化”現(xiàn)象。

1.2 PCA降維簡介

PCA方法是一種基于二階統(tǒng)計的數(shù)據(jù)分析方法。該方法通過對各個變量之間相關(guān)關(guān)系進行分析，用一組較少的、互不相關(guān)的、盡可能多地保留原來復(fù)雜變量所反映信息的主元變量來代替原來較多的變量，從而實現(xiàn)對數(shù)據(jù)的降維。這種方法可以有效地找出數(shù)據(jù)中最“主要”的元素和結(jié)構(gòu)，去除噪音和冗余，將原有的復(fù)雜數(shù)據(jù)降維，揭示隱藏在復(fù)雜數(shù)據(jù)背后的簡單結(jié)構(gòu)。它的優(yōu)點是簡單，而且無參數(shù)限制，可以方便的應(yīng)用于各個場合［9－10］。

2 基于NMF和PCA的狀態(tài)識別

傅里葉變換對非平穩(wěn)與非線性信號進行頻譜分析往往得不到理想的效果。因此，為精確的描述這些非穩(wěn)定信號的局部信息，時頻分析方法得到廣泛的重視和研究。時頻方法得到的結(jié)果往往可以用圖像的形式(如輪廓圖像、灰度圖像、浮雕圖像)表達出來。因此故障診斷實質(zhì)上是對圖像的分類:如果圖像分類正確，則可正確的識別故障。關(guān)于對圖像的分類，傳統(tǒng)的方法是人工判斷，這將在故障診斷中引入主觀因素。而用計算機自動處理(如智能識別系統(tǒng))時，往往由于圖像數(shù)據(jù)的尺寸和維數(shù)太大而無法進行處理。因此，對圖像在其信息盡可能完整的條件下進行降維是非常必要的。

由于圖像像素的非負性，傳統(tǒng)的降維方法如主成分分析(PCA)、矢量量化(VQ)等方法的應(yīng)用效果往往不是很理想。而加入非負限定條件的NMF能很好的解決上述問題。根據(jù)設(shè)備故障診斷的實際特點，本文采用Hilbert時頻譜構(gòu)建時頻圖像，應(yīng)用非負矩陣分解與主元分析方法進行特征提取和狀態(tài)識別。本方法的流程圖如圖1所示。

圖1 基于NMF和PCA的設(shè)備狀態(tài)識別方法流程圖Fig.1 Flow chart for classification based on NMF and PCA

此方法能夠有效對包含設(shè)備故障信息的振動信號進行表征。采用NMF方法能夠獲取其主要特征向量從而提高識別的可靠性，PCA能夠進行有效降維，便于實際設(shè)備故障診斷的應(yīng)用。

3 實驗測試

3.1 實驗數(shù)據(jù)

本文采用的數(shù)據(jù)來自于大連理工大學(xué)振動工程研究所的多功能故障診斷試驗臺進行的滾動軸承故障模擬實驗。模擬的工作狀態(tài)主要有正常、外環(huán)磨損、內(nèi)環(huán)磨損和滾動體磨損四種。每個狀態(tài)共采集了125組數(shù)據(jù)，一組數(shù)據(jù)包含4個周期，即每個狀態(tài)包含了500個周期的數(shù)據(jù)。滾動軸承的型號為N205，轉(zhuǎn)速為1 200 r/min。關(guān)于本實驗對振動信號采集的要求、獲得單周期信號Yi(t)的方法以及循環(huán)平均以得到去一階信號yi(t)的過程見文獻［11］。

3.2 Hilbert時頻譜分析

Hilbert時頻譜是利用瞬時頻率進行信號時頻分析的方法［12］。將濾去一階的信號yi(t)進行Hilbert時頻譜分析。為了削弱噪聲的干擾、突出故障信息，可以對時頻譜進行平均，本文采用5幅Hilbert圖像平均的方法。這樣每種情況將得到100組Hilbert時頻譜圖像。圖2為平均后不同狀態(tài)下的Hilbert時頻譜。

圖2 不同狀態(tài)的Hilbert時頻圖像Fig.2 Different Hilbert Time-frequency image

由圖2可知，各工作狀態(tài)下的Hilbert時頻圖像并無明顯可區(qū)分特點，因此單單依靠Hilbert時頻譜是無法進行狀態(tài)識別的。鑒于Hilbert時頻圖像含有豐富的狀態(tài)信息，為提取有用信息，需對其進行進一步的圖像處理。降低;r值選取過大時，則可能帶入多余的冗余信息和噪聲。當r值為45時，平均識別成功率達到最大［13］。因此本文r選擇45。則每種情況將得到100組維數(shù)為45的特征向量ei。

4 基于NMF和PCA的狀態(tài)識別

4.1 特征向量的提取

在圖像處理之前，要對Hilbert時頻圖像(為201×479的矩陣)進行真彩化和歸一化處理(這樣可以減小特征向量不同量級之間差異的影響，提高訓(xùn)練和預(yù)測的準確性)。

將真彩化和歸一化處理后的矩陣(201×479)表示成一個96279維的列向量。每種狀態(tài)(正常、外環(huán)磨損、內(nèi)環(huán)磨損和滾動體磨損)都有100個圖像數(shù)據(jù)。由100組96 279維的列向量構(gòu)成96279×100矩陣V作為非負矩陣。經(jīng)非負矩陣分解后，得到矩陣因子W和H，滿足V≈WH。W的列向量可認為是組成原圖像矩陣V的基本圖像。將代表信號圖像的96 279維列向量vi與W進行線性組合，即ei=WTvi。這樣，r維向量ei即可認為是該幅信號圖像vi的特征向量。其中r的值滿足式(2)，即 r＜99.896。

非負矩陣分解r值的選取至關(guān)重要。r為描述原始圖像矩陣基本圖像的個數(shù)以及得到的特征向量的維數(shù)。r值選取過小時，得到的基本圖像不足以描述原始圖像矩陣，因而向基投影后部分信息丟失，導(dǎo)致識別率

4.2 利用歐式距離分離各種狀態(tài)

歐式距離是圖像識別中經(jīng)常采用的方法，因此本文首先采用歐式距離對不同狀態(tài)下的故障進行判別。從上述四種情況的100個特征向量中隨機各選一個向量作為參考向量，分別命名為 ball，inner，normal，outer。將特征向量組矩陣NMF(NMF矩陣由45×400構(gòu)成，其中1～100列、101～200列、201～300列、301～400列分別是滾動體磨損、內(nèi)環(huán)磨損、正常、外環(huán)磨損的特征向量)的每一列與上述四個參考向量分別作歐式距離。當對應(yīng)情況的歐式距離最小時則識別成功，圖3為識別結(jié)果。

圖3 歐式距離識別結(jié)果Fig.3 Classification result by using euclidean distance

從圖3可以看出，除對正常情況識別成功以外，其余都不理想。這說明雖然特征向量的維數(shù)較大時，可以較全面的反應(yīng)原始圖像的特征信息，但維數(shù)的增加并不利于利用特征向量進行狀態(tài)識別。所以在狀態(tài)識別之前對特征向量進行降維變得非常必要。

4.3 基于PCA的降維處理

由于PCA的諸多優(yōu)點(如前所述)，本文采用PCA方法對特征向量進行降維處理。經(jīng)過PCA處理，將NMF處理得到的400組維數(shù)為45的特征向量ei降維成維數(shù)為3的特征向量ei'。通過上述的數(shù)據(jù)處理，得到正常、外環(huán)磨損、內(nèi)環(huán)磨損和滾動體磨損四種情況分別對應(yīng)的特征向量組。每種情況對應(yīng)100組維數(shù)為3的特征向量。將降維后的特征向量組矩陣命名為PCA(PCA矩陣由3×400構(gòu)成，其中1～100列、101～200列、201～300列、301～400列分別是滾動體磨損、內(nèi)環(huán)磨損、正常、外環(huán)磨損降維后的特征向量)。本文將一個三維特征向量當作三維空間中的一個點并利用三維坐標系來分離不同的狀態(tài)，圖4、圖5為特征向量在三維坐標系上的反映。可以看出，正常、外環(huán)磨損、內(nèi)環(huán)磨損和滾動體磨損四種狀態(tài)被識別的非常清楚。正常、外環(huán)磨損和滾動體磨損的分布區(qū)域依次鏈狀分布，而內(nèi)環(huán)磨損的點則成環(huán)狀圍繞在鏈的周圍。因而解決了只采用Hilbert時頻譜識別故障狀態(tài)不理想的問題，也解釋了歐式距離無法識別狀態(tài)的原因。

因此經(jīng)過主元分析之后可以有效獲取識別向量作為參數(shù)進行設(shè)備的狀態(tài)識別，此方法可以為以后的狀態(tài)辨識研究打下堅實的基礎(chǔ)。

5 結(jié)論

本文研究了基于NMF的時頻圖像特征提取與PCA降維進行振動信號特征提取與狀態(tài)識別的方法，并采用滾動軸承信號的故障診斷進行方法的有效性驗證。通過對比，此方法在識別滾動體磨損等沖擊現(xiàn)象不明顯的故障狀態(tài)時比Hilbert時頻譜更為有效。此方法的研究將有助于設(shè)備故障診斷和預(yù)知維修的發(fā)展。

［1］ Newton R E.Fragility assessment theory and practice［R］.Monterey Research Laboratory， Inc. Monterey，Califonia，1968.

［2］ Lee D D，Seung H S.Learning the parts of objects by nonnegative matrix factorization ［J］.Nature，1999，(401):788－791.

［3］ Paatero P，Tapper U.Positive matrix factorization:a nonnegative factor model with optimal utilization of error estimates of data values［J］.Environmetrics，1994，5(2):111－126.

［4］ 汪 鵬.非負矩陣分解:數(shù)學(xué)的奇妙力量［J］.計算機教育，2004(10):38－40.

［5］ 蔡 蕾，朱永生.基于稀疏性非負矩陣分解和支持向量機的時頻圖像識別［J］.自動化學(xué)報，2009，35(10):1272－1277.

［6］ Lee D， Seung H. Algorithmsfornon-negative matrix factorization［M］.USA:NeuralInformationProcessing Systems(NIPS).Denver，CO，2000，7.

［7］ Lee D D，Seung H S.Algorithms for non-negative matrix factorization［C］.Advances in Neural Information Processing Systems.Cambridge:MTT Press，2000:556－562.

［8］ 李 樂，章毓晉.非負矩陣分解算法綜述［J］.電子學(xué)報，2008，36(4):737－743.

［9］ Jolliffe I T.Principal component analysis(2nded)［M］.New York:Springer，2002.

［10］ He Q B，Yan R Q，Kong F R，et al.Machine condition monitoring using principal component representations［J］.Mechanical Systems and Signal Processing，2009，23(2):446－466.

［11］ 李宏坤，周 帥，黃文宗，基于時頻圖像特征提取的狀態(tài)識別方法研究與應(yīng)用［J］.振動與沖擊，2010，29(7):184－188.

［12］ Huang N E，Shen Z，Long S R，et al.The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis［J］.Proceedings of The Royal Society of London Series A-Mathematical Physical and Engineering Sciences，1998，454:903－ 995.

［13］ 宦若虹，楊汝良.基于小波域NMF特征提取的SAR圖像目標識別方法［J］.電子與信息學(xué)報，2009，31(3):588－591.