基于多尺度基本熵的CFS聚類滾動軸承故障診斷

海波，姚海龍，孟叢叢

（蘭州城市學(xué)院 電子與信息工程學(xué)院，蘭州 730070）

使用軸承振動信號進(jìn)行故障診斷是常用的故障診斷方式，振動信號的特征提取和故障識別則是滾動軸承故障診斷的重要步驟［1-3］。

在特征提取方面，近似熵（Approximate Entropy，AE）與樣本熵（Sample Entropy，SE）對信號數(shù)據(jù)的長度較為敏感，且在熵值的計算過程中需對信號進(jìn)行2次重構(gòu)，算法較為復(fù)雜［4-7］；排列熵（Permutation Entropy，PE）雖然僅需對信號進(jìn)行1次重構(gòu)，但計算過程中仍需要復(fù)雜的排序運算［8-10］；基本尺度熵（Base-Scale Entropy，BSE）計算簡單、效率高，且具有一定的抗干擾能力，在齒輪故障診斷中有一定的應(yīng)用［11-13］；然而，上述方法只能檢測時間序列信號在單一尺度上的復(fù)雜性，而多尺度熵（Multiscale Entropy，ME）能衡量時間序列信號在不同尺度下的復(fù)雜性［14-15］，其中多尺度排列熵（Multiscale Permutation Entropy，MPE）模型在滾動軸承信號特征提取中的效果較好［16-17］。

在故障識別方面，GK（Gustafson-Kessel）和GG（Gath-Geva）聚類算法比較常用，但均需要事先劃分聚類數(shù)目［18-22］；CFS（Clustering by Fast Search）聚類算法［23］簡單明快，其根據(jù)數(shù)據(jù)點的局部密度和截斷距離對具有跳躍性的數(shù)據(jù)點進(jìn)行聚類中心點的選定，不需要聚類數(shù)目的劃分。

綜上所述，雖然ME可獲得多尺度的信息，但卻存在數(shù)據(jù)維數(shù)高，無法可視化等問題，可采用主成分分析法（Principal Component Analysis，PCA）進(jìn)行數(shù)據(jù)降維，從而作為CFS的輸入。因此，將多尺度基本熵（Multiscale Base-scale Entropy，MBSE）與CFS聚類算法結(jié)合，嘗試用于滾動軸承的故障診斷。

1 多尺度基本熵

1.1 基本尺度熵

BSE的計算步驟如下：

1）對于一個給定數(shù)據(jù)長度的時間序列u()i（1≤i≤N），數(shù)據(jù)集 Xmi（1≤i≤N-m+1）需進(jìn)行如下重構(gòu)

3）對得到的BS值進(jìn)行符號劃分，將轉(zhuǎn)變?yōu)榉栃蛄?Si（X（i））＝｛s（i），…，s（i+m-1）｝，s∈ A（1，2，3，4），轉(zhuǎn)換方式為

4）統(tǒng)計符號序列Si的概率分布情況P（π），由于Si（X（i））由m維時間序列Xmi形成，且存在4種符號，所以共有4m種不同組合狀態(tài)π，則統(tǒng)計概率 P（π）為

5）計算BSE，其計算式為

1.2 多尺度基本熵

BSE確定的是時間序列在單一尺度上的復(fù)雜度和無規(guī)則程度，MBSE反映的則是時間序列在不同尺度下的復(fù)雜性程度，其定義為時間序列在不同尺度下的BSE，其中多尺度通過粗粒化過程獲得，具體的計算步驟如下：

1）對于原始數(shù)據(jù)Xmi，建立新的粗粒向量為

式中：τ為尺度因子（τ＝1，2，…），若 τ＝1則為原始數(shù)據(jù)；M為樣本總數(shù)。對于非零原始序列Xmi，其被分為τ個每段長M／τ的粗粒序列yj（τ）。

2）分別計算τ個粗粒序列為yj（τ）的 BSE，形成MBSE。

2 主成分分析

MBSE模型提取的信號維度較高，存在數(shù)據(jù)冗余、計算效率低、無法可視化等問題，而PCA可以抽取數(shù)據(jù)的主要成分，對特征向量進(jìn)行降維以作為CFS聚類算法的輸入。

對于一個給定的MBSE或MPE矩陣X＝［x1，x2，…，xN］，其中N為樣本數(shù)量，n為樣本的維數(shù)，且滿足n＜N。其計算步驟為：

3）使用奇異值方法分解協(xié)方差矩陣，則Rx＝ATΛA，ATA＝AAT＝I。A＝［A1，A2，…，An］T為矩陣X對應(yīng)的特征值矩陣，其中Λ＝diag（λ1，λ2，…，λn），是矩陣A按照降序排列的對角矩陣。

4）使用累積貢獻(xiàn)率θk確定X的主要成分個數(shù)，即

式中：θk為從第1到第k個主要成分的百分比。

3 CFS聚類算法

CFS聚類算法的核心思想在于對聚類中心的刻畫上，聚類中心點的特點為：1）自身密度大，即被密度均不超過自己的近鄰數(shù)據(jù)點包圍；2）與其他密度更大數(shù)據(jù)點之間的距離更大。

CFS聚類的詳細(xì)計算過程［23］如下：

1）設(shè)具有N個數(shù)據(jù)點的數(shù)據(jù)集合X＝｛x1，x2，…，xN｝，該集合中的任何數(shù)據(jù)點xi與xj之間的距離為

式中：當(dāng)xi具有最大局部密度時，δi表示數(shù)據(jù)集合X中與數(shù)據(jù)點xi距離最大的數(shù)據(jù)點與xi之間的距離；否則δi表示在所有局部密度大于xi的數(shù)據(jù)點中與xi距離最小數(shù)據(jù)點之間的距離。

4）計算每個數(shù)據(jù)點的γ值并進(jìn)行降序排列

5）根據(jù)γ確定其聚類中心，γi越大表示成為聚類中心點的可能性越大。故對γ進(jìn)行降序排列，依次選取γ較大且具有跳躍性的數(shù)據(jù)點作為聚類中心點，即聚類中心點過渡到非聚類中心點時γ值具有明顯的跳躍性。

6）統(tǒng)計每個數(shù)據(jù)點與聚類中心的距離小于截斷距離dc的個數(shù)（即ρi），從而實現(xiàn)分類。

4 試驗分析

4.1 數(shù)據(jù)來源

試驗軸承為6205-2RS JEM-SKF深溝球軸承，使用電火花技術(shù)分別在軸承內(nèi)圈、外圈、球上設(shè)置單點故障［24］。本試驗選用的故障直徑為0.177 8 mm，電動機功率為1 496W，轉(zhuǎn)速為1 750 r／min，采樣頻率為12 kHz。在此情況下采集到正常（NR），內(nèi)圈故障（IRF），外圈故障（ORF）和球故障（BF）4種狀態(tài)的振動信號。

4.2 參數(shù)設(shè)置

BSE：嵌入維數(shù) m的取值范圍一般為3～7［11-13］，且要滿足 4m≤ N條件（N為序列的長度）；參數(shù) a一般取值為0.1～0.4［11-13］；根據(jù)軸承振動信號的特性，選擇m＝4，a取值為0.2或0.3。

MBSE：尺度因子 τ一般設(shè)置為 20［15-17］，振動信號的長度統(tǒng)一取為N＝2 048。

CFS：需選取一個階段距離參數(shù)dc，使得每個數(shù)據(jù)點的平均近鄰數(shù)據(jù)點（即距離不超過dc的數(shù)據(jù)點）約占數(shù)據(jù)點總數(shù)的1%～2%。dc設(shè)置太大，會使每個數(shù)據(jù)點xi的局部密度數(shù)值ρi較大，導(dǎo)致區(qū)分度不高；若dc設(shè)置過小，則會導(dǎo)致同一個簇類被拆分成多個［23］。綜合考慮，設(shè)置dc為比例k／M＝1.5%情況下的數(shù)值。

4.3 試驗流程

試驗流程如圖1所示，主要步驟為：

圖1 試驗流程Fig.1 Experimental flow chart

1）預(yù)先設(shè)置模型的參數(shù)；

2）使用MBSE模型對采集的振動信號進(jìn)行多尺度分解，得到不同尺度因子下的MBSE值；

3）統(tǒng)計MBSE的計算總時間和平均時間，進(jìn)行計算效率對比；

4）使用PCA方法進(jìn)行數(shù)據(jù)降維，選擇（9）式中的參數(shù)k＝2時的主要成分作為CFS聚類模型的輸入；

5）與MPE模型進(jìn)行對比分析。

4.4 對比分析

對采集到4種狀態(tài)的實際振動信號進(jìn)行編號：NR為 1～50，IRF為51～100，BF為 101～150，ORF為151～200。每種狀態(tài)下振動信號某一樣本的幅值如圖2所示。從圖中可以看出：NR和BF的信號沒有明顯的周期性，這是由于NR和BF振動信號的隨機性較強，故兩者的自相似性較低，兩者難以區(qū)分；IRF和ORF信號具有較為明顯的振動規(guī)律，這是由于兩者的自相似性較高，特別是內(nèi)圈固定、外圈轉(zhuǎn)動時的振動規(guī)律更加明顯，但這2種信號的振動規(guī)律相似，也難以區(qū)分。

圖2 軸承不同故障狀態(tài)下的振動信號Fig.2 Vibration signals of rolling bearing under different faults

為準(zhǔn)確進(jìn)行故障識別，需進(jìn)一步提取信號的特征向量，根據(jù)預(yù)先設(shè)置的參數(shù)，使用MPE／MBSE模型對實際振動信號進(jìn)行多尺度分解，得到的不同尺度因子下的熵值，如圖3所示。

圖3 多尺度分解結(jié)果Fig.3 Multi-scale decomposition results

從圖3中可以看出，MBSE與MPE值整體上呈現(xiàn)減小的趨勢，因為振動信號在分解之后，頻率依此從高到低，故呈現(xiàn)減小趨勢，表明其信息主要包含在尺度因子較小的MBSE與MPE值中。經(jīng)過統(tǒng)計得到：MBSE模型計算用時22.219 8 s，MPE模型計算用時38.246 7 s，即MBSE模型的計算效率較高。主要原因是由于PE模型在統(tǒng)計每種符號序列出現(xiàn)的概率前需要對N-m+1個m維矢量進(jìn)行升序排序［9-10］，但 BSE模型在統(tǒng)計概率P(π)前只需要計算相鄰點數(shù)間隔的差值方均根值，不需要進(jìn)行復(fù)雜的排序運算［11-13］；另外，在劃分到多尺度上時，需要計算尺度從1～20情況下的PE／BSE值，會加大MPE／MBSE模型的計算時間落差。

由于上述熵值特征矩陣的維度為20，在進(jìn)行故障識別時會造成數(shù)據(jù)不易可視化和信息冗余等問題。故使用PCA模型進(jìn)行降維，根據(jù)（9）式計算所得的累計貢獻(xiàn)率結(jié)果見表1，表中λ1與λ2代表原始數(shù)據(jù)經(jīng)過PCA降維后第1～2個特征數(shù)值。由表可知，前2個特征數(shù)值之和約占總特征數(shù)值的80%左右，這是由于特征數(shù)值在計算過程中按降序排列，故λ1與λ2基本包含了信號的大部分信息，可選作故障特征向量。

表1 累計貢獻(xiàn)率θ隨k的數(shù)值變化表Tab.1 The value of θ changes with k

隨后使用CFS模型對故障特征向量進(jìn)行聚類，依據(jù)（11）～（13）式計算得到200個樣本的的聚類中心點的局部密度ρ、距離δ以及γ，結(jié)果如圖4和圖5所示。

圖4 局部密度ρ與距離δFig.4 Local density ρ and distance δ

圖5 樣本點的γ值Fig.5 Sample point of γ

從圖中可以看出：1）聚類中心點的局部密度ρ與距離δ的數(shù)值明顯大于非聚類中心點；2）聚類中心點的γ整體上明顯大于非聚類中心點。

其中，γ最小為0.818 4，最大為8.851 4，而非聚類中心點的γ數(shù)值幾乎為0，其原因如下：

1）對于聚類中心點而言，利用（10）式計算任意樣本之間的距離時，聚類中心點之間的距離最大，而（11）式為ρ的減函數(shù)，即計算所有樣本的局部密度ρ時，聚類中心點之間的距離是最大。設(shè)置的階段距離為1.5%，因此（11）式中的 -（dij／dc）2較小，局部密度數(shù)值ρ偏大，非聚類中心點的局部密度ρ與距離δ數(shù)值接近于0，與圖4的結(jié)果相符。

2）利用（13）式計算γ時，由于聚類中心點的局部密度ρ與距離δ遠(yuǎn)大于非聚類中心點，故得到的γ也大于非聚類中心點。而且圖中的聚類中心點具有明顯的跳躍性，而非聚類中心點則無此特點，說明應(yīng)該選擇γ數(shù)值較大且有跳躍性這些數(shù)據(jù)點作為聚類中心點。最后的2維聚類結(jié)果如圖6所示，從圖中可以看出，4種故障狀態(tài)得到了較好的區(qū)分，各樣本點之間幾乎沒有重疊現(xiàn)象，聚類中心點被周圍樣本點圍繞的密度較大，中心點相互之間的距離也較大。且樣本點與其他中心點之間的距離比與自身中心點之間的距離大，統(tǒng)計數(shù)據(jù)也表明，CFS模型對試驗信號的區(qū)分程度近似100%，具有較好的區(qū)分效果。

圖6 2維聚類結(jié)果Fig.6 The 2D clustering results

5 結(jié)束語

提出了一種MBSE與CFS相結(jié)合的故障診斷方法，并將其應(yīng)用于滾動軸承的故障診斷。與MPE模型相比，使用MBSE模型對軸承振動信號進(jìn)行特征提取的計算效率更高，而通過選用不同故障尺寸的軸承樣本數(shù)據(jù)進(jìn)行試驗，表明經(jīng)過PCA降維處理后，CFS聚類算法可對不同故障類型進(jìn)行有效的區(qū)分。

但也存在不足之處，如CFS聚類算法在選擇聚類中心點時選擇局部密度ρ與距離δ的乘積作為標(biāo)準(zhǔn)，雖然綜合考慮了這2個參數(shù)，但兩者的數(shù)值可能處于不同數(shù)量級，是否能通過其他方式綜合考慮這2個參數(shù)，更好的選擇聚類中心抉擇則需要進(jìn)行更深入的研究。