一種基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法

1、權(quán)利要求書發(fā)明名稱：一種基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法發(fā)明人：劉宇（身份證號碼、張凡、陳初杰、李彥鋒、楊圓鑒、米金華、黃洪鐘1. 基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法，如圖1所示，其步驟包括：步驟1：采集被監(jiān)測對象的原始振動數(shù)據(jù)，并劃定訓(xùn)練樣本集和測試樣本集；步驟2：利用集合經(jīng)驗(yàn)?zāi)B(tài)分解方法（Ensemble Empirical Mode Decomposition，EEMD）和希爾伯特變換等信號處理方法從各樣本集中提取指定時(shí)頻域特征；步驟3：利用基于距離的特征評估方法從步驟2中得到的特征集中選擇出少數(shù)敏感特征，并且計(jì)算這些

2、敏感特征的分類貢獻(xiàn)率；步驟4：利用核密度估計(jì)計(jì)算訓(xùn)練樣本中不同故障類別樣本集關(guān)于步驟3中提取的敏感特征的概率密度函數(shù)，并計(jì)算加入一個(gè)未知故障類別待測樣本后各類樣本集新的概率密度函數(shù)；步驟5：求出在選定的特征描述下，訓(xùn)練樣本中各類故障樣本集原始概率密度函數(shù)和加入一個(gè)待測樣本后的新概率密度函數(shù)兩者的K-L散度值；步驟6：計(jì)算集成K-L散度，并通過集成K-L散度的大小判斷待測樣本的故障類別。2. 根據(jù)權(quán)利要求1所示的基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法，其特征在于，所述步驟3不但從原始特征集中選擇出少量敏感特征，而且將選擇出來的敏感特征按照其貢獻(xiàn)率區(qū)別對待，使特征選擇更加合理。對于選擇出

3、的第j個(gè)特征，其貢獻(xiàn)率計(jì)算公式為：，（1）其中，為在第j個(gè)特征描述下，各類樣本集類間距離與類內(nèi)距離的比值，n為選擇的敏感特征個(gè)數(shù)。3. 根據(jù)權(quán)利要求1所示的基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法，其特征在于，步驟4中利用核密度估計(jì)方法得到每類故障樣本集的概率密度函數(shù)，在任意點(diǎn)x處的總體密度函數(shù) 的核密度估計(jì)定義為：（2）其中， 稱為核函數(shù)(kernel function)，h為窗寬。相比傳統(tǒng)的分類器算法，這種算法不但能夠考慮到每個(gè)樣本所帶來的個(gè)體信息，更重要的是，它很好的考慮到了各類樣本集中樣本的聚集情況；換句話說，它很好的考慮到了樣本間的關(guān)聯(lián)統(tǒng)計(jì)信息，這對于處理樣本分布的隨機(jī)性是至

4、關(guān)重要的，因此可以更好的實(shí)現(xiàn)正確分類。4. 根據(jù)權(quán)利要求1所示的基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法，其特征在于，步驟5中使用對稱化的K-L散度來量化兩個(gè)分布的相似性。原始的K-L散度定義不滿足距離的對稱性，于是在本發(fā)明中使用了如下的對稱化K-L散度：（3）使用該K-L散度來量化分布相似性，可以方便地實(shí)現(xiàn)算法的程序化，提高算法的效率。除此之外，在步驟6中定義了集成K-L散度指標(biāo)，從而將所有敏感特征描述下的散度值歸一化，方便分類判斷，這樣做既考慮了全面性，又考慮了可行性。說明書一種基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法技術(shù)領(lǐng)域本發(fā)明屬于機(jī)械裝備故障智能診斷領(lǐng)域，具體涉及基于

5、核密度估計(jì)和K-L 散度等統(tǒng)計(jì)學(xué)工具的故障診斷方法。背景技術(shù)旋轉(zhuǎn)機(jī)械廣泛應(yīng)用于工業(yè)生產(chǎn)實(shí)踐中，例如風(fēng)力發(fā)電機(jī)、數(shù)控機(jī)床、航空航天發(fā)動機(jī)等關(guān)系到國防民生的重要領(lǐng)域。在生產(chǎn)工作中，滾動軸承、齒輪等旋轉(zhuǎn)機(jī)械中的關(guān)鍵部件由于需要承受交變機(jī)械應(yīng)力和偶然沖擊，加上本身固有的制造誤差，經(jīng)常會產(chǎn)生一些早期缺陷，例如輕度磨損、點(diǎn)蝕等。這些缺陷如果不及時(shí)診斷發(fā)現(xiàn)，就會不斷惡化，最終導(dǎo)致系統(tǒng)失效，帶來很大的財(cái)產(chǎn)損失，甚至對國防和人身安全帶來巨大威脅。對于現(xiàn)代化大型復(fù)雜旋轉(zhuǎn)機(jī)械設(shè)備而言，盡管可以通過改善設(shè)計(jì)、制造工藝來提高零部件的質(zhì)量，但仍難以確保不出故障。所以，有必要利用先進(jìn)的傳感和監(jiān)測技術(shù)，對關(guān)鍵零部件及系統(tǒng)進(jìn)行

6、有計(jì)劃、有組織、有針對性的狀態(tài)監(jiān)測和故障診斷，盡早發(fā)現(xiàn)設(shè)備運(yùn)行過程中的各種隱患，從而防止巨額財(cái)產(chǎn)損失和災(zāi)難性事故的發(fā)生。數(shù)據(jù)驅(qū)動方法是近年來逐漸興起的一種故障診斷技術(shù)，計(jì)算機(jī)技術(shù)的快速發(fā)展使得大數(shù)據(jù)并行高速計(jì)算變得非常容易，推動了依靠大量數(shù)據(jù)分析的故障診斷技術(shù)的發(fā)展。從應(yīng)用的角度看，基于數(shù)據(jù)驅(qū)動的故障診斷方法較基于模型的方法更為切實(shí)可行，這是由于數(shù)據(jù)采集通常要比精確建立物理模型更加容易。除此之外，數(shù)據(jù)驅(qū)動的故障診斷方法還有兩個(gè)明顯的優(yōu)點(diǎn)：一是該類方法更容易實(shí)現(xiàn)自動診斷，這與現(xiàn)代工業(yè)的智能化發(fā)展是切合的；二是該類方法不需要太多參數(shù)設(shè)置和專家經(jīng)驗(yàn)知識。一般來說，一種數(shù)據(jù)驅(qū)動的故障診斷方法應(yīng)包括數(shù)據(jù)

7、獲取、特征提取、特征降維、分類器設(shè)計(jì)和結(jié)果輸出等五個(gè)步驟，其中分類器設(shè)計(jì)和選擇是該類方法的關(guān)鍵?，F(xiàn)有的數(shù)據(jù)驅(qū)動的故障診斷方法大多都是在樣本的特征空間尋找一個(gè)最優(yōu)分類超曲面，從而將不同類型的故障樣本分開。例如，基于支持向量機(jī)（Support Vector Machine，SVM）的故障診斷方法，基于BP（Back Propagation）神經(jīng)網(wǎng)絡(luò)的故障診斷方法等。然而，由于噪聲、測量誤差等對有效振動信號的污染，使得分類問題存在一定的不確定性，進(jìn)而導(dǎo)致了分類錯誤的現(xiàn)象，即很難找到一個(gè)恰當(dāng)?shù)某鎸⑺袠颖救繗w類正確。以上常用方法產(chǎn)生分類錯誤的現(xiàn)象的原因與它們的分類原理是分不開的。傳統(tǒng)的智能故障診

8、斷方法往往忽略了樣本間的統(tǒng)計(jì)信息和關(guān)聯(lián)信息，而統(tǒng)計(jì)信息對于隨機(jī)信號處理是極為關(guān)鍵的，即對于正確分類很有幫助。目前，國內(nèi)外從樣本統(tǒng)計(jì)角度開展的智能故障診斷方法的研究或報(bào)道還非常少。發(fā)明內(nèi)容本發(fā)明的目的是為了從原始樣本中提取更加全面、有效的統(tǒng)計(jì)學(xué)信息，從而提高分類器的準(zhǔn)確率和推廣能力，提出一種基于核密度估計(jì)（KDE，Kernel Density Estimation）和K-L散度（Kullback-Leibler Divergence）兩種統(tǒng)計(jì)學(xué)工具的智能故障診斷方法。本發(fā)明的技術(shù)方案是：一種基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法，具體包括如下步驟：步驟1：采集被監(jiān)測對象的原始振動數(shù)據(jù)，

9、并劃定訓(xùn)練樣本集和測試樣本集；步驟2： 利用集合經(jīng)驗(yàn)?zāi)B(tài)分解方法（Ensemble Empirical Mode Decomposition，EEMD）和希爾伯特變換等信號處理方法，從原始振動數(shù)據(jù)中提取指定時(shí)頻域特征；步驟3： 利用基于距離的特征評估方法，從步驟2中得到的特征集中選擇出少數(shù)敏感特征，并且計(jì)算這些敏感特征的分類貢獻(xiàn)率；改進(jìn)的基于距離的特征評估方法包括以下四步：（1）計(jì)算第j個(gè)特征C個(gè)類的類內(nèi)距離的平均值：（1）（2）式中，Mc表示第c類的樣本個(gè)數(shù)，J表示特征個(gè)數(shù)，C表示類別個(gè)數(shù)，qm,c,j表示第c類第m個(gè)樣本的第j個(gè)特征的特征值。（2）計(jì)算第j個(gè)特征C個(gè)類的類間距離的平均值：（

10、3）（4）式中，uc,j、ue,j分別表示第c和第e個(gè)類第j個(gè)特征的平均值。（3）數(shù)據(jù)集A的類間距離與類內(nèi)距離的比值：（5）越大表示第j個(gè)特征對于類別越敏感，更符合評估原則，也是更應(yīng)該被選擇使用的特征。 （4）計(jì)算選擇的前n個(gè)特征的分類貢獻(xiàn)率：對于第j個(gè)特征，定義：， （6）為其分類貢獻(xiàn)率。當(dāng)完成特征選擇之后，得到的敏感特征將作為后續(xù)分類器的輸入，基于核密度估計(jì)和K-L散度的分類器是本發(fā)明的核心內(nèi)容，包含步驟4到步驟6。步驟4： 利用核密度估計(jì)計(jì)算訓(xùn)練樣本中不同故障類別樣本集關(guān)于步驟3中提取的敏感特征的概率密度函數(shù)，并計(jì)算加入一個(gè)未知故障類別待測樣本后各類樣本集新的概率密度函數(shù)；核密度估計(jì)是概

11、率論中一種用來估計(jì)隨機(jī)變量概率密度函數(shù)的非參數(shù)方法。設(shè)是取自一元連續(xù)總體X的樣本，在任意點(diǎn)x處的總體密度函數(shù) 的核密度估計(jì)定義為：（7）其中， 稱為核函數(shù)（Kernel Function），h為窗寬。選用高斯函數(shù)作為核密度估計(jì)方法的核函數(shù)。h是核密度估計(jì)方法中唯一要優(yōu)化的參數(shù)，根據(jù)最小化平均集成平方誤差方法，最優(yōu)化h的取值由（8）式得到：（8）式中，是樣本標(biāo)準(zhǔn)差，n為各類樣本集所包含的樣本數(shù)。步驟5：求出在選定的特征描述下，訓(xùn)練樣本中各類故障樣本集原始概率密度函數(shù)，以及加入一個(gè)待測樣本后的新概率密度函數(shù)兩者的K-L散度值；K-L散度也稱為相對熵或者信息增益。在概率論和信息論中，K-L散度用來計(jì)

12、算兩個(gè)分布的分對稱性或者差異性。K-L散度值越小，表明兩個(gè)分布越相似。原始的K-L散度定義式為：（9）由于其不滿足距離的對稱性，于是在本發(fā)明中定義了下式的對稱化K-L散度：（10）由步驟4得到各類樣本集的原始密度函數(shù)和加入測試樣本后的新密度函數(shù)之后，就可以根據(jù)（10）式計(jì)算出對應(yīng)的K-L散度值，即量化原始密度函數(shù)和新密度函數(shù)的差異性。 步驟6：計(jì)算集成K-L散度，并通過集成K-L散度的大小判斷待測樣本的故障類別。由步驟5可知，我們計(jì)算得到了前n個(gè)特征描述下的K-L散度值，為了便于判斷，在這里，我們定義了一個(gè)集成K-L散度值（）：（11）式中，為步驟3得到的前n個(gè)特征分類貢獻(xiàn)率。（11）式中通過

13、將第j個(gè)特征描述下的K-L散度值與第j個(gè)特征的分類貢獻(xiàn)率加權(quán)，可以得到各類樣本集的集成K-L散度值。越小，原始分布和加入待測樣本后的新概率密度函數(shù)越相似；反之，當(dāng)集比較大時(shí)，原始概率密度函數(shù)與新概率密度函數(shù)差別越大。換句話說，待測樣本應(yīng)該歸類于K-L散度值最小的類別。附圖說明圖1為本發(fā)明的基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械故障診斷方法的流程示意圖圖2為本發(fā)明實(shí)例所針對的來自凱斯西儲大學(xué)的滾動軸承故障診斷平臺圖3數(shù)據(jù)集A所有特征類間/類內(nèi)距離比值圖4基于核密度估計(jì)和K-L散度的分類器原理示意圖圖5將本發(fā)明用于表1中數(shù)據(jù)集A的實(shí)驗(yàn)結(jié)果：四類樣本集原始概率密度曲線及對應(yīng)的加入一個(gè)正常樣本后新的概率

14、密度函數(shù)圖6為不同特征數(shù)量下，本發(fā)明方法在數(shù)據(jù)集A-E上的分類準(zhǔn)確率圖7為不同訓(xùn)練樣本下，本發(fā)明方法在數(shù)據(jù)集A-E上的分類準(zhǔn)確率具體實(shí)施方式下面結(jié)合附圖和具體的實(shí)例對本發(fā)明做進(jìn)一步的闡述。附圖1所示為基于核密度估計(jì)和K-L散度的旋轉(zhuǎn)機(jī)械智能故障診斷方法的流程圖，如圖1所示，其步驟包括：步驟1：采集被監(jiān)測對象的原始振動數(shù)據(jù)，并劃定訓(xùn)練樣本集和測試樣本集。本發(fā)明以凱斯西儲大學(xué)的滾動軸承故障診斷平臺為例，具體的實(shí)驗(yàn)參數(shù)如下：如附圖2所示，實(shí)驗(yàn)平臺包括一個(gè)2馬力的電機(jī)（左側(cè)），一個(gè)功率計(jì)（右側(cè)）和電子控制設(shè)備（沒有顯示）。被測軸承起到支撐電機(jī)軸的作用。利用電火花加工技術(shù)在驅(qū)動端軸承上人為制造了單點(diǎn)故障

15、，故障直徑分別為0.007、0.014、0.021和0.028英寸。這些故障單獨(dú)分布在軸承內(nèi)圈、滾動體和外圈上。電機(jī)轉(zhuǎn)軸上加載了一個(gè)沖擊力，使用兩個(gè)加速度傳感器測量振動，一個(gè)傳感器布置在電機(jī)殼體上，另一個(gè)傳感器布置在驅(qū)動端軸承外圈上。采樣頻率為12KHz，單采樣樣本長度為12000。為了驗(yàn)證上文所提方法的有效性，我們進(jìn)行了兩種實(shí)驗(yàn)，第一種是不同故障類型判別實(shí)驗(yàn)，第二種是同故障損傷程度判別實(shí)驗(yàn)，詳細(xì)信息如表1和表2所示。數(shù)據(jù)集A包括軸承4種不同的運(yùn)行狀態(tài)（正常、外圈故障、內(nèi)圈故障、滾動體故障）的280個(gè)數(shù)據(jù)樣本，損傷尺寸為0.007 英寸。數(shù)據(jù)集A分為兩部分，140個(gè)樣本用來訓(xùn)練，140個(gè)用來測

16、試，很明顯這是一個(gè)四種不同健康狀態(tài)的模式識別問題。數(shù)據(jù)集B同樣包含280個(gè)數(shù)據(jù)樣本，其中訓(xùn)練樣本與測試樣本仍然各占一半，包含內(nèi)圈和滾動體兩種故障類型。數(shù)據(jù)集B又包含B1和 B2兩個(gè)數(shù)據(jù)子集，每個(gè)包含140個(gè)樣本。數(shù)據(jù)子集B1包含70個(gè)訓(xùn)練樣本，每個(gè)樣本的損傷程度為0.007英寸，同時(shí)包含損傷程度為0.021英寸的70個(gè)測試樣本。數(shù)據(jù)子集B2與B1類似，只是將B1的訓(xùn)練樣本和測試樣本進(jìn)行交換，即B1的測試樣本作為B2的訓(xùn)練樣本。使用數(shù)據(jù)集B的作用是驗(yàn)證該方法對于同種故障下、不同于訓(xùn)練樣本損傷程度的測試樣本的檢測能力，即驗(yàn)證該方法的穩(wěn)健性。數(shù)據(jù)集C、D和E各包含210個(gè)樣本，訓(xùn)練和測試樣本各占一半

17、。數(shù)據(jù)集C、D和E的故障類型分別為內(nèi)圈故障、滾動體故障和外圈故障。每種故障損傷程度分為三個(gè)層次，0.007、0.021和0.028英寸。使用數(shù)據(jù)C、D和E的目的是驗(yàn)證該方法對于同種故障不同損傷程度的分辨能力。表1 故障類型判別實(shí)驗(yàn)數(shù)據(jù)集數(shù)據(jù)集訓(xùn)練樣本數(shù)測試樣本數(shù)損傷尺寸（英寸）狀態(tài)A35350/0正常35350.007/0.007內(nèi)圈故障35350.007/0.007滾動體故障35350.007/0.007外圈故障BB135350.007/0.021內(nèi)圈故障35350.007/0.021滾動體故障B235350.021/0.007內(nèi)圈故障35350.021/0.007滾動體故障表2 故障損傷程

18、度判別實(shí)驗(yàn)數(shù)據(jù)集數(shù)據(jù)集訓(xùn)練樣本數(shù)測試樣本數(shù)損傷程度（英寸）狀態(tài)C35350.007內(nèi)圈故障35350.01435350.021D35350.007滾動體故障35350.01435350.021E35350.007外圈故障35350.01435350.021為了精簡篇幅，這里僅以數(shù)據(jù)集A為例敘述本發(fā)明方法的后續(xù)步驟。步驟2： 利用EEMD和希爾伯特變換等信號處理方法從原始振動數(shù)據(jù)中提取指定時(shí)頻域特征。利用EEMD對數(shù)據(jù)集A所對應(yīng)的振動信號樣本進(jìn)行分解，白噪聲初始值為0.3，集成數(shù)為100，于是可以得到每個(gè)樣本的前四階IMF（固有模態(tài)分量）。接下來，分別從數(shù)據(jù)集A原始振動信號和四階IMF分量中提取

19、制定的9個(gè)時(shí)域特征和10個(gè)頻域特征，從而得到每類樣本集所對應(yīng)的95個(gè)原始特征。所提取的時(shí)域特征和頻域特征種類如表3和表4所示。表3 時(shí)域特征表()振動信號的離散時(shí)間序列；是數(shù)據(jù)點(diǎn)數(shù)。表4 頻域特征表是時(shí)間序列的頻譜， ()； 是頻譜譜線條數(shù)； 是第k條譜線。步驟3：利用基于距離的特征評估方法從步驟2中得到的特征集中選擇出前10個(gè)敏感特征，并且計(jì)算這些敏感特征的分類貢獻(xiàn)率，以供后續(xù)步驟使用。應(yīng)用基于距離的特征評估方法包括以下四步：（1）計(jì)算數(shù)據(jù)集A中四類樣本集在第j個(gè)特征描述下的類內(nèi)距離的平均值：數(shù)據(jù)集A中共包含四類故障，各類樣本數(shù)為35，總的特征數(shù)為95，所以Mc=35，C=4，J=95；將這

20、些參數(shù)帶入（1）和（2）式，有：式中，qm,c,j表示第c類第m個(gè)樣本的第j個(gè)特征的特征值。（2）計(jì)算數(shù)據(jù)集A中四類樣本集在第j個(gè)特征描述下的類間距離的平均值：同上步驟，將相關(guān)參數(shù)帶入（3）和（4）式，得：式中，uc,j、ue,j分別表示第c和第e個(gè)類第j個(gè)特征的平均值。（3）數(shù)據(jù)集A的類間距離與類內(nèi)距離的比值：，95個(gè)特征的值如附圖3所示。（4）在本實(shí)驗(yàn)中，我們選取了值排名前10的特征作為敏感特征，并計(jì)算得到了這10個(gè)特征的分類貢獻(xiàn)率。這10個(gè)特征的編號和值如表5所示，進(jìn)一步，計(jì)算選擇的這10個(gè)敏感特征的分類貢獻(xiàn)率。對于第j個(gè)特征，有：， 表5 敏感特征的類間、類內(nèi)距離比值特征編號15142

21、0221316191056 1.000.950.860.860.850.850.790.760.690.66 附圖4中以兩種故障為例，介紹了本發(fā)明基于核密度估計(jì)和K-L散度的分類器原理，應(yīng)用于此實(shí)例，具體如以下步驟4到6所述。步驟4：利用核密度估計(jì)計(jì)算訓(xùn)練樣本中不同故障類別樣本集關(guān)于步驟3中提取的敏感特征的概率密度函數(shù)和加入一個(gè)未知故障類別待測樣本后各類樣本集新的概率密度函數(shù)。針對數(shù)據(jù)集A，各類故障樣本集訓(xùn)練樣本數(shù)是35，帶入（8）式有：則可以計(jì)算出核密度估計(jì)的最佳帶寬h。 基于上述方法，可得到四類樣本集的關(guān)于第j特征的原始核密度函數(shù)（）以及加入一個(gè)相同待測樣本后四類樣本集新的核密度函數(shù)（），示例如附圖5所示。步驟5： 求出在選定的特征描述下，訓(xùn)練樣本中各類故障樣本集原始概率密度函數(shù)，以及加入一個(gè)待測樣本后的新概率密度函數(shù)兩者的K-L散度值。由步驟4得到數(shù)據(jù)集A中四類樣本集的原始密度函數(shù)和加入測試樣本后的新密度函數(shù)之后，根據(jù)對稱化K-L散度定義，（10）求得在步驟3中選擇的10個(gè)敏感特征描述下，數(shù)據(jù)集A中各類樣本集原始密度函數(shù)與新密度函數(shù)的K-L散度值。步驟6：計(jì)算集成K-L散度，并通過集成K-L散度的大小判斷待測樣本的故障類別。由步驟5可知，我們計(jì)算得到了前10個(gè)特征描述下的K-L散度值，根據(jù)式（11）中集成K-L散度值的定義，有：，式中