基于有監督增量式局部線性嵌入的故障辨識

李 鋒，田大慶，王家序，楊榮松

目前的旋轉機械故障診斷方法一般是先利用信號處理技術及衍生的特征構造方法對故障信號進行故障特征選擇與優化，再通過人工智能算法或決策融合機制對故障特征進行模式識別［1－2］。這種診斷模式多采用單一/單域的故障特征提取方式，該特征提取方式作用域少、泛化性差，難以全面、準確地捕捉非線性、強耦合、特征表現不確定的故障特征［3］;該診斷模式還常要借助人為分析才能完成故障特征的選擇與優化，故障特征優選質量和模式識別精度取決于用戶的專業知識和現場經驗，精度低、效率低、可靠性差［4］。

為改善現有故障診斷方法的識別精度和計算效率，本文首先構造20個時域和頻域特征參數，從時域和頻域兩個作用域來全面準確地挖掘復雜旋轉機械系統不同部位、不同類型、不同程度故障的特征信息。但高維多征兆域特征集不可避免地會摻雜一些冗余信息和干擾成分，造成不同故障特征集之間具有某種統計相關性，這種不該具有的相關性容易扭曲特征的分布結構甚至嚴重惡化分類性能，因此需要采用維數化簡方法代替人工的特征優選方式來對高維時頻域特征集中狀態敏感的特征進行二次提取，以獲取維數低、敏感性高、獨立不相關的主要特征矢量。流形學習比傳統約簡方法更能體現事物本質，在提取出主要變量的同時還獲得了原始觀測空間的真實結構分布［5］，為對高維、非線性的旋轉機械故障特征進行自動特征約簡和高精度的模式分類提供了更好的解決思路。目前最典型的流形學習方法包括局部保持映射(Locality Preserving Projection，LPP)［6］、鄰域保持嵌入(Neighborhood Preserving Embedding，NPE)［7］、線性局部切空間排列(Linear Local Tangent Space Alignment，LLTSA)［8］、局部線性嵌入(Locally Linear Embedding，LLE)［9］等。但這些方法都是無監督方法，忽視了大量訓練樣本的類信息，因而有損其自身的分類特性;此外，這些方法都需要通過丟棄已獲得的低維流形結構，重新合并原訓練樣本和新樣本并重構權值矩陣及投影模型來實現新增樣本的推廣，無法高效處理實際大批量數據，故以上方法對以快速的解耦和分類為目的的應用來說是非優化的［10］。本文在LLE基礎上研究了有監督增量式局部線性嵌入(Supervised Incremental Locally Linear Embedding，SILLE)這一新流形學習理論，實現了故障特征流形的解耦與分類及新增樣本的增量處理，得到了更高的故障辨識精度和計算效率。為直觀表達辨識結果，本文最后借助新型模式識別技術Morlet小波支持向量機(Morlet Wavelet Support Vector Machine，MWSVM)［4］建立SILLE降維后的特征矢量與故障模式之間的映射關系。

本文提出了基于有監督增量式局部線性嵌入(SILLE)的故障辨識方法:“時頻域特征集－SILLEMWSVM”，用于典型旋轉機械—深溝球軸承的故障診斷和某型號空間軸承的壽命狀態辨識，達到了較高的故障診斷精度、壽命狀態識別精度和計算效率。

1 有監督增量式局部線性嵌入

令 X={x1，x2，…，xN}和 Xnew={xnew1，xnew2，…，xnewM}分別為高維數據空間?D中的N個訓練樣本和M個新增測試樣本。這些樣本都存在本征維數為d(一般d?D)的非線性流形。SILLE的目的就是通過將D維數據映射到?d空間來尋找X和Xnew的低維嵌入式，即d維非線性流形。這里將嵌入空間Rd中的N個訓練樣本表示為Y={y1，y2，…，yN}，M 個新增測試樣本表示為 Ynew={ynew1，ynew2，…，ynewM}。

1.1 SILLE理論

首先對每一個訓練用振動信號樣本xi和新增待測振動信號樣本xnewi(i=1，2，…，M)進行去直流分量和歸一化處理，使它們具有零均值，初步消除噪聲干擾。

SILLE再結合類標簽信息尋找每個訓練樣本xi∈X的k個同類最鄰近點。這里采用放大異類樣本間的歐幾里得距離，并保持同類樣本間的距離不變的方法來拉大異類樣本之間的距離，實現同類樣本的相對靠攏，即:

以保證最鄰近點以同類點居多。其中，S=‖xi－xj‖為未考慮類標簽信息時的原始歐幾里得距離，max(S)=maxi，j‖xi－ xj‖為樣本間的最大距離，S'為融入類標簽信息的距離。若xi和xj與屬于異類，則δ(xi，xj)=0;否則，δ(xi，xj)=1。式(1)中，α 控制類信息的融入程度。當α=0時，同于無監督LLE;當α=1時，就是全監督LLE。若α在0與1之間變動，就形成了局部有監督LLE。

由于SILLE結合訓練樣本的類標簽信息和類信息融入程度控制參數α來調節樣本間的局部鄰域結構(即形成了有監督的學習機制)，因此可以有效分離高維、非線性、強耦合(強相關性)的異類故障特征，強化了同類樣本的聚集性和異類樣本的互斥性。相比無監督(即流形局部鄰域結構完全取決于原始樣本的幾何關系，無類判別信息指導)的LLE，SILLE的解耦和分類性能顯然更為優越。

通過最小化xi的局部重構誤差來求解xi的最近鄰點的重構權值。即:

該局部重構誤差為xi與其重構值距離的平方，并服從約束:

為計算重構權值，引入k×k維局部協方差矩陣Qi(k為近鄰點數目)，即Qjsi=(xi－xj)T(xi－xs)。使用逆局部協方差矩陣和拉格朗日乘數法來求解以上有約束最值問題，獲得最優權值如下［9］:

由wij可構造N×N維稀疏權值矩陣W=［wij］N×N。

SILLE基于所得到的權值矩陣W來計算可保持X的局部幾何結構的最佳低維嵌入式Y。這對應于最小化以下代價函數:

基于 M=［Mij］N×N，式(4)可表示成二次方程式形式，即 φ(Y)=根據Rayleigh-Ritz理論，求式(4)所示最小值問題的解yi可以通過尋求矩陣M的d個最小非零特征值(即將M的d+1個最小特征值中最趨近于0的特征值剔除之后所保留下的d個最小特征值)所對應的特征向量得到。剔除M的零特征值所對應的特征向量目的是保證嵌入式Y={y1，y2，…，yN}具有零均值屬性，以消除嵌入式中隨機尺度因子波動所帶來的干擾。非零特征值所對應的特征向量經過正交化處理之后其對應元素的累加之和必然為零，因此容易實現零均值的目的。

SILLE最后采用局部線性投影算法來將高維新增樣本映射到訓練樣本的嵌入空間。實現流程如下:

(1)對于新增高維樣本xnewi，根據歐幾里得距離尺度在訓練樣本集中尋找其k個最近鄰點x1，x2，…，xk。

(2)尋找滿足Y=AX的映射矩陣A(此時精確解可能不存在)，其中 A=［a1，a2，…，ad］T∈?d×D，X=［x1，x2，…，xk］，Y=［y1，y2，…，yk］為 xnewi的近鄰點 X相對應的嵌入式。A的每一個基向量aj可以通過求解以下線性最小二乘回歸問題得到:

若XT為列滿秩矩陣，則aj可以求解如下:

若XT不是列滿秩矩陣(例如近鄰點數目k小于數據維數D)，則問題(6)就成為不適定問題，因為此時可能存在無窮多解。為解決不適定問題，使局部線性投影算法能更好地推廣應用于新增樣本，使用嶺回歸法來求解 aj，即求解以下正則化最小二乘目標函數［10－11］:

其中:β＞0為正則化參數。則aj求解如下:

式中:I表示D×D單位矩陣。嶺回歸法就是要將xnewi的每一個近鄰點線性投影到d維嵌入空間，同時將噪聲從該低維空間中剔除。故β值由噪聲量值σ2確定。使用特征分析法來估計β值如下:

式中:λi表示矩陣XXT的非零特征值(λi按從大到小的順序排列)。

顯然，SILLE中的局部線性投影算法使新樣本加入時，無需丟棄已獲得的訓練樣本d維嵌入式，無需重復計算近鄰點、重構權值矩陣和投影模型，而能有效利用訓練樣本的低維嵌入式來對新增樣本進行增量處理，也即SILLE建立了類似人工神經網絡和支持向量機的獨立的學習訓練模型，當有新的待測樣本時，SILLE可直接由基于訓練樣本的學習訓練模型來對新增樣本進行分類，而不用重新對訓練樣本進行學習訓練;而LLE沒有提供明確的將訓練樣本和新增樣本從高維空間向低維空間化簡的映射功能，LLE只是粗略地將矩陣M的d個最小非零特征值所對應的特征向量視為訓練樣本的d維嵌入式(即化簡后的d維特征向量)，若存在新增測試樣本，LLE需要丟棄已得到的訓練樣本d維嵌入式，將原訓練樣本和新增樣本合并以重新計算權值矩陣W、投影模型M和d維嵌入式，因此LLE未建立獨立的訓練模型，無法對新增大批量實測數據進行高效處理。SILLE的批處理能力明顯優于LLE，SILLE可以提高所提故障辨識方法的計算速率。

(3)xnewi的d維嵌入坐標ynewi計算如下:

ynewi=Axnewi，也即 Ynew=AXnew。

1.2 約簡維數、近鄰點數目和距離參數估計

①約簡維數d的估計方法如下:首先對每一個訓練樣本xi的局部協方差矩陣Qi(m，n)=(xi－xm)T×(xi－xn)進行特征分析;然后為每一個訓練樣本xi計的約簡維數 di，其中 λj為Qi的非零特征值;最后選擇di中的最大值作為總體上可以接受的d值。②近鄰點數目k和距離參數α可以由十折交叉驗證法來確定:首先以1為步長在區間［2，20］內隨機選取一個數作為k值，同時以0.05為步長在區間［0，1］內隨機選取一個數作為α值。然后對選取的各對參數(k，α)執行十折交叉驗證操作。即:將原始訓練樣本X={x1，x2，…，xN}隨機地等分為10個子集。在每一折(次)故障辨識運算中，將其中一個子集作為測試集，而剩余子集共同作為訓練集。10折(次)運算下來，每個子集都已做過一次測試集。將這10折(次)平均分類結果記錄下來。最后將得到最佳分類結果的(k，α)作為最優參數。

2 故障辨識方法“時頻域特征集－SILLEMWSVM”

故障辨識方法“時頻域特征集－SILLE－MWSVM”實現流程如圖1，步驟如下:

①對每個訓練樣本和測試樣本信號構造10個時域特征參數和10個頻域特征參數。②將10個時域特征參數和10個頻域特征參數組合得到20維時頻域特征集。③將作為訓練樣本的20維時頻域特征集輸入SILLE進行學習，求得訓練樣本的d維嵌入式(特征矢量)。1≤d＜k，最優化簡維數d、近鄰點數目k和距離參數α的估計方法如1.2節所述。④訓練好的SILLE再用局部線性投影法對測試樣本進行維數化簡，得到測試樣本的d維特征矢量。⑤將訓練樣本的d維特征矢量輸入MWSVM對其進行訓練。⑥用訓練好的MWSVM判別測試樣本d維特征矢量的故障位置、類型或程度。

圖1 辨識方法“時頻域特征集－SILLE－MWSVM”流程Fig.1 Implementation process of the method“Time-frequency domain feature set-SILLE-MWSVM”

顯然，SILLE在本故障辨識方法中起到承上啟下的作用，它自動銜接時頻域特征集和MWSVM，并具有優良的維數化簡、模式分類和增量處理特性，是所提故障辨識方法實現較高辨識精度和計算效率的關鍵技術。

當旋轉機械出現故障時，時域信號的幅值和概率分布將會發生變化;信號中的頻率成分、不同頻率成分的能量，以及頻譜的主能量譜峰位置也將發生改變。因此，通過描述信號時域波形和頻譜中能量的大小和分布等，可以反映振動信號的時域和頻域信息，從而指示故障的出現和不同故障之間的差異性［3］。因此為獲取更多早期故障信息，這里構造了10個時域特征參數和10個頻域特征參數(如表1所示)，以較全面準確地反映旋轉機械的運行狀態。

表1中，時域特征參數c1和c3～c5表征時域信號的幅值和能量大小;c2和c6～c10反映時域信號的時間序列分布情況。頻域特征參數c11表征頻域振動能量的大小;c12～c13、c15和c18～c20表征頻譜的分散或集中程度;c14和c16～c17反映主頻帶位置的變化。通過以上10個時域特征參數和10個頻域特征參數的組合來全面準確地描述旋轉機械不同部位、不同類型、不同程度故障的特征信息。

表1 時域、頻域特征參數Tab.1 Characteristic parameters in time domain and frequency domain

所用到的Morlet小波支持向量機(MWSVM)是一種新型支持向量機(SVM)，其MWSVM小波核函數通過平移伸縮即能生成二次可積空間L2(R)上的一組完備正交基，可逼近L2(R)空間上任意函數，故MWSVM具有比一般SVM更好的自適應分類決策力［5］。

3 實例分析

3.1 深溝球軸承故障診斷

本文以深溝球軸承正常狀態以及不同部位、不同程度故障的同步診斷實例驗證本故障辨識方法的有效性。實驗采用6205－2RS型深溝球軸承，軸承內徑25 mm，外徑52 mm，厚度15 mm。軸承實驗由電動機、扭矩傳感器/譯碼器、測力計和電器控制裝置組成，由電機帶動輸入軸，轉速控制在1 772 r/min，輸出軸帶動負載。在3個軸承的外圈、內圈、滾動體上分別加工寬0.178 mm，深0.28 mm 的小槽模擬軸承外圈、內圈、滾動體局部輕微裂紋，在另外2個軸承的內圈分別加工寬0.356 mm和0.533 mm(深度尺寸不變)的小槽模擬內圈中度裂紋和嚴重裂紋，在另1個軸承的滾動體上加工寬0.533 mm(深度尺寸不變)的小槽模擬滾動體嚴重裂紋。通過壓電式加速度傳感器、電荷放大器、數據采集器分別采集以上6類故障和正常狀態的振動數據各50組，每組數據長度為0.1 s，從中隨機抽取20組用于樣本訓練，另外30組作為測試樣本，采樣頻率48 kHz。

對每組訓練樣本和測試樣本截取4 096個數據點用于測試樣本的故障辨識，數據截取長度的選定遵循既要基本覆蓋各類故障的特征頻帶(即振動周期)又不增加辨識方法過多計算量的原則。辨識流程如第2節所述。時域、頻域特征參數總個數為20，也即SILLE的輸入特征維數為20。MWSVM的參數設置如下:懲罰因子 γ =1，核參數 c=1，ω0=1.75［4］。故障模式與MWSVM期望輸出的對應關系設定為:正常狀態→1、外圈輕度裂紋→2、內圈輕度裂紋→3、滾動體輕度裂紋→4、內圈中度裂紋→5、內圈嚴重裂紋→6、滾動體嚴重裂紋→7。以下從四個層次評估本故障辨識方法的性能:

(1)在本辨識方法的維數化簡和模式識別兩個環節不做任何變動的前提下，將時頻域特征集的辨識精度與時域特征參數、頻域特征參數等單一/單域方法進行對比，對比結果如表2所示。顯然，時頻域特征集所達到的測試樣本辨識精度比二者的單獨作用都更好，表明時頻域特征集全面挖掘軸承不同部位、不同程度故障特征信息的能力要優于單一/單域特征提取方法。

(2)將SILLE的維數化簡效果與LPP、LLTSA、LLE對比，對比結果如表3。

顯然，采用有監督學習機制(即由訓練樣本的類標簽信息來優化流形局部鄰域結構和重構權值矩陣)的SILLE，其對軸承七種故障模式(包括正常狀態)的辨識結果明顯好于無監督的LPP、LLTSA和LLE的辨識效果。表4是本辨識方法“時頻域特征集 －SILLEMWSVM”的故障辨識結果詳表。

表2 三種特征提取方法的故障辨識精度對比%Tab.2 Comparison of fault identification accuracy achieved by 3 feature extraction methods%

表3 SILLE維數化簡后的故障辨識精度和其它三種流形學習算法以及不進行維數化簡的辨識精度對比%Tab.3 Comparison of identification accuracy achieved by dimension reduction with SILLE and another three manifold learning algorithms and no dimension reduction%

(3)再將SILLE維數化簡后的故障辨識精度與不進行維數化簡，將高維時頻域特征集直接輸入MWSVM的辨識結果進行對比，結果見表3。由于高維時頻域特征集不可避免地會摻雜一些冗余信息和干擾成分，造成不同故障的高維時頻域特征集之間具有某種統計相關性，因此若不進行維數化簡處理，直接將高維時頻域特征集輸入MWSVM進行故障辨識，正如表3第6行結果所示，異類特征集之間的統計相關性造成了異類故障特征分布的混疊與耦合，進而引起故障辨識精度的下降。而經SILLE維數化簡后的辨識精度處于高水平，表明SILLE在維數化簡的同時，在基于類標簽的有監督機制的調節下消除了不同故障信號之間的相關性，保留了高維故障信息中的低維不相關主成分，提高了對正常狀態和6類不同層次故障的區分度。

(4)將原測試樣本作為第二批新增樣本再分別輸入到已訓練好的基于SILLE的辨識方法中和基于LPP、LLTSA、LLE的其它辨識方法(如表5所示)中，對比這四種方法的計算速率以檢驗SILLE的增量處理性能，對比結果如表5。顯然，由于SILLE能有效利用前面的處理結果，不必對訓練樣本的嵌入式(包括參數d、k、α的估計過程在內)做重復計算，只單獨對第二批測試樣本做投影計算，故基于SILLE的辨識方法對第二批樣本的分析效率大幅提高，由對第一批樣本的處理時間64 s降至現在的6 s，這表明了SILLE增量處理的快速性;而 LPP、LLTSA和 LLE因不具有類似 SILLE和MWSVM的增量處理能力，故基于LPP、LLTSA、LLE的辨識方法對第二批樣本的處理時間依然較長，效率較低。

表4 故障辨識方法“時頻域特征集－SILLE－MWSVM”故障識別結果詳表Tab.4 Concrete fault identification results of the method“Time-frequency domain feature set-SILLE-MWSVM”

表5 四種故障辨識方法在相同測試樣本數條件下的計算時間對比sTab.5 Computation time comparison of 4 fault identification methods under the condition of same test sample numbers s

以上軸承故障辨識實例充分驗證了時頻域特征集用于故障特征的全面準確挖掘的有效性和SILLE用于自動維數化簡、模式分類與增量處理的有效性，也即驗證了故障辨識方法“時頻域特征集－SILLE－MWSVM”具有較高的故障辨識精度和計算效率。

3.2 軸承壽命狀態識別

通過對某型號空間軸承不同的壽命階段(即對早期故障出現之前的正常狀態按相同時間間隔劃分階段)進行識別來進一步驗證本故障辨識方法的泛化性。本實驗中同時監測了同時開始運轉的3個同型號的空間軸承(分別編號為1、2、3)，3個軸承的運行時轉速恒為1500 r/min，并都加載3 kg軸向載荷，采集振動數據時轉速恒為1000 r/min，采樣頻率25.6 kHz。所監測到的3個軸承的壽命狀態有T1(運行不到1天)、T2(已運行15天)、T3(已運行30天)、T4(已運行45天)、T5(已運行60天)、T6(已運行75天)、T7(已運行90天)等。采集以上7種壽命狀態的振動數據樣本各25組(其中1號軸承8組、2號軸承8組、3號軸承9組)，將樣本打亂后從中隨機抽取10組用于樣本訓練(結果顯示含1號軸承3組、2號軸承3組、3號軸承4組)，另外15組作為測試樣本，訓練和測試樣本長度都為4 096個數據點。壽命階段模式與MWSVM的期望輸出的對應關系設為:T1→1、T2→2、T3→3、T4→4、T5→5、T6→6、T7→7。本故障辨識方法的壽命狀態識別結果如表6所示。

表6 故障辨識方法對7個壽命階段的識別精度%Tab.6 Recognition accuracy of the proposed fault identification method for seven life states

結果表明所提故障辨識方法對某空間軸承7個壽命階段的識別結果正確，并達到較好的精度，這說明該故障辨識方法泛化性較好，可推廣應用于旋轉機械壽命狀態識別或早期故障辨識，這為旋轉機械的可靠性評估和剩余壽命預測提供了可借鑒的理論與方法。

4 結論

提出了基于有監督增量式局部線性嵌入(SILLE)的故障辨識方法:“時頻域特征集－SILLE－MWSVM”。①時頻域特征集較為全面準確地表征了不同故障的特征信息;②SILLE利用類標簽調節流形局部幾何結構和設計重構權值矩陣，實現了有監督流形解耦與分類機制，并采用局部線性投影計算新增故障樣本的嵌入映射，發揮了自動簡化高維時頻域特征集、高精度區分故障類別、增量處理新樣本等關鍵作用;③所提故障辨識方法集成了時頻域特征集在故障特征的全面提取、SILLE在信息化簡和MWSVM在模式識別上的優勢，具有較高的故障辨識精度、壽命識別精度和計算效率。

［1］ Zhou G N，Hou S L，Bo L，et al.Fault diagnosis of rotor system based on wavelet energy immune recognition［J］.Electronics Optics and Control，2010，17(6):81 －84.

［2］李 鋒，湯寶平，董紹江.基于正交鄰域保持嵌入特征約簡的故障診斷模型［J］.儀器儀表學報，2011，32(3):621－627.LI Feng， TANG Bao-ping， DONG Shao-jiang. Fault diagnosis model based on feature compression with orthogonal neighborhood preserving embedding［J］.Chinese Journal of Scientific Instrument，2011，32(3):621 －627.

［3］雷亞國，何正嘉，訾艷陽.基于混合智能新模型的故障診斷［J］.機械工程學報，2008，44(7):112－117.LEI Ya-guo，HE Zheng-jia，ZI Yan-yang.Fault diagnosis based on novel hybrid intelligent model［J］.Chinese Journal of Mechanical Engineering，2008，44(7):112－117.

［4］ Tang B P，Li F，Qin Y.Fault diagnosis model based on feature compression with orthogonal locality preserving projection［J］.Chinese Journal of Mechanical Engineering，2011，24(5):897－904.

［5］Huang H，Li J W，Liu J M.Enhanced semi-supervised local fisher discriminant analysis for face recognition［J］.Future Generation Computer Systems，2012(28):244 －253.

［6］ Cai D， He Xiao-fei， Han J W， et al. Orthogonal laplacianfaces for face recognition［J］.IEEE Trans.Image Process，2006，15(11):3608 －3614.

［7］ Liu X M， Yin J W， Feng Z L，et al.Orthogonal neighborhood preserving embedding for face recognition［C］//2007 IEEE International Conference on Image Processing，ICIP 2007，New York，USA，2007:133－136.

［8］ Zhang T H，Yang J，Zhao D L，et al.Linear local tangent space alignment and application to face recognition［J］.Neurocomputing，2007，70:1547－1553.

［9］ Samue K，Martin D L.Face detection in gray scale images using locally linear embeddings［J］.Comput.Vision Image Understanding，2006:1－20.

［10］ Li B W，Zhang Y.Supervised locally linear embedding projection(SLLEP)for machinery fault diagnosis［J］.Mechanical Systems and Signal Processing，2011，25:3125－3134.

［11］ Tikhonov A N.Solutions of Ill-Posed problems［M］.Wiley，New York:1997.