王道嶸，梁濤，王建輝，姜文

(1.河北工業(yè)大學(xué) 人工智能與數(shù)據(jù)科學(xué)學(xué)院，天津 300130；2.河北建投能源投資股份有限公司，石家莊 050011)

軸承作為一種機(jī)械標(biāo)準(zhǔn)件，廣泛應(yīng)用于各種旋轉(zhuǎn)機(jī)械設(shè)備，軸承健康度的預(yù)測對于各種旋轉(zhuǎn)機(jī)械的運(yùn)維具有重要的意義[1]。軸承健康度預(yù)測技術(shù)可以分為基于物理模型的方法與數(shù)據(jù)驅(qū)動(dòng)的方法兩類[2]：由于設(shè)備的結(jié)構(gòu)、機(jī)理、運(yùn)行環(huán)境等錯(cuò)綜復(fù)雜，基于物理模型的方法難以建立起適用性強(qiáng)的模型；數(shù)據(jù)驅(qū)動(dòng)的方法通過分析、處理軸承振動(dòng)信號數(shù)據(jù)建立振動(dòng)信號與軸承健康度的映射關(guān)系模型，從而預(yù)測軸承的健康度，是目前的主流方法。隨著機(jī)器學(xué)習(xí)的發(fā)展，更使得基于數(shù)據(jù)驅(qū)動(dòng)的健康度預(yù)測成為國內(nèi)外研究的重點(diǎn)方向[3]。

從軸承振動(dòng)信號中提取出退化特征信息并將其構(gòu)建為用于描述軸承健康度的指標(biāo)，稱為健康指標(biāo)(Health Indicator，HI)，其構(gòu)建的優(yōu)劣將直接影響模型預(yù)測的準(zhǔn)確性[4]。很多學(xué)者對健康指標(biāo)的構(gòu)建進(jìn)行了大量的研究。文獻(xiàn)[5]使用dispersion-entropy描述設(shè)備的退化過程；文獻(xiàn)[6]用變分自編碼器對軸承原始信號進(jìn)行處理，提取并學(xué)習(xí)軸承的故障特征，用于軸承的故障診斷；文獻(xiàn)[7]利用WP-EMD 提取軸承振動(dòng)信號特征，利用自組織競爭網(wǎng)絡(luò)進(jìn)行融合并生成健康指標(biāo)；文獻(xiàn)[8]計(jì)算了多個(gè)時(shí)、頻域的特征指標(biāo)，選擇與剩余使用壽命相關(guān)性強(qiáng)的幾個(gè)特征輸入循環(huán)神經(jīng)網(wǎng)絡(luò)構(gòu)建健康指標(biāo)；文獻(xiàn)[9]用t-SNE方法處理軸承振動(dòng)信號來提取相關(guān)特征，得到了性能更好的健康指標(biāo)：上述研究采用的特征信息類型大都比較單一，沒有綜合考慮各種類型的特征；而且大多數(shù)學(xué)者都選擇構(gòu)建一個(gè)一維的單調(diào)健康指標(biāo)描述當(dāng)前的健康度，導(dǎo)致所構(gòu)建健康指標(biāo)中包含的健康度信息不夠充分，也導(dǎo)致其提取的特征只能識別出健康度，并不能同時(shí)提取工況等信息。

針對以上問題，本文提出了一種多維健康指標(biāo)的構(gòu)建方法，通過蜉蝣優(yōu)化算法(Mayfly Optimization Algorithm，MOA)優(yōu)化參數(shù)后的變分模態(tài)分解(Variational Mode Decomposition，VMD)對軸承振動(dòng)信號進(jìn)行去噪，針對去噪信號選取合適的特征組成特征向量并降維后得到三維健康指標(biāo)，驗(yàn)證健康指標(biāo)的可靠性后將其送入支持向量機(jī)(Support Vector Machine，SVM)分類器與長短期記憶(Long Short-Term Memory，LSTM)神經(jīng)網(wǎng)絡(luò)進(jìn)行工況識別與健康度預(yù)測。

1 振動(dòng)信號去噪

1.1 變分模態(tài)分解

變分模態(tài)分解能夠自適應(yīng)地實(shí)現(xiàn)信號的頻域剖分及各分量的有效分離，且可以自由確定模態(tài)分解的個(gè)數(shù)[10]， 將變分模態(tài)分解后的模態(tài)分量定義為一個(gè)幅頻調(diào)制信號，分解個(gè)數(shù)為K，則第k個(gè)分量的表達(dá)式為

uk(t)=Ak(t)cos[φk(t)]；k∈{1,…,K} ，

(1)

式中：uk(t)為諧波信號；Ak(t)為uk(t)的瞬時(shí)幅值；φk(t)為相位，其一階導(dǎo)數(shù)表示瞬時(shí)頻率。

變分模態(tài)分解的約束模型公式為

式中：x為原始信號，是各分量uk的累加；ωk為各分量的中心頻率；δ(t)為脈沖函數(shù)。

為得到變分問題的最優(yōu)解，引入增廣的拉格朗日函數(shù)，即

L({uk},{ωk},λ)=

(3)

式中：λ(t)為拉格朗日乘子；α為二次懲罰因子。用交替方向乘子算法迭代更新un+1，wn+1，λn+1，可求得(3)式中的“鞍點(diǎn)”，即(2)式的最優(yōu)解。

1.2 蜉蝣算法優(yōu)化變分模態(tài)分解參數(shù)

在變分模態(tài)分解處理過程中，分解個(gè)數(shù)K和懲罰參數(shù)α均會(huì)極大地影響其分解效果，而其他參數(shù)對分解效果的影響則較小[11]，因此，設(shè)置tau=0,init=1,DC=0,ε=1.0×10-7，分析變分模態(tài)分解前如何選取適當(dāng)?shù)腒和α。

蜉蝣優(yōu)化算法是一種新的仿真優(yōu)化算法，用于解決復(fù)雜的函數(shù)優(yōu)化問題。蜉蝣優(yōu)化算法由雌性蜉蝣群體和雄性蜉蝣群體組成，受蜉蝣動(dòng)物的交配行為啟發(fā)，將雄性蜉蝣的最優(yōu)個(gè)體和雌性蜉蝣的最優(yōu)個(gè)體進(jìn)行交配，得到一個(gè)最優(yōu)子代。同理，雄性次優(yōu)個(gè)體與雌性次優(yōu)個(gè)體進(jìn)行交配得到次優(yōu)的個(gè)體。這一過程與適者生存的規(guī)律一樣，逐步淘汰適應(yīng)度較差的個(gè)體[12]。

最初，隨機(jī)產(chǎn)生2組蜉蝣，分別代表雄性和雌性種群。也就是說，每個(gè)蜉蝣被隨機(jī)放置在問題空間中，作為由d維向量表示的候選解x= (x1,…,xd),并根據(jù)預(yù)先確定的目標(biāo)函數(shù)f(x)對其性能進(jìn)行評價(jià)。蜉蝣的速度v=(v1,…,vd)定義為其位置的變化，每個(gè)蜉蝣的飛行方向是個(gè)體和社會(huì)飛行經(jīng)驗(yàn)的動(dòng)態(tài)交互作用。每個(gè)蜉蝣都會(huì)調(diào)整自己的軌跡，使其朝向目前為止的個(gè)人最佳位置(pbest)以及迄今為止群中任何蜉蝣所獲得的最佳位置(gbest)。雌、雄性蜉蝣具有不同的運(yùn)動(dòng)特點(diǎn)：雄性蜉蝣成群的聚集，意味著每只雄性蜉蝣的位置都是根據(jù)自己和周圍蜉蝣的經(jīng)驗(yàn)來調(diào)整，而雌性蜉蝣則不會(huì)成群聚集，它們會(huì)飛到雄性群體中來繁殖。蜉蝣交配行為則用交叉算子來表示，基于個(gè)體的適應(yīng)度函數(shù)來選擇親本，最好的雌性與最好的雄性繁殖，次好的雌性與次好的雄性繁殖。交叉的結(jié)果是產(chǎn)生2個(gè)后代。其尋優(yōu)步驟見表1。

表1 蜉蝣優(yōu)化算法的尋優(yōu)步驟

為了驗(yàn)證蜉蝣優(yōu)化算法的收斂性和優(yōu)化性能，將(4)式作為適應(yīng)度函數(shù)，選擇遺傳算法(Genetic Algorithm，GA)、粒子群優(yōu)化(Particle Swarm Optimization，PSO)算法進(jìn)行性能對比。

(4)

設(shè)置遺傳算法的迭代次數(shù)為100，種群規(guī)模為100，交叉概率為1，變異概率為0.01；粒子群優(yōu)化算法迭代次數(shù)為100，種群規(guī)模為100；設(shè)置蜉蝣優(yōu)化算法迭代次數(shù)為100，種群數(shù)量為100。3種算法的優(yōu)化迭代曲線如圖1所示：粒子群優(yōu)化算法和遺傳算法都陷入了局部最優(yōu)解，蜉蝣優(yōu)化算法的迭代速度優(yōu)于其他2種算法且其得到的最優(yōu)值更加接近理論最優(yōu)值。說明相對于遺傳算法與粒子群優(yōu)化算法，蜉蝣優(yōu)化算法具有更快的收斂速度和更強(qiáng)的全局搜索能力。

圖1 不同算法的優(yōu)化迭代曲線

用蜉蝣優(yōu)化算法對變分模態(tài)分解的參數(shù)進(jìn)行優(yōu)化時(shí)需要選定合適的適應(yīng)度函數(shù)。軸承振動(dòng)信號經(jīng)過變分模態(tài)分解后各分量的包絡(luò)熵值Ep可反應(yīng)出該分量的稀疏特性[13]。如果分量中包含的噪聲較多，會(huì)掩蓋信號的沖擊特征，則該分量的稀疏性較弱，包絡(luò)熵較大；反之，則包絡(luò)熵較小。選擇每次分解后所有分量中包絡(luò)熵值中最小的一個(gè)作為局部最小包絡(luò)熵minEp，該最小熵值對應(yīng)的分量有著豐富的特征信息。包絡(luò)熵Ep的計(jì)算公式為

(5)

式中：ej為a(j)量綱一化的結(jié)果；a(j)為x(j)經(jīng)希爾伯特解調(diào)得到的包絡(luò)信號。

將局部最小包絡(luò)熵作為參數(shù)優(yōu)化的適應(yīng)度函數(shù)，整個(gè)搜索過程就是要找到全局最小包絡(luò)熵以及對應(yīng)的最佳分量所在的最優(yōu)參數(shù)組合[K,α]。

1.3 相關(guān)性

相關(guān)系數(shù)是反映變量之間相關(guān)關(guān)系密切程度的指標(biāo)。這里采用皮爾遜相關(guān)系數(shù)來分析K個(gè)分量與原始信號的相關(guān)性。設(shè)樣本X和樣本Y，則兩者的相關(guān)系數(shù)為

(6)

式中：r為相關(guān)系數(shù)；Cov(X,Y)為樣本X與樣本Y的協(xié)方差；D(X)，D(Y)分別為樣本X，Y的方差。r越大，代表樣本間的相關(guān)性越高。計(jì)算各個(gè)分量與原始振動(dòng)信號的相關(guān)系數(shù)，挑選最能代表原始信號的部分分量重構(gòu)信號以達(dá)到降噪目的[14]。

2 特征選取

隨著軸承的退化，普通時(shí)域下的簡單特征會(huì)呈現(xiàn)出明顯變化，但這些簡單特征的變化往往呈現(xiàn)為“階段性”，單獨(dú)采用這類特征會(huì)使得退化特征信息的提取效果較差；奇異值分解可以將包含信號特征信息的矩陣分解到不同的子空間中，是一種能夠在擾動(dòng)和噪聲下保持信號特征相對穩(wěn)定的特征提取方法[14]；熵類特征一般用來表征信號中所蘊(yùn)含的各種能量的大小，樣本熵可以反映出時(shí)間序列的復(fù)雜性，序列的復(fù)雜性越高，樣本熵的值就越大[15]，能量熵則會(huì)隨著振動(dòng)信號的能量分布而發(fā)生變化[16]：綜上，選擇重構(gòu)信號的最大值、標(biāo)準(zhǔn)值、峭度以及主要分量的奇異值、樣本熵、能量熵作為提取的多維特征。

t分布隨機(jī)鄰近嵌入(t-SNE)算法是一種深度學(xué)習(xí)的非線性流行學(xué)習(xí)算法，對高維非線性數(shù)據(jù)集有優(yōu)異的降維效果[15]，其核心是引入自由度為1的分布函數(shù)代替高斯分布，構(gòu)建高維空間數(shù)據(jù)樣本的概率分布，并在低維子空間構(gòu)建對應(yīng)樣本的概率分布，采用KL散度(相對熵)衡量高低維空間概率分布的相似程度，從而實(shí)現(xiàn)降維。

3 軸承多維健康指標(biāo)構(gòu)建方法

基于上述理論，軸承多維健康指標(biāo)構(gòu)建方法的流程如圖2所示，具體步驟如下：

1)獲取軸承的振動(dòng)信號。

2)利用蜉蝣優(yōu)化算法優(yōu)化變分模態(tài)分解參數(shù)[K0,α0]，用尋優(yōu)得到的最優(yōu)參數(shù)分解振動(dòng)信號，得到K0個(gè)分量。

3)計(jì)算各分量與原始信號的相關(guān)性，得到相關(guān)性最高的3個(gè)分量進(jìn)行信號重構(gòu)。

4)計(jì)算重構(gòu)信號的最大值、平均值、峭度以及3個(gè)分量的奇異值、樣本熵、能量熵，得到12維特征向量。

5)用t-SNE降維方法將12維特征向量降至3維。

6)將訓(xùn)練集的特征向量輸入用蜉蝣優(yōu)化算法優(yōu)化參數(shù)的SVM分類器和LSTM神經(jīng)網(wǎng)絡(luò)進(jìn)行訓(xùn)練，分別得到工況分類模型與健康度預(yù)測模型。

7)將測試集的特征向量輸入2個(gè)模型，驗(yàn)證模型的可靠性。

圖2 軸承多維健康指標(biāo)構(gòu)建方法的流程圖

4 試驗(yàn)驗(yàn)證

4.1 數(shù)據(jù)描述

選用XJTU-SY數(shù)據(jù)集驗(yàn)證方法的有效性，該數(shù)據(jù)集通過維持固定的徑向力和轉(zhuǎn)速保持高負(fù)載來加速軸承的退化，從而獲取測試軸承全壽命周期的振動(dòng)數(shù)據(jù)，試驗(yàn)軸承型號為LDK UER204，相關(guān)參數(shù)見表2[17]。

表 2 LDK UER204軸承的參數(shù)

設(shè)置采樣頻率為 25.6 kHz，采樣間隔為 1 min。2個(gè)加速度傳感器分別用于測量軸承水平與垂直方向上的信號。忽略出現(xiàn)的混合故障，以及與其他數(shù)據(jù)有明顯差異的數(shù)據(jù)。數(shù)據(jù)共有3類工況，每類工況下有4套軸承。工況數(shù)據(jù)和軸承壽命數(shù)據(jù)見表3。

表3 試驗(yàn)軸承工況及其壽命數(shù)據(jù)

Bearing1_1水平、垂直方向的全壽命周期振動(dòng)信號如圖3所示，對比表明水平方向振動(dòng)信號的振幅明顯大于垂直方向，更有利于對信號進(jìn)行試驗(yàn)研究，故只選取了水平方向信號進(jìn)行后續(xù)分析。

圖3 水平與垂直方向的信號對比

4.2 數(shù)據(jù)處理

軸承原始振動(dòng)信號往往包含著大量噪聲，嚴(yán)重影響特征信息的提取，因此在特征提取之前對信號進(jìn)行一次“去噪”很有必要。每進(jìn)行一次采樣得到32 768個(gè)點(diǎn)，截取前10 000個(gè)樣本點(diǎn)作為試驗(yàn)數(shù)據(jù)，數(shù)據(jù)采樣示意圖如圖4所示。

圖4 數(shù)據(jù)采樣示意圖

利用蜉蝣優(yōu)化算法對變分模態(tài)分解的參數(shù)組合[K0,α0]進(jìn)行尋優(yōu)，尋優(yōu)目標(biāo)是使得分解后各分量的最小包絡(luò)熵Ep值最小。這里以一個(gè)采樣點(diǎn)的數(shù)據(jù)為例(一個(gè)采樣點(diǎn)指的是進(jìn)行一次數(shù)據(jù)采樣得到的10 000個(gè)數(shù)據(jù)點(diǎn))，得到的最優(yōu)參數(shù)組合為α=4 000,K=10，對每個(gè)采樣點(diǎn)的10 000個(gè)數(shù)據(jù)進(jìn)行變分模態(tài)分解，得到10個(gè)分量。

每個(gè)分量與原始信號的相關(guān)系數(shù)如圖5所示，選擇與原始信號相關(guān)性最高的3個(gè)分量重構(gòu)信號。Bearing1_1的振動(dòng)信號及其重構(gòu)信號如圖6所示，可以看出重構(gòu)信號邊緣的高頻噪聲明顯減少，從去噪信號中更容易提取軸承的相關(guān)特征信息。

圖5 各分量與原始信號的相關(guān)程度

圖6 Bearing1_1數(shù)據(jù)去噪前、后的時(shí)域圖

4.3 特征提取

計(jì)算重構(gòu)信號每個(gè)采樣點(diǎn)的最大值、平均值、峭度以及相關(guān)性最高3個(gè)分量的奇異值、樣本熵、能量熵，組成12維特征量，再利用t-SNE降維降至3維，進(jìn)行特征量可視化。

4.3.1 工況特征

取12套軸承的全壽命周期數(shù)據(jù)(Bearing3_1與Bearing3_3數(shù)據(jù)量過大，分別為2 538與1 515個(gè)數(shù)據(jù)點(diǎn)，故每隔10個(gè)點(diǎn)取樣一次，最后取得253與151個(gè)數(shù)據(jù)點(diǎn))，求得3維特征進(jìn)行可視化(圖7)，紅色、綠色、藍(lán)色的點(diǎn)分別用來代表工況1，2，3下的數(shù)據(jù)。從圖7可以看出，3種工況的數(shù)據(jù)呈不同的曲線分布，不同工況部分有較好的分離效果，重疊部分較少。

圖7 不同工況下軸承數(shù)據(jù)特征提取的可視化

4.3.2 健康度特征

將工況1下4套軸承所提取的3維特征進(jìn)行可視化，結(jié)果如圖8所示，圖下方的色譜條從藍(lán)色到紅色代表軸承數(shù)據(jù)從正常逐漸到故障，不同形狀代表工況1下4種不同的軸承數(shù)據(jù)。

(a)本文所提方法

圖8a是通過本文所提方法得到的可視化視圖，可以看出振動(dòng)信號從正常狀態(tài)至損壞狀態(tài)在三維空間中幾乎呈一條曲線分布， 具有很好的聚集性與連續(xù)性，特征點(diǎn)所處空間的位置能很大程度地反映軸承健康狀況。 圖8b由相關(guān)性最強(qiáng)的3個(gè)分量的最大值、標(biāo)準(zhǔn)值和峭度所組成的9維特征量降維所得，圖8c則由這3個(gè)分量的奇異值、能量熵和樣本熵所組成的9維特征量降維得到，可以看出這2種情況均不能很好表征軸承的退化情況。

4.4 故障分類與剩余壽命預(yù)測

4.4.1 故障分類

取4.3.1節(jié)得到的3維特征數(shù)據(jù)進(jìn)行故障分類試驗(yàn)，12套軸承共計(jì)2 977個(gè)數(shù)據(jù)點(diǎn)，將數(shù)據(jù)順序打亂，劃分前2 477個(gè)數(shù)據(jù)點(diǎn)為訓(xùn)練集，剩余500個(gè)數(shù)據(jù)點(diǎn)作為測試集，利用支持向量機(jī)對故障進(jìn)行分類處理，每個(gè)數(shù)據(jù)點(diǎn)的3維特征作為輸入，工況標(biāo)簽作為輸出。

支持向量機(jī)進(jìn)行分類預(yù)測時(shí)需要調(diào)節(jié)懲罰參數(shù)c與核函數(shù)參數(shù)g以得到最佳參數(shù)組合[c0,g0]，從而使分類預(yù)測效果最佳。同樣，利用蜉蝣優(yōu)化算法對支持向量機(jī)的這2個(gè)參數(shù)進(jìn)行優(yōu)化，適應(yīng)度函數(shù)為分類準(zhǔn)確率的相反數(shù)，當(dāng)其最小時(shí)，即分類準(zhǔn)確率最高。通過尋優(yōu)得到最佳參數(shù)組合c=100，g=0.1，分類結(jié)果如圖9所示，測試集中，499個(gè)數(shù)據(jù)點(diǎn)實(shí)現(xiàn)了正確分類，準(zhǔn)確率為99.8%，效果理想。

圖9 軸承故障的支持向量機(jī)分類結(jié)果

4.4.2 健康度預(yù)測

將軸承的健康度標(biāo)簽Y定義為

(7)

式中：Xnow為軸承當(dāng)前剩余壽命；Xall為軸承總壽命。初始健康標(biāo)簽Y為0，隨著軸承運(yùn)行狀態(tài)逐漸惡化，Y最后變?yōu)?，這樣就把軸承的健康度標(biāo)簽歸至0與1之間。

長短期記憶神經(jīng)網(wǎng)絡(luò)作為循環(huán)神經(jīng)網(wǎng)絡(luò)的變體，能夠解決循環(huán)神經(jīng)網(wǎng)絡(luò)長序列訓(xùn)練過程中的梯度消失和梯度爆炸問題，可以很好地捕捉時(shí)間序列上的特征，是用于剩余壽命預(yù)測的一種經(jīng)典神經(jīng)網(wǎng)絡(luò)模型。本文選擇長短期記憶神經(jīng)網(wǎng)絡(luò)預(yù)測軸承健康度，用提取的3維特征作為輸入，對應(yīng)的健康度標(biāo)簽作為輸出。

取第1個(gè)工況下4套軸承的的數(shù)據(jù)，共計(jì)564個(gè)采樣點(diǎn)，所提取的特征數(shù)據(jù)為R564*3，劃分Bearing1_1的123個(gè)數(shù)據(jù)為測試集，其他軸承的441個(gè)數(shù)據(jù)為訓(xùn)練集，訓(xùn)練長短期記憶神經(jīng)網(wǎng)絡(luò)并進(jìn)行測試。試驗(yàn)結(jié)果如圖10所示，僅通過當(dāng)前一個(gè)時(shí)間點(diǎn)振動(dòng)信號提取的特征，即可得到相對不錯(cuò)的健康度預(yù)測效果，驗(yàn)證了所構(gòu)建健康指標(biāo)的合理性。

圖10 工況1下試驗(yàn)軸承長短期記憶神經(jīng)網(wǎng)絡(luò)的 健康度預(yù)測結(jié)果

5 結(jié)束語

利用蜉蝣優(yōu)化算法優(yōu)化變分模態(tài)分解的參數(shù)，以各分量的最小包絡(luò)熵作為適應(yīng)度函數(shù)，可以獲得比遺傳算法和粒子群優(yōu)化算法更好的準(zhǔn)確性與收斂性；以最大值、平均值、峭度以及主要分量的奇異值、樣本熵、能量熵作為特征量組合，將其降維后含有軸承的工況與健康度信息在三維空間中直觀展現(xiàn)出來，通過三維健康指標(biāo)在空間中呈現(xiàn)一條代表軸承健康度逐漸變化的過程曲線，作為神經(jīng)網(wǎng)絡(luò)的輸入時(shí)能提供更多的特征量信息。