本發(fā)明涉及一種基于數(shù)學(xué)形態(tài)學(xué)和IFOA-SVR的滾動(dòng)軸承可靠度預(yù)測方法，涉及滾動(dòng)軸承可靠度預(yù)測技術(shù)領(lǐng)域。

背景技術(shù)：

滾動(dòng)軸承是旋轉(zhuǎn)機(jī)械中的關(guān)鍵部件，一旦出現(xiàn)故障將會(huì)造成大量的經(jīng)濟(jì)損失甚至危害人的生命安全^[1-2]。因此，準(zhǔn)確預(yù)測滾動(dòng)軸承在下一階段的工作狀態(tài)是合理制定機(jī)械設(shè)備維修計(jì)劃的前提基礎(chǔ)^[3-4]。

目前，滾動(dòng)軸承振動(dòng)信號(hào)的特征提取方法研究受到學(xué)者的廣泛關(guān)注。文獻(xiàn)[5]提出了基于形態(tài)分量分析和包絡(luò)譜的滾動(dòng)軸承故障診斷方法，該方法可有效提取滾動(dòng)軸承振動(dòng)信號(hào)中的故障特征。文獻(xiàn)[6]提出局部均值分解與形態(tài)學(xué)分形維數(shù)相結(jié)合的滾動(dòng)軸承故障診斷方法，可有效地對滾動(dòng)軸承進(jìn)行故障診斷。本發(fā)明將Hilbert包絡(luò)解調(diào)與形態(tài)學(xué)分形維數(shù)結(jié)合，對滾動(dòng)軸承全壽命數(shù)據(jù)提取出數(shù)學(xué)形態(tài)學(xué)分形維數(shù)，并將其作為性能退化狀態(tài)特征。

在建立支持向量回歸(Support vector regression,SVR)模型方面，文獻(xiàn)[7]利用粒子群算法(Particle swarm optimization,PSO)優(yōu)化最小二乘SVR對滾動(dòng)軸承的性能退化趨勢進(jìn)行了20點(diǎn)預(yù)測，平均絕對百分誤差(Mean absolute percentage error,MAPE)為15.82％。文獻(xiàn)[8]利用交叉驗(yàn)證法對C，g尋優(yōu)建立SVR模型，對滾動(dòng)軸承故障趨勢進(jìn)行了7點(diǎn)預(yù)測，MAPE為0.2％。文獻(xiàn)[9]通過經(jīng)驗(yàn)選取C，g參數(shù)建立SVR模型，并對電廠汽輪發(fā)電機(jī)組轉(zhuǎn)子的振動(dòng)幅值進(jìn)行了30點(diǎn)預(yù)測，平均絕對誤差(Mean absolute error,MAE)為0.6424，MAPE為1.48％。文獻(xiàn)[10]利用PSO對SVR模型中的C，g以及ε尋優(yōu)，預(yù)測航空飛行器的剩余壽命。文獻(xiàn)[7]到文獻(xiàn)[9]預(yù)測步長較短，文獻(xiàn)[8]和[9]只對SVR模型中的C和g進(jìn)行選取，并沒有提出對不敏感誤差ε的取值問題。同時(shí)，憑經(jīng)驗(yàn)選取C和g有很大的盲目性，通過交叉對比測試確定C和g的值則需要進(jìn)行大量的實(shí)驗(yàn)，而PSO相對復(fù)雜，并且容易陷入局部最優(yōu)，導(dǎo)致對模型參數(shù)的賦值不夠準(zhǔn)確。文獻(xiàn)[11]提出一種果蠅優(yōu)化算法(Fruit fly optimization algorithm,FOA)。文獻(xiàn)[12]提出遞減步長果蠅優(yōu)化算法(Diminishing step fruit fly optimization algorithm,DS-FOA)實(shí)現(xiàn)了全局搜索能力和局部尋優(yōu)能力的平衡。但DS-FOA在二維空間進(jìn)行搜索，不能真實(shí)反映果蠅的覓食行為。SVR模型內(nèi)部參數(shù)的選取會(huì)對SVR的預(yù)測精度以及預(yù)測步長產(chǎn)生較大影響，本發(fā)明采用在三維空間搜索的改進(jìn)果蠅優(yōu)化算法(Improved fruit fly optimization algorithm,IFOA)對SVR模型中的3個(gè)參數(shù)C、g以及ε同時(shí)進(jìn)行參數(shù)尋優(yōu)。并且在進(jìn)行預(yù)測時(shí)，將MAE、均方根誤差(Root mean square error,RMSE)、歸一化均方誤差(Normalized mean square error,NMSE)以及MAPE的和作為適應(yīng)度函數(shù)，相比于任一誤差作為適應(yīng)度函數(shù)具有更強(qiáng)的尋優(yōu)能力。

在可靠度模型方面，可靠度理論已被應(yīng)用于機(jī)械工程以及航空航天等各個(gè)領(lǐng)域^[13]。比例故障率模型的優(yōu)勢在于將軸承的性能退化特征與可靠度理論相結(jié)合，這樣，軸承的可靠度評估就能夠在設(shè)備狀態(tài)監(jiān)測數(shù)據(jù)的基礎(chǔ)上得到更新。文獻(xiàn)[14]提出將峭度、均方根作為特征指標(biāo)，利用威布爾比例故障率模型(Weibull proportional hazard model,WPHM)作為可靠度評估模型，并利用fminsearch優(yōu)化函數(shù)求解極大似然方程組確定WPHM的待定參數(shù)，有效地對鐵路機(jī)車輪對軸承進(jìn)行可靠度評估。文獻(xiàn)[15]在估計(jì)WPHM的待定參數(shù)時(shí)，利用牛頓迭代法對極大似然方程組進(jìn)行求解。但上述兩種求解方法都需要根據(jù)經(jīng)驗(yàn)設(shè)定搜索初值，并且運(yùn)算時(shí)間較長。本發(fā)明將極大似然方程組中每個(gè)方程的絕對值的和作為IFOA的適應(yīng)度函數(shù)，向適應(yīng)度函數(shù)最小的方向優(yōu)化，可以快速求出方程組的解，確定WPHM的待定參數(shù)。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明為了保證滾動(dòng)軸承運(yùn)行狀態(tài)可靠度的預(yù)測精度同時(shí)增加預(yù)測步長，提出一種數(shù)學(xué)形態(tài)學(xué)分形維數(shù)結(jié)合改進(jìn)果蠅優(yōu)化算法-支持向量回歸(Improved fruit fly optimization algorithm-support vector regression,IFOA-SVR)的滾動(dòng)軸承可靠度預(yù)測方法，即提出一種基于數(shù)學(xué)形態(tài)學(xué)和IFOA-SVR的滾動(dòng)軸承可靠度預(yù)測方法。

本發(fā)明為解決上述技術(shù)問題采取的技術(shù)方案是：

一種基于數(shù)學(xué)形態(tài)學(xué)和IFOA-SVR的滾動(dòng)軸承可靠度預(yù)測方法，所述方法的實(shí)現(xiàn)過程為：

步驟一、求出包絡(luò)信號(hào)的數(shù)學(xué)形態(tài)學(xué)分形維數(shù)D_M；

步驟二、從獲取的數(shù)學(xué)形態(tài)學(xué)分形維數(shù)D_M選取訓(xùn)練樣本對，基于所述訓(xùn)練樣本對利用IFOA對SVR模型中的參數(shù)C、g以及ε同時(shí)進(jìn)行尋優(yōu)，建立IFOA-SVR預(yù)測模型：

構(gòu)建IFOA的過程：

(1)初始化算法參數(shù)：設(shè)置果蠅種群規(guī)模Sizepop，最大覓食代數(shù)Maxgen，并隨機(jī)初始化果蠅群體位置坐標(biāo)(X₀,Y₀,Z₀)；

(2)果蠅個(gè)體利用嗅覺隨機(jī)搜索的方向和距離可以通過式(12)獲得

式中，i＝1,2,...,Sizepop，L₀為初始步長值，gen為當(dāng)前覓食代數(shù)；

(3)由于無法確定食物源的具體位置，所以需要通過式(13)估計(jì)第i個(gè)果蠅個(gè)體的當(dāng)前位置與坐標(biāo)原點(diǎn)間的距離Dist_i，之后計(jì)算出味道濃度判定值S_i

S_i＝1/Dist_i (14)

(4)將S_i代入味道濃度判定函數(shù)，計(jì)算出果蠅個(gè)體當(dāng)前位置的味道濃度

Smell_i＝function(S_i) (15)

(5)當(dāng)前果蠅群體中具有最高味道濃度的個(gè)體，可由式(16)獲得

[bestSmell,bestIndex]＝max(Smell_i) (16)

式中，bestSmell表示果蠅群體中具有最高味道濃度的個(gè)體的味道濃度值，bestIndex表示果蠅群體中具有最高味道濃度的個(gè)體的位置；

(6)保留果蠅群體中最佳味道濃度值和與其對應(yīng)的個(gè)體坐標(biāo)，同時(shí)果蠅群體利用自身的視覺對食物源進(jìn)行定位，然后飛往食物源所在的位置；

(7)進(jìn)入迭代尋優(yōu)過程，重復(fù)步驟(2)-(5)，并判斷當(dāng)前味最高道濃度是否好于前一迭代味道濃度，且gen<Maxgen；若成立，則執(zhí)行步驟(6)；

建立IFOA-SVR模型，其過程為：

(1)初始化IFOA參數(shù)，包括Sizepop、Maxgen以及果蠅個(gè)體初始位置；對SVR中的3個(gè)參數(shù)：懲罰系數(shù)C、核函數(shù)參數(shù)g以及不敏感誤差ε進(jìn)行尋優(yōu)，所以初始坐標(biāo)為(X₀¹,Y₀¹,Z₀¹)，(X₀²,Y₀²,Z₀²)和(X₀³,Y₀³,Z₀³)；

(2)附與每個(gè)果蠅個(gè)體隨機(jī)方向和飛行距離，并用遞減搜索步長代替固定步長，得到(X_i¹,Y_i¹,Z_i¹)，(X_i²,Y_i²,Z_i²)，(X_i³,Y_i³,Z_i³)，計(jì)算當(dāng)前果蠅個(gè)體與原點(diǎn)間距離的倒數(shù)，得到味道濃度判定值S_i¹，S_i²和S_i³；

(3)確定SVR中參數(shù)C、g和ε的范圍，即C∈[2^-14,2¹⁴]，g∈[2^-14,2¹⁴]，ε∈[2^-14,2¹⁴]；

(4)將所述訓(xùn)練樣本對輸入到SVR模型中，對模型進(jìn)行訓(xùn)練，將平均絕對誤差MAE、均方根誤差RMSE、歸一化均方誤差NMSE以及平均絕對百分誤差MAPE的和作為適應(yīng)度函數(shù)，即

(5)找到適應(yīng)度函數(shù)Fitness的最小值對應(yīng)的果蠅個(gè)體，開始迭代尋優(yōu)，并判斷最小Fitness是否低于前一代最小Fitness；如低于，則保留最小Fitness值及其對應(yīng)的坐標(biāo)，并將其賦給初始坐標(biāo)；如高于，則返回步驟(2)；

(6)找到C、g和ε的最佳值，建立IFOA-SVR預(yù)測模型；

步驟三、從所述性能退化狀態(tài)特征中截取一定長度Z作為WPHM的響應(yīng)協(xié)變量，結(jié)合極大似然估計(jì)獲得似然函數(shù)方程組，把方程組中每個(gè)方程絕對值的和作為IFOA的適應(yīng)度函數(shù)，求出方程組的解，確定WPHM的待定參數(shù)，進(jìn)而得到可靠度模型R(t,Z)，過程如下：

基于WPHM建立滾動(dòng)軸承性能退化狀態(tài)特征與可靠度之間的數(shù)學(xué)關(guān)系，WPHM的表達(dá)式為

式中β為形狀參數(shù)，η為尺度參數(shù)，μ為協(xié)變量回歸參數(shù)，t表示時(shí)間，Z為截取一定長度的性能退化狀態(tài)特征，即從數(shù)學(xué)形態(tài)學(xué)分形維數(shù)D_M中截取一定長度數(shù)據(jù)；

h(t,Z)和可靠度函數(shù)R(t,Z)之間的關(guān)系為

從而，R(t,Z)可表示為

利用極大似然估計(jì)法得到關(guān)于待定參數(shù)β、η、μ的方程組并利用IFOA對該方程組進(jìn)行求解，確定WPHM的表達(dá)式h(t,Z)的待定參數(shù)β、η、μ的最終取值；將參數(shù)β、η、μ代入式(21)即可建立可靠度模型；

步驟四、將性能退化狀態(tài)特征的最后一部分?jǐn)?shù)據(jù)作為IFOA-SVR預(yù)測模型的輸入，采用長期迭代預(yù)測法獲取性能退化狀態(tài)特征預(yù)測結(jié)果；

步驟五、將IFOA-SVR預(yù)測模型得出的預(yù)測結(jié)果作為響應(yīng)協(xié)變量嵌入到可靠度模型R(t,Z)中，即可計(jì)算出所述性能退化狀態(tài)特征預(yù)測結(jié)果所對應(yīng)的可靠度，實(shí)現(xiàn)滾動(dòng)軸承可靠度的預(yù)測。

在步驟一中，通過Hilbert變換計(jì)算出狀態(tài)監(jiān)測數(shù)據(jù)的包絡(luò)信號(hào)；

步驟一一、獲取滾動(dòng)軸承的狀態(tài)監(jiān)測數(shù)據(jù)：滾動(dòng)軸承的振動(dòng)信號(hào)；

步驟一二、通過Hilbert變換計(jì)算出狀態(tài)監(jiān)測數(shù)據(jù)的包絡(luò)信號(hào)；

步驟一三、利用式(4)求出每次采集振動(dòng)信號(hào)的包絡(luò)信號(hào)的數(shù)學(xué)形態(tài)學(xué)分形維數(shù)D_M，將D_M作為滾動(dòng)軸承性能退化狀態(tài)特征；

式中：D_M為信號(hào)的Minkowski-Bouligand維數(shù)即數(shù)學(xué)形態(tài)學(xué)分形維數(shù)，c為常數(shù)，對和log(1/λ)進(jìn)行最小二乘線性擬合即可得到對信號(hào)數(shù)學(xué)形態(tài)學(xué)分形維數(shù)D_M；λ＝1,2,…,λ_max，λ_max為最大分析尺度；

A_Se(λ)為定義在分析尺度λ下f(n)關(guān)于Se(n)膨脹和腐蝕運(yùn)算的覆蓋面積，f(n)表示滾動(dòng)軸承振動(dòng)信號(hào)的包絡(luò)信號(hào)；Se(n)為結(jié)構(gòu)元素，為一維離散向量；Se(n)的長度及λ_max的值通過控制變量法求得。

本發(fā)明的有益效果是：

本發(fā)明提出的一種基于數(shù)學(xué)形態(tài)學(xué)分形維數(shù)結(jié)合改進(jìn)果蠅算法-支持向量回歸(Improved fruit fly optimization algorithm-support vector regression,IFOA-SVR)的滾動(dòng)軸承可靠度預(yù)測方法，充分發(fā)揮IFOA-SVR預(yù)測模型的優(yōu)勢，對提取的滾動(dòng)軸承狀態(tài)特征數(shù)據(jù)進(jìn)行預(yù)測。將預(yù)測結(jié)果嵌入到已建立的可靠度模型中進(jìn)行滾動(dòng)軸承可靠度預(yù)測，在保證滾動(dòng)軸承運(yùn)行狀態(tài)可靠度預(yù)測精度的同時(shí)增加了預(yù)測步長。

該方法首先提取振動(dòng)信號(hào)的包絡(luò)信號(hào)，計(jì)算該包絡(luò)信號(hào)的數(shù)學(xué)形態(tài)學(xué)分形維數(shù)，將其作為滾動(dòng)軸承性能退化狀態(tài)特征；其次，利用IFOA對SVR中的參數(shù)C，g以及ε同時(shí)進(jìn)行尋優(yōu)，建立IFOA-SVR預(yù)測模型。同時(shí)，利用極大似然估計(jì)結(jié)合IFOA建立威布爾比例故障率模型(Weibull proportional hazard model,WPHM)，進(jìn)而得到可靠度模型；最后，將退化狀態(tài)特征作為IFOA-SVR預(yù)測模型的輸入，采用長期迭代預(yù)測法獲取特征預(yù)測結(jié)果，并將該結(jié)果嵌入到可靠度模型中，從而預(yù)測出軸承運(yùn)行狀態(tài)的可靠度。實(shí)驗(yàn)結(jié)果表明，利用所提方法對滾動(dòng)軸承可靠度進(jìn)行預(yù)測，能在保證預(yù)測精度的前提下增加預(yù)測步長。

具體而言，本發(fā)明將Hilbert包絡(luò)解調(diào)與形態(tài)學(xué)分形維數(shù)結(jié)合，對滾動(dòng)軸承全壽命數(shù)據(jù)提取出數(shù)學(xué)形態(tài)學(xué)分形維數(shù)，并將其作為性能退化狀態(tài)特征。本發(fā)明采用在三維空間搜索的改進(jìn)果蠅優(yōu)化算法(Improved fruit fly optimization algorithm,IFOA)對SVR模型中的3個(gè)參數(shù)C、g以及ε同時(shí)進(jìn)行參數(shù)尋優(yōu)。并且在進(jìn)行預(yù)測時(shí)，將MAE、均方根誤差(Root mean square error,RMSE)、歸一化均方誤差(Normalized mean square error,NMSE)以及MAPE的和作為適應(yīng)度函數(shù)，相比于任一誤差作為適應(yīng)度函數(shù)具有更強(qiáng)的尋優(yōu)能力。本發(fā)明將極大似然方程組中每個(gè)方程的絕對值的和作為IFOA的適應(yīng)度函數(shù)，向適應(yīng)度函數(shù)最小的方向優(yōu)化，可以快速求出方程組的解，確定WPHM的待定參數(shù)。

附圖說明

圖1是本發(fā)明中的果蠅三維空間覓食示意圖，圖2是本發(fā)明的可靠度預(yù)測方法流程圖；

圖3為本發(fā)明方法在應(yīng)用中的第860次采集數(shù)據(jù)的時(shí)域波形圖，圖4是本發(fā)明方法在應(yīng)用中對應(yīng)的第860次采集數(shù)據(jù)的包絡(luò)信號(hào)波形圖，圖5為本發(fā)明在應(yīng)用中第1點(diǎn)至第860點(diǎn)的可靠度曲線圖，圖6為本發(fā)明在應(yīng)用中GA優(yōu)化SVR模型得出的預(yù)測結(jié)果圖，圖7為本發(fā)明在應(yīng)用中PSO優(yōu)化SVR模型得出的預(yù)測結(jié)果圖，圖8為本發(fā)明在應(yīng)用中DS-FOA優(yōu)化SVR模型得出的預(yù)測結(jié)果圖，圖9為本發(fā)明在應(yīng)用中IFOA優(yōu)化SVR模型得出的預(yù)測結(jié)果圖，圖10為本發(fā)明在應(yīng)用中50點(diǎn)預(yù)測數(shù)據(jù)的可靠度趨勢曲線圖；

圖11為由包絡(luò)求分形維數(shù)與直接求分形維數(shù)的對比圖，本發(fā)明采用由由包絡(luò)求分形維數(shù)。

具體實(shí)施方式

具體實(shí)施方式一：如圖1至11所示，本實(shí)施方式針對所述的基于數(shù)學(xué)形態(tài)學(xué)和IFOA-SVR的滾動(dòng)軸承可靠度預(yù)測方法的實(shí)現(xiàn)過程具體描述如下：

1基于數(shù)學(xué)形態(tài)學(xué)的分形維數(shù)

數(shù)學(xué)形態(tài)學(xué)包括兩種基本算子，分別是膨脹運(yùn)算和腐蝕運(yùn)算。設(shè)原始信號(hào)f(n)和結(jié)構(gòu)元素Se(n)分別為定義在集合F＝{0,1,…,N-1}和集合G＝{0,1,…,M-1}上的兩個(gè)一維離散函數(shù)，且N≥M。在每一個(gè)分析尺度λ下，令Se(n)對f(n)進(jìn)行一次膨脹和腐蝕運(yùn)算，即：

式中：表示膨脹運(yùn)算，Θ表示腐蝕運(yùn)算，λ＝1,2,…,λ_max，λ_max為最大分析尺度。

定義在尺度λ下f(n)關(guān)于Se(n)膨脹和腐蝕運(yùn)算的覆蓋面積A_Se(λ)為

Maragos證明，A_Se(λ)滿足式(4)

式中D_M為信號(hào)的Minkowski-Bouligand維數(shù)即數(shù)學(xué)形態(tài)學(xué)分形維數(shù)，c為常數(shù)，對和log(1/λ)進(jìn)行最小二乘線性擬合即可得到對信號(hào)數(shù)學(xué)形態(tài)學(xué)分形維數(shù)的估計(jì)。

2支持向量回歸機(jī)

SVR算法一般用于解決數(shù)據(jù)回歸預(yù)測問題，通過建立一個(gè)最優(yōu)超平面，使?fàn)顟B(tài)空間內(nèi)的各個(gè)數(shù)據(jù)點(diǎn)距離該超平面最近，并將該超平面作為回歸模型進(jìn)行預(yù)測。

訓(xùn)練樣本對(x₁,y₁)，(x₂,y₂)，…，(x_r,y_r)，x_i∈Rⁿ，y_i∈R，x_i為輸入樣本，y_i為輸出樣本，i＝1,2,...,r。SVR的回歸函數(shù)表達(dá)式為

f(x)＝<ω·x_i>+b (5)

式中<ω·x_i>為ω和x_i的內(nèi)積。系數(shù)ω和b通過求解式(6)的最小值優(yōu)化問題獲得

由于擬合曲線必然會(huì)存在誤差，但要把誤差控制在一定的可容許范圍內(nèi)，ε為不敏感誤差，式(6)中L_ε為不敏感損失函數(shù)，其表達(dá)式為

對于回歸錯(cuò)誤的數(shù)據(jù)點(diǎn)，需要引入松弛變量ξ和ξ^*。將不敏感損失函數(shù)Lε代入式(6)可得

滿足約束條件

求解上式時(shí)，一般采用對偶理論將其轉(zhuǎn)化為二次規(guī)劃問題。對于非線性數(shù)據(jù)，引入非線性映射函數(shù)Φ，建立拉格朗日方程，通過化簡，可得式(8)的對偶式

式中α，α^*為拉格朗日乘子。

令K(x_i,x_j)＝Φ(x_i)^TΦ(x_j)為特征空間的內(nèi)積，稱K(x_i,x_j)為核函數(shù)，選擇徑向基函數(shù)作為核函數(shù)。根據(jù)KKT定理，可求得變量α_i，α_i^*，b，最終可得支持向量機(jī)的回歸函數(shù)為

3基于IFOA的SVR尋優(yōu)

3.1 IFOA算法

FOA是一種基于果蠅覓食行為推演出尋求全局優(yōu)化的算法。利用在三維空間搜索的IFOA對SVR模型進(jìn)行優(yōu)化，可以有效地獲取模型中的最優(yōu)參數(shù)。

如圖1所示為果蠅三維空間覓食示意圖。

IFOA步驟：

(1)初始化算法參數(shù)。設(shè)置果蠅種群規(guī)模Sizepop，最大覓食代數(shù)Maxgen，并隨機(jī)初始化果蠅群體位置坐標(biāo)(X₀,Y₀,Z₀)。

(2)果蠅個(gè)體利用嗅覺隨機(jī)搜索的方向和距離可以通過式(12)獲得

式中，i＝1,2,...,Sizepop，L₀為初始步長值，gen為當(dāng)前覓食代數(shù)。

(3)由于無法確定食物源的具體位置，所以需要通過式(13)估計(jì)第i個(gè)果蠅個(gè)體的當(dāng)前位置與坐標(biāo)原點(diǎn)間的距離Dist_i，之后計(jì)算出味道濃度判定值S_i

S_i＝1/Dist_i (14)

(4)將S_i代入味道濃度判定函數(shù)，計(jì)算出果蠅個(gè)體當(dāng)前位置的味道濃度

Smell_i＝function(S_i) (15)

(5)當(dāng)前果蠅群體中具有最高味道濃度的個(gè)體，可以由式(16)獲得

[bestSmell,bestIndex]＝max(Smell_i) (16)

(6)保留果蠅群體中最佳味道濃度值和與其對應(yīng)的個(gè)體坐標(biāo)，同時(shí)果蠅群體利用自身的視覺對食物源進(jìn)行定位，然后飛往食物源所在的位置。

(7)進(jìn)入迭代尋優(yōu)過程，重復(fù)步驟(2)-(5)，并判斷當(dāng)前味最高道濃度是否好于前一迭代味道濃度，且gen<Maxgen；若成立，則執(zhí)行步驟(6)。

3.2 IFOA-SVR模型

(1)初始化IFOA參數(shù)。包括Sizepop、Maxgen以及果蠅個(gè)體初始位置。由于本發(fā)明需要對SVR中的3個(gè)參數(shù)C、g以及ε進(jìn)行尋優(yōu)，所以初始坐標(biāo)為(X₀¹,Y₀¹,Z₀¹)，(X₀²,Y₀²,Z₀²)和(X₀³,Y₀³,Z₀³)。

(2)附與每個(gè)果蠅個(gè)體隨機(jī)方向和飛行距離，并用遞減搜索步長代替固定步長，得到(X_i¹,Y_i¹,Z_i¹)，(X_i²,Y_i²,Z_i²)，(X_i³,Y_i³,Z_i³)。計(jì)算當(dāng)前果蠅個(gè)體與原點(diǎn)間距離的倒數(shù)，得到味道濃度判定值S_i¹，S_i²和S_i³。

(3)確定SVR中參數(shù)C、g和ε的范圍，即C∈[2^-14,2¹⁴]，g∈[2^-14,2¹⁴]，ε∈[2^-14,2¹⁴]。

(4)將訓(xùn)練樣本對輸入到SVR模型中，對模型進(jìn)行訓(xùn)練，將MAE、RMSE、NMSE以及MAPE的和作為適應(yīng)度函數(shù)，即

(5)找到Fitness的最小值對應(yīng)的果蠅個(gè)體，開始迭代尋優(yōu)，并判斷最小Fitness是否低于前一代最小Fitness。如低于，則保留最小Fitness值及其對應(yīng)的坐標(biāo)，并將其賦給初始坐標(biāo)。如高于，則返回步驟(2)。

(6)找到C、g和ε的最佳值，建立IFOA-SVR預(yù)測模型。

4滾動(dòng)軸承運(yùn)行狀態(tài)可靠度

WPHM建立了設(shè)備運(yùn)行狀態(tài)特征指標(biāo)與可靠度之間的數(shù)學(xué)關(guān)系，表達(dá)式為

式中β為形狀參數(shù)，η為尺度參數(shù)，μ為協(xié)變量回歸參數(shù)，Z為振動(dòng)特征指標(biāo)數(shù)據(jù)。

h(t,Z)和可靠度函數(shù)R(t,Z)之間的關(guān)系為

從而，R(t,Z)可表示為

利用極大似然估計(jì)法對待定參數(shù)β、η、μ進(jìn)行估計(jì)，確定待定參數(shù)最終取值。

5滾動(dòng)軸承可靠度預(yù)測方法及流程

在WPHM參數(shù)確定的情況下，建立可靠度模型，將特征預(yù)測結(jié)果嵌入已建立的可靠度模型中，即可預(yù)測未來一段時(shí)間內(nèi)軸承運(yùn)行的可靠度，具體的預(yù)測流程見圖2。

(1)通過Hilbert變換計(jì)算出狀態(tài)監(jiān)測數(shù)據(jù)的包絡(luò)信號(hào)，利用式(4)，求出包絡(luò)信號(hào)的數(shù)學(xué)形態(tài)學(xué)分形維數(shù)D_M。

(2)選取適當(dāng)?shù)挠?xùn)練樣本對，并利用IFOA對SVR模型中的參數(shù)C、g以及ε同時(shí)進(jìn)行尋優(yōu)，建立IFOA-SVR預(yù)測模型。

(3)選取一定長度的振動(dòng)特征指標(biāo)數(shù)據(jù)作為WPHM的響應(yīng)協(xié)變量，結(jié)合極大似然估計(jì)獲得似然函數(shù)方程組，把方程組中每個(gè)方程絕對值的和作為IFOA的適應(yīng)度函數(shù)，求出方程組的解，確定WPHM的待定參數(shù)，進(jìn)而得到可靠度模型。

(4)將振動(dòng)特征指標(biāo)數(shù)據(jù)作為IFOA-SVR預(yù)測模型的輸入，采用長期迭代預(yù)測法^[16]獲取振動(dòng)特征預(yù)測結(jié)果。

(5)將IFOA-SVR預(yù)測模型得出的預(yù)測結(jié)果作為響應(yīng)協(xié)變量嵌入到可靠度模型中，即可計(jì)算出這些特征預(yù)測結(jié)果所對應(yīng)的可靠度，實(shí)現(xiàn)滾動(dòng)軸承可靠度的預(yù)測。

針對本明提出的方法進(jìn)行具體應(yīng)用與分析

滾動(dòng)軸承全壽命數(shù)據(jù)來自于Cincinnati大學(xué)^[17]，每隔10min進(jìn)行一次數(shù)據(jù)采集，共采集984次，每次采樣長度為20480點(diǎn)，采樣頻率為20KHz。圖3為第860次采集數(shù)據(jù)的時(shí)域波形，圖4為其對應(yīng)的包絡(luò)信號(hào)波形。

在每次采集的數(shù)據(jù)中，取中間的4096點(diǎn)進(jìn)行Hilbert變換，求出包絡(luò)信號(hào)，進(jìn)而計(jì)算出數(shù)學(xué)形態(tài)學(xué)分形維數(shù)，并作為一個(gè)特征數(shù)據(jù)點(diǎn)。通過控制變量法實(shí)驗(yàn)對比證明，當(dāng)λ_max＝11，Se(n)＝(0,0,0,0,0,0,0,0,0,0,0,0)時(shí)，利用式(4)求出包絡(luò)信號(hào)的數(shù)學(xué)形態(tài)學(xué)分形維數(shù)能夠作為刻畫滾動(dòng)軸承性能退化過程的特征。

將滾動(dòng)軸承性能退化狀態(tài)特征作為響應(yīng)協(xié)變量，取600點(diǎn)到860點(diǎn)之間的特征指標(biāo)值，采用極大似然估計(jì)，利用IFOA求得WPHM模型的待定參數(shù)，如表1所示。

表1 WPHM待定參數(shù)的估計(jì)結(jié)果

將待定參數(shù)代入式(19)，即可確定當(dāng)數(shù)學(xué)形態(tài)學(xué)分形維數(shù)作為響應(yīng)協(xié)變量時(shí)的WPHM，進(jìn)而建立可靠度模型，按照式(22)可以計(jì)算出滾動(dòng)軸承在任意時(shí)刻下的可靠度

式中m＝1,2,…,984，984表示數(shù)據(jù)采集總次數(shù)；

每個(gè)數(shù)據(jù)點(diǎn)等于10分鐘，數(shù)據(jù)點(diǎn)與時(shí)間存在對應(yīng)關(guān)系。

將第1點(diǎn)至第860點(diǎn)之間對應(yīng)的特征指標(biāo)值嵌入可靠度模型即可計(jì)算出在這段時(shí)間內(nèi)軸承運(yùn)行的可靠度，如圖5所示。

由圖5可知，滾動(dòng)軸承在525點(diǎn)左右出現(xiàn)早期故障，在525點(diǎn)以前軸承對應(yīng)的可靠度為1。之后可靠度開始逐漸減小，到860點(diǎn)時(shí)滾動(dòng)軸承可靠度已經(jīng)降到0.25附近，此時(shí)軸承處于故障程度比較嚴(yán)重的狀態(tài)。

利用IFOA對SVR中的參數(shù)C，g以及ε進(jìn)行尋優(yōu)，L₀＝30，Sizepop＝20，Maxgen＝100。訓(xùn)練數(shù)據(jù)取第680點(diǎn)到第860點(diǎn)，輸入向量維數(shù)為20，建立IFOA-SVR預(yù)測模型。表2為利用不同的優(yōu)化算法對SVR進(jìn)行尋優(yōu)得到的參數(shù)值。

表2利用不同算法尋優(yōu)得到的3個(gè)參數(shù)值

預(yù)測數(shù)據(jù)從第861點(diǎn)開始，共進(jìn)行了50次迭代運(yùn)算，通過長期迭代預(yù)測法得到最終預(yù)測結(jié)果。如果只優(yōu)化SVR模型中的C和g，預(yù)測誤差會(huì)很大并且預(yù)測步長較短。

利用GA對SVR模型中的C、g以及ε同時(shí)進(jìn)行尋優(yōu)，預(yù)測結(jié)果如圖6所示。

利用PSO對SVR模型中的C、g以及ε同時(shí)進(jìn)行尋優(yōu)，預(yù)測結(jié)果如圖7所示。

利用DS-FOA對SVR模型中的C、g以及ε同時(shí)進(jìn)行尋優(yōu)，預(yù)測結(jié)果如圖8所示。

利用IFOA對SVR模型中的C、g以及ε同時(shí)進(jìn)行尋優(yōu)，并將MAE、RMSE、NMSE以及MAPE的和作為適應(yīng)度函數(shù)，預(yù)測結(jié)果如圖9所示。

實(shí)際退化曲線在890點(diǎn)以后呈上升趨勢，由圖6可以看出，預(yù)測曲線在890點(diǎn)以后，明顯偏離實(shí)際退化曲線，預(yù)測效果不好；由圖7可以看出，預(yù)測曲線在890點(diǎn)以后上升趨勢不明顯；由圖8可以看出，預(yù)測曲線與實(shí)際退化曲線較為接近，但是在893至901點(diǎn)，與實(shí)際退化曲線相差較大；對比圖6至圖9可以看出，IFOA優(yōu)化SVR所得出的預(yù)測曲線與實(shí)際性能退化曲線更為接近，預(yù)測效果更好。

為了定量地評估預(yù)測結(jié)果，分別計(jì)算所預(yù)測50個(gè)點(diǎn)的MAE、RMSE、NMSE和MAPE，計(jì)算結(jié)果如表3所示。

表3 50點(diǎn)預(yù)測誤差分析

由表3可知，在本發(fā)明對SVR中3個(gè)參數(shù)均優(yōu)化的情況下，進(jìn)行50點(diǎn)預(yù)測時(shí)，與GA、PSO以及DS-FOA優(yōu)化的SVR預(yù)測模型相比，IFOA-SVR模型的4類預(yù)測誤差均最小，其4類誤差和也最小。證明IFOA具有更好的空間尋優(yōu)能力，能夠準(zhǔn)確地找到SVR模型中參數(shù)的最佳值，體現(xiàn)了IFOA-SVR模型在滾動(dòng)軸承性能退化趨勢預(yù)測方面的優(yōu)勢。

與文獻(xiàn)[7]到[9]對比分析可知，本發(fā)明所得實(shí)驗(yàn)結(jié)果MAPE低于文獻(xiàn)[7]，預(yù)測步長增加了30點(diǎn)；與文獻(xiàn)[8]相比，MAPE在同一數(shù)量級內(nèi)稍高，但預(yù)測步長增加了43點(diǎn)；相比于文獻(xiàn)[9]，在降低了MAE和MAPE的同時(shí)預(yù)測步長增加了20點(diǎn)。

將預(yù)測所得50點(diǎn)數(shù)據(jù)嵌入可靠度模型中，可以計(jì)算出這50點(diǎn)的可靠度值，可靠度曲線如圖10所示。

從圖10中可以看出，第861點(diǎn)的可靠度值在0.25附近，并且這50點(diǎn)的可靠度呈下降趨勢，到910點(diǎn)時(shí)可靠度已經(jīng)下降到0.08左右。這表明軸承運(yùn)行狀態(tài)很差，已經(jīng)出現(xiàn)非常嚴(yán)重的故障，應(yīng)該做好及時(shí)更換軸承的準(zhǔn)備。

通過對本發(fā)明方法的應(yīng)用得出如下結(jié)論：

(1)提出Hilbert包絡(luò)解調(diào)與數(shù)學(xué)形態(tài)學(xué)分形維數(shù)相結(jié)合的特征提取方法。通過控制變量法得出最大分析尺度λ_max以及結(jié)構(gòu)元素Se這兩個(gè)參數(shù)的最佳設(shè)置，計(jì)算出包絡(luò)信號(hào)的出數(shù)學(xué)形態(tài)學(xué)分形維數(shù)并作為刻畫滾動(dòng)軸承性能退化過程的特征。

(2)將MAE、RMSE、NMSE以及MAPE的和作為IFOA的適應(yīng)度函數(shù)，對SVR模型中的3個(gè)參數(shù)進(jìn)行尋優(yōu)，建立了IFOA-SVR模型。利用該模型對滾動(dòng)軸承的性能退化趨勢進(jìn)行了50點(diǎn)預(yù)測，與GA，PSO和DS-FOA優(yōu)化的SVR模型相比預(yù)測誤差最小。同時(shí)，與相關(guān)文獻(xiàn)進(jìn)行了橫向?qū)Ρ?，在保證預(yù)測精度的前提下，增加了預(yù)測步長。

(3)利用IFOA求解WPHM的待定參數(shù)，建立了可靠度模型。將IFOA-SVR模型所得的特征預(yù)測結(jié)果嵌入到可靠度模型中，得到可靠度預(yù)測曲線，能夠預(yù)測滾動(dòng)軸承在未來50點(diǎn)的運(yùn)行狀態(tài)可靠度。

本發(fā)明中涉及的參考文獻(xiàn)：

[1]何正嘉,曹宏瑞,訾艷陽,等.機(jī)械設(shè)備運(yùn)行可靠性評估的發(fā)展與思考[J].機(jī)械工程學(xué)報(bào),2014,50(2)：171-186.

HE Zhengjia,CAO Hongrui,ZI Yanyang,et al.Developments and Thoughts on Operational Reliability Assessment of Mechanical Equipment[J].Journal of Mechanical Engineering,

2014,50(2)：171-186.

[2]WANG Zhonglai,HUANG Hongzhong,LI Yanfeng.An approach to reliability assessment under degradation and shock process[J].IEEE Transactions on Reliability,2011,60(4)：852-863.

[3]何群,李磊,江國乾,等.基于PCA和多變量極限學(xué)習(xí)機(jī)的軸承剩余壽命預(yù)測[J].中國機(jī)械工程,2014,25(7)：984-988+989.

HE Qun,LI Lei,JIANG Guoqian,et al.Residual Life Predictions for Bearings Based on PCA and MELM[J].China Mechanical Engineering,2014,25(7)：984-988+989.

[4]申中杰,陳雪峰,何正嘉,等.基于相對特征和多變量支持向量機(jī)的滾動(dòng)軸承剩余壽命預(yù)測[J].機(jī)械工程學(xué)報(bào),2013,49(2)：183-189.

SHEN Zhongjie,CHEN Xuefeng,HE Zhengjia,et al.Remaining Life Predictions of Rolling Bearing Based on Relative Features and Multivariable Support Vector Machine[J].Journal of Mechanical Engineering,2013,49(2)：183-189.

[5]陳向民,于德介,李蓉.基于形態(tài)分量分析和包絡(luò)譜的軸承故障診斷[J].中國機(jī)械工程,2014,25(8)：1047-1052+1053.

CHEN Xiangmin,YU Dejie,LI Rong.Fault Diagnosis of Rolling Bearings Based on Morphological Componenet Analysis and Envelope Spectrum[J].China Mechanical Engineering,2014,25(8)：1047-1052+1053.

[6]張亢,程軍圣,楊宇.基于局部均值分解與形態(tài)學(xué)分形維數(shù)的滾動(dòng)軸承故障診斷方法[J].振動(dòng)與沖擊,2013,32(9)：90-94.

ZHANG Kang,CHENG Junsheng,YANG Yu.Roller bearing fault diagnosis based on local mean decomposition and morphological fractal dimension[J].Journal of Vibration and Shock,2013,32(9)：90-94.

[7]陳昌,湯寶平,呂中亮.基于威布爾分布及最小二乘支持向量機(jī)的滾動(dòng)軸承退化趨勢預(yù)測[J].振動(dòng)與沖擊,2014,33(20)：52-56.

CHEN Chang,TANG Baoping,LU Zhongliang.Degradation trend prediction of rolling bearings based on Weibull distribution and least squares support vector machine[J].Journal of Vibration and Shock,2014,33(20)：52-56.

[8]宋梅村,蔡琦.基于支持向量回歸的設(shè)備故障趨勢預(yù)測[J].原子能科學(xué)技術(shù),2011,45(8)：972-976.

SONG Meicun,CAI Qi.Fault Trend Prediction of Device Based on Support Vector Regression[J].Atomic Energy Science and Technology,2011,45(8)：972-976.

[9]朱霄珣.基于支持向量機(jī)的旋轉(zhuǎn)機(jī)械故障診斷與預(yù)測方法研究[D].北京：華北電力大學(xué),2013.

ZHU Xiaoxun.Research on Rotating Machine Fault Diagnosis and State Prediction Method Based on Support Vector Machine[D].Beijing：North China Electric Power University,2013.

[10]P.J.García Nieto,E.García-Gonzalo,F.Sánchez Lasheras et al.Hybrid PSO-SVM-based method for forecasting of the remaining useful life for aircraft engines and evaluation of its reliability[J].Reliability Engineering and System Safety,2015,138：219-231.

[11]PAN W T.A new fruit fly optimization algorithm:Taking the financial distress model as an example[J].Knowledge-Based Systems,2012,26(2)：69-74.

[12]寧劍平,王冰,李洪儒,等.遞減步長果蠅優(yōu)化算法及應(yīng)用[J].深圳大學(xué)學(xué)報(bào)理工版,2014,31(4)：367-373.

NING Jianping,WANG Bing,LI Hongru,et al.Research on and application of diminishing step fruit fly optimization algorithm[J].Journal of Shenzhen University Science and Engineering,2014,31(4)：367-373.

[13]張義民,孫志禮.機(jī)械產(chǎn)品的可靠性大綱[J].機(jī)械工程學(xué)報(bào),2014,50(14)：14-20.

ZHANG Yimin,SUN Zhili.The Reliability Syllabus of Mechanical Products[J].Journal of Mechanical Engineering,2014,50(14)：14-20.

[14]丁鋒,何正嘉,訾艷陽,等.基于設(shè)備狀態(tài)振動(dòng)特征的比例故障率模型可靠性評估[J].機(jī)械工程學(xué)報(bào),2009,45(12)：89-94.

DING Feng,HE Zhengjia,ZI Yanyang,et al.Reliability Assessment Based on Equipment Condition Vibration Feature Using Proportional Hazards Model[J].Journal of Mechanical Engineering,2009,45(12)：89-94.

[15]陳昌.基于狀態(tài)振動(dòng)特征的空間滾動(dòng)軸承可靠性評估方法研究[D].重慶：重慶大學(xué),2014.

CHEN Chang.Reliability Assessment Method for Space Rolling Bearing Based on Condition Vibration Feature[D].Chongqing：Chongqing University,2014.

[16]董紹江.基于優(yōu)化支持向量機(jī)的空間滾動(dòng)軸承壽命預(yù)測方法研究[D].重慶：重慶大學(xué),2012.

DONG Shaojiang.Research on Space Bearing Life Prediction Method Based on Optimized Support Vector Machine[D].Chongqing：Chongqing University,2012.

[17]Rexnord technical services,bearing data set,IMS,University of Cincinnati,NASA Ames Prognostics Data Repository,NASA Ames,Moffett Field,CA,2007[2014-04].http://ti.arc.nasa.gov/project/prognostics-data-repository.