孟 蕾 許愛(ài)強(qiáng) 董 超

（1.海軍航空工程學(xué)院飛行器工程系 煙臺(tái) 264001）（2.海軍航空工程學(xué)院科研部 煙臺(tái) 264001）（3.海軍航空工程學(xué)院基礎(chǔ)部 煙臺(tái) 264001）

1 引言

航空機(jī)電系統(tǒng)發(fā)生故障模式有兩種，退化故障和突發(fā)故障，其最終故障是兩種故障模式競(jìng)爭(zhēng)的結(jié)果。目前國(guó)內(nèi)外文獻(xiàn)研究競(jìng)爭(zhēng)故障問(wèn)題，都是在突發(fā)故障與退化故障兩者不相連的情況下，著重于從單邊退化角度展開(kāi)論述，是退化量、突發(fā)故障單一化的集中體現(xiàn)，退化量的多元化始終未能在學(xué)界占主導(dǎo)。對(duì)于這樣的問(wèn)題，本文主要是通過(guò)航空機(jī)電系統(tǒng)提出了退化特性，構(gòu)建了多元退化量的航空機(jī)電系統(tǒng)的競(jìng)爭(zhēng)和預(yù)測(cè)模型。

本文中，以最小二乘支持向量機(jī)預(yù)測(cè)法為基本計(jì)算工具，簡(jiǎn)要的對(duì)退化參數(shù)做了基本測(cè)算，如非線性、小樣本等；退化量、故障間存在著非常特殊的聯(lián)系，為尋找隱含性質(zhì)，引入了“位置-尺度模型”，得出突發(fā)故障與退化量的相關(guān)參數(shù)，進(jìn)而根據(jù)競(jìng)爭(zhēng)故障預(yù)測(cè)模型得到了航空機(jī)電系統(tǒng)未來(lái)一段時(shí)間內(nèi)的競(jìng)爭(zhēng)故障概率。立足案例，在比較研究之后，認(rèn)為本文在方法論上具有優(yōu)勢(shì)，滿(mǎn)足精度要求。

2 構(gòu)建航空機(jī)電系統(tǒng)競(jìng)爭(zhēng)故障預(yù)測(cè)模型

假設(shè)存在的性能特征參數(shù)數(shù)量是n，xi(t)、Ui和Li分別代表i(i=1，…，n)個(gè)性能特征參數(shù)的退化量與上、下閾值。測(cè)試起始時(shí)間為t，則有G(x，βi)（其中βi=(β1，β2，…，βk)，為該分布的參數(shù)向量），對(duì)應(yīng)的概率密度函數(shù)g(x，βi) 可表示為

工業(yè)部門(mén)對(duì)故障的概念進(jìn)行了界定，對(duì)性能特征參數(shù)的退化量的標(biāo)準(zhǔn)值和閾值范圍，有著嚴(yán)格的退化量限制，如果大于了規(guī)定的標(biāo)準(zhǔn)，則故障難免。假設(shè)第i個(gè)性能特征參數(shù)在發(fā)生了退化故障，則可依據(jù)下式計(jì)算其可能性：

因?yàn)楹娇諜C(jī)電系統(tǒng)本身就具有n個(gè)特定性能的參數(shù)，任一參數(shù)出現(xiàn)了退化故障，將會(huì)導(dǎo)致整個(gè)系統(tǒng)會(huì)產(chǎn)生故障。所以，可以把系統(tǒng)的故障視為由n個(gè)性能特征參數(shù)的退化故障通過(guò)競(jìng)爭(zhēng)而產(chǎn)生的。假設(shè)所有的性能和特征參數(shù)之間保持相對(duì)獨(dú)立，則在t時(shí)刻航空機(jī)電系統(tǒng)的退化故障的概率可以通過(guò)下面的公式進(jìn)行估算。

突發(fā)故障這一指標(biāo)不但和時(shí)刻t等有關(guān)聯(lián)，而且與t時(shí)刻某性能特征參數(shù)的退化量xi(t) 等指標(biāo)有關(guān)聯(lián)。在航空機(jī)電系統(tǒng)中，若突發(fā)故障必然由退化量引起，故障出現(xiàn)時(shí)間是Th，且在產(chǎn)生關(guān)系的臨界時(shí)間Th會(huì)伴隨失效率（λi(t，x)）的同時(shí)出現(xiàn)，則Th的可靠性關(guān)系式為

航空機(jī)電系統(tǒng)在t時(shí)刻的突發(fā)故障概率可表示為

對(duì)于競(jìng)爭(zhēng)故障出現(xiàn)的概率，其算式是

3 航空機(jī)電系統(tǒng)競(jìng)爭(zhēng)故障的預(yù)測(cè)流程

航空機(jī)電系統(tǒng)進(jìn)行競(jìng)爭(zhēng)故障預(yù)測(cè)的具體流程如下：首先確定性能退化參數(shù)與突發(fā)故障參數(shù)的分布類(lèi)型；然后預(yù)測(cè)退化的未知參數(shù)，確定下一階段性能退化參數(shù)的分布函數(shù)；在考慮突發(fā)故障與退化量的相關(guān)度的情況下，求出突發(fā)故障預(yù)測(cè)模型的未知參數(shù)；最后對(duì)航空機(jī)電系統(tǒng)退化與突發(fā)競(jìng)爭(zhēng)故障概率進(jìn)行預(yù)測(cè)。

3.1 分布函數(shù)的確定

性能退化參數(shù)分布類(lèi)型的檢驗(yàn)，本文采取Копмогоров檢驗(yàn)方法。假設(shè)其分布函數(shù)G(x，β)的分布類(lèi)型與時(shí)間t無(wú)關(guān)，只是參數(shù)向量β隨時(shí)間t的變化而改變，則G(x，β)即可看作β與時(shí)間t的函數(shù)，可表示為G(x，β(t))。在確定G(x，β(t))的分布種類(lèi)前，不妨借鑒工程實(shí)踐經(jīng)驗(yàn)，并且針對(duì)假設(shè)進(jìn)行分布擬合檢驗(yàn)，最終確定出G(x，β(t))的分布類(lèi)型。綜合考慮航空機(jī)電系統(tǒng)的雙向退化與單邊退化。通常情況下，都會(huì)將正態(tài)分布N(x;μ，σ2)作為理想化選擇，相應(yīng)的密度函數(shù)如式（6）所示：

式中，μ和σ2分別表示t時(shí)刻性能退化數(shù)據(jù)的均值和方差，則性能退化數(shù)據(jù)的分布參數(shù)向量β=(μ，σ2)，易知，β與時(shí)間t相關(guān)。

對(duì)航空機(jī)電系統(tǒng)性能退化數(shù)據(jù)分布類(lèi)型的檢驗(yàn)就變?yōu)橐粋€(gè)正態(tài)性檢驗(yàn)過(guò)程，可利用樣本(X1，X2，…，Xn)對(duì)總體的分布是否服從N(μ，σ2)進(jìn)行檢驗(yàn)。 樣本從小到大排列成：X(1)≤X(2)…≤X(n)。利用式（7）計(jì)算統(tǒng)計(jì)量W，式中ak(W)的值可查表獲得［1］。

然后在給定顯著性水平α與樣本容量n下，通過(guò)查表可得Wn，α。若W＜Wα，則拒絕H0，否則接受H0。

針對(duì)突發(fā)故障數(shù)據(jù)，本文根據(jù)工程經(jīng)驗(yàn)假設(shè)其服從威布爾分布。檢驗(yàn)方法仍采用Копмогоро對(duì)假設(shè)的航空機(jī)電系統(tǒng)突發(fā)故障數(shù)據(jù)分布類(lèi)型進(jìn)行驗(yàn)證，驗(yàn)證過(guò)程與性能退化數(shù)據(jù)分布類(lèi)型的檢驗(yàn)類(lèi)似。

3.2 參數(shù)的預(yù)測(cè)

3.2.1 退化故障參數(shù)預(yù)測(cè)

3.1節(jié)得出性能退化數(shù)據(jù)的分布類(lèi)型，從而得出性能退化分布參數(shù)。根據(jù)歷史性能退化數(shù)據(jù)求解分布參數(shù)的歷史數(shù)據(jù)值，求解將來(lái)參數(shù)，確定分布函數(shù)，對(duì)退化故障、發(fā)生幾率做預(yù)測(cè)。

航空機(jī)電系統(tǒng)的性能退化數(shù)據(jù)的獲取主要基于狀態(tài)監(jiān)測(cè)系統(tǒng)或者是加速退化實(shí)驗(yàn)，本文研究的性能退化數(shù)據(jù)主要是來(lái)源于狀態(tài)監(jiān)測(cè)系統(tǒng)獲取的數(shù)據(jù)。通過(guò)測(cè)試取得了相應(yīng)的性能數(shù)據(jù)是和時(shí)間相關(guān)的指標(biāo)，分布參數(shù)也將會(huì)隨著時(shí)間的變化而調(diào)整。為此，可從時(shí)間序列角度審視。根據(jù)上文相關(guān)論述，可以很清晰得知：從時(shí)間序列角度展開(kāi)預(yù)測(cè)性能退化故障活動(dòng)將顯得非常的簡(jiǎn)便與合理。

在上述領(lǐng)域中，其主要的方法為神經(jīng)網(wǎng)絡(luò)、時(shí)間序列分析、支持向量機(jī)等［2～3］。由于航空機(jī)電系統(tǒng)性能大多數(shù)分布參數(shù)都集中體現(xiàn)出了非線性性質(zhì)，這決定了只有應(yīng)用統(tǒng)計(jì)學(xué)原理方可得出理想結(jié)果。支付向量機(jī)（SupportVector Machine，SVM）在處理小樣本、高維數(shù)的非線性問(wèn)題方面具有很強(qiáng)的適用性，但是SVM算法隨著樣本量的增大，運(yùn)算速度會(huì)降低［4］。最小二乘支持向量機(jī)（LV-SVM）算法既節(jié)約了計(jì)算時(shí)間又降低了計(jì)算難度。因此，本文選用LV-SVM預(yù)測(cè)模型對(duì)航空機(jī)電系統(tǒng)各性能退化數(shù)據(jù)的分布參數(shù)時(shí)間序列進(jìn)行預(yù)測(cè)［5］。

上式（8）中，αi為拉格朗日乘子，K為核函數(shù)，b?R為常數(shù)。由于，則預(yù)測(cè)模型可表示為

3.2.2 突發(fā)故障參數(shù)預(yù)測(cè)

本文將突發(fā)故障時(shí)間當(dāng)作響應(yīng)變量，在回歸函數(shù)中，退化量并不需要被排除。這說(shuō)明退化量、故障出現(xiàn)時(shí)間的研究能夠被放置于回歸模型中討論。借助回歸分析，還需要尋找理想的模型——位置-尺度模型［6］，而且它在退化量x等參數(shù)的基礎(chǔ)上對(duì)Y=lnT的分布進(jìn)行了描述，模型可以使用下面的關(guān)系式進(jìn)行表述：

當(dāng)中，e的分布和x沒(méi)有聯(lián)系。H(e)-e的可靠度函數(shù)。μ(x)-位置參數(shù)，σ＞0-不變的尺度參數(shù)，x，Y的可靠度函數(shù)為

針對(duì)其中的突發(fā)故障時(shí)間Th，使用Yh=lnTh的分布，則Th=exрYh的可靠度函數(shù)可表示為

式中，α(x)=exр(μ(x))，δ=1/σ，S(w)=H(ln(w))。利用該模型，式（6）可進(jìn)一步表示為

由3.1節(jié)分析可知，航空機(jī)電突發(fā)故障時(shí)間Th服從威布爾分布［7］?，F(xiàn)對(duì)給定退化量x下Yh=lnTh的分布進(jìn)行分析，由于Th吻合于威布爾分布，則Yh=lnTh服可以吻合于極值分布［8］，其密度函數(shù)使用以下的公式進(jìn)行計(jì)算：

式中，位置參數(shù)μ(x)=lnη(x)，尺度參數(shù)σ=1/m，若采用標(biāo)準(zhǔn)極值分布［9］進(jìn)行表示，則Yh如式（15）所示：

式中，e的密度函數(shù)是exр(e-exр(e))。μ(x)以線性為主，令μ(x)=γ1+γ2·x，σ=γ3，則式（16）可表示為

在綜合研究了諸多的機(jī)電系統(tǒng)后，發(fā)現(xiàn)對(duì)象系統(tǒng)的故障并不是單一存在的，而是事件、時(shí)間的綜合表現(xiàn)，即tj(j=1，2，…，K)（退化量）必然在某個(gè)tj(j=1，2，…，K)（突發(fā)故障時(shí)間）中出現(xiàn)，此類(lèi)關(guān)系被稱(chēng)為似然關(guān)系［10-12］，由式（17）函數(shù)來(lái)闡述突發(fā)故障數(shù)據(jù)：

4 實(shí)例分析

某海軍航空飛行團(tuán)生產(chǎn)的某型機(jī)電系統(tǒng)綜合ATS測(cè)試數(shù)據(jù)為例，隨機(jī)抽取10組機(jī)電系統(tǒng)某特征參數(shù)作為樣本進(jìn)行故障預(yù)測(cè)。其測(cè)試信息進(jìn)行測(cè)試的時(shí)間從2010年持續(xù)到2016年，主要是針對(duì)性能退化數(shù)據(jù)和突發(fā)故障數(shù)據(jù)進(jìn)行測(cè)試，因此可根據(jù)2010年～2015年的測(cè)試信息對(duì)機(jī)電系統(tǒng)2016年的故障概率進(jìn)行預(yù)測(cè)，而且把結(jié)果和機(jī)電系統(tǒng)的真實(shí)情況進(jìn)行對(duì)比，以以便可以檢驗(yàn)本次系統(tǒng)設(shè)計(jì)的有效性和針對(duì)性。

結(jié)合在工程中的相關(guān)經(jīng)驗(yàn)，可假定航空機(jī)電系統(tǒng)的一些性能特征的參數(shù)退化中的數(shù)據(jù)可以吻合于正態(tài)分布。以提取的十組機(jī)電系統(tǒng)中的特征性能參數(shù)2010年的相關(guān)的數(shù)據(jù)為實(shí)例，并最終得到以下退化數(shù)據(jù)（抽取而來(lái)）：9.93，10.61，11.00，10.43，12.10，10.22，10.53，10.55，10.67，11.12。按照分布擬合法對(duì)抽取值開(kāi)展驗(yàn)算，并最終得到分布狀態(tài)。

在集中整理了測(cè)試信息后，設(shè)該性能退化數(shù)據(jù)為X，則檢驗(yàn)X服從正態(tài)分布與否就變?yōu)闄z驗(yàn)假設(shè)H0:X～N(μ，σ2)是否成立。按升序排列性能退化數(shù)據(jù)，并對(duì)其進(jìn)行正態(tài)性檢驗(yàn)，為方便計(jì)算，列表如下。

表1 性能退化數(shù)據(jù)的正態(tài)性檢驗(yàn)

由表1可計(jì)算得到

將上述計(jì)算結(jié)果代入式（7），可得：W=0.8907。

若取顯著性水平α=0.05，則通過(guò)查閱W檢驗(yàn)統(tǒng)計(jì)量W的α分位數(shù)表Wα可得W10，0.05=0.842。由于W＜W10，0.05，所以在顯著性水平α=0.05的情況下，可以接受原假設(shè)H0，即可確定2010年測(cè)試時(shí)某的一些性能參數(shù)的特征和正態(tài)分布相吻合。

按照提取的十組機(jī)電系統(tǒng)某的一些性能參數(shù)，針對(duì)2010～2015的性能參數(shù)進(jìn)行退化數(shù)據(jù)的正態(tài)性檢驗(yàn)，結(jié)果表明所有的年份的數(shù)據(jù)都吻合正態(tài)分布。此外，正態(tài)性測(cè)試也適合剩余數(shù)據(jù)，相關(guān)操作方式如上。通過(guò)測(cè)算，發(fā)現(xiàn)如果結(jié)果出現(xiàn)了α=0.05的情況，其余性能特征參數(shù)的性能中的退化數(shù)據(jù)也吻合于正態(tài)分布。同理，結(jié)合工程中的經(jīng)驗(yàn)，假定機(jī)電系統(tǒng)的突發(fā)故障數(shù)據(jù)仍然使用Копмогоров方法對(duì)其中的分布類(lèi)型進(jìn)行證實(shí)，可得機(jī)電系統(tǒng)的突發(fā)故障參數(shù)服從威布爾分布。

由3.2節(jié)分析結(jié)果得出，航空機(jī)電系統(tǒng)性能退化參數(shù)可以吻合于正態(tài)分布，所以在針對(duì)性能參數(shù)進(jìn)行預(yù)測(cè)的過(guò)程中，也只是需要針對(duì)預(yù)測(cè)均值和方差等指標(biāo)進(jìn)行測(cè)試。

經(jīng)過(guò)對(duì)數(shù)據(jù)的分析處理，從2007～2015年中，提取了十組機(jī)電系統(tǒng)某性能特征參數(shù)，并對(duì)此類(lèi)參數(shù)進(jìn)行特殊處理，得到方差、均值。其時(shí)間序列為：（10.7160，10.7530，10.7960，10.8000，10.8160，10.8240，10.9340，11.0020，11.0870）和（0.3534，0.3771，0.6871，0.5370，0.5003，0.4050，1.0234，0.5728，0.5980）。首先預(yù)測(cè)均值μ，選取m=1，相空間重構(gòu)均值序列的前5個(gè)數(shù)據(jù)，可得樣本對(duì)Xμ=[10.7160，10.7530，10.7960，10.8000，10.8160，10.8240]T、Yμ=[10.7530，10.7960，10.8000，10.8160，10.8240，10.9340]T。選取RBF函數(shù)作為核函數(shù)，設(shè)定參數(shù)γ和δ的搜索區(qū)域分別為[0.1，500]和[0.01，100]，采用交叉檢驗(yàn)和網(wǎng)格搜索算法進(jìn)行求解，可得最優(yōu)參數(shù)δ2=3.2463，γ=116.4963。使用LS-SVM訓(xùn)練對(duì)樣本執(zhí)行輸入和輸出等操作，可以取得均值μ的預(yù)測(cè)的回歸直線，如下圖。得到訓(xùn)練完畢的LS-SVM回歸函數(shù)之后，即可開(kāi)始預(yù)測(cè)輸入樣本10.7530， 10.7960， 10.8000， 10.8160， 10.8240，10.9340]T。預(yù)測(cè)結(jié)果顯示實(shí)測(cè)值與預(yù)測(cè)值基本吻合：同理，可以預(yù)測(cè)性能退化數(shù)據(jù)的方差時(shí)間序列，結(jié)果顯示絕對(duì)誤差小于-0.05，相對(duì)誤差小于0.009。

圖1 均值預(yù)測(cè)回歸曲線

分析2007年～2015年的故障數(shù)據(jù)可以發(fā)現(xiàn)，提取了十組機(jī)電系統(tǒng)的記錄，持續(xù)的時(shí)間為2008年～2015年，其故障的總數(shù)分別為2,2,2,2,1,1,1,1。因?yàn)閷?duì)于機(jī)電系統(tǒng)進(jìn)行了實(shí)時(shí)的測(cè)試，所以需要以天為單位進(jìn)行測(cè)試，所以提取的十組數(shù)據(jù)的突發(fā)時(shí)間可以表示為具體的天數(shù)。以某退化特征參數(shù)的退化量為例，與突發(fā)故障時(shí)間的對(duì)應(yīng)關(guān)系可得，同時(shí)可以得出突發(fā)故障與退化量的相關(guān)參數(shù)的估計(jì)值。明確了數(shù)據(jù)的分類(lèi)模型和相應(yīng)的參數(shù)，即可根據(jù)式（13）對(duì)機(jī)電系統(tǒng)在2015年和2016年的故障概率進(jìn)行預(yù)測(cè)。結(jié)果發(fā)現(xiàn)，實(shí)際進(jìn)行測(cè)評(píng)的結(jié)果和預(yù)測(cè)值是基本上一致的，所以建立的故障預(yù)測(cè)模型是合理而有效的。

5 結(jié)語(yǔ)

文章中，借助于航空機(jī)電系統(tǒng)中關(guān)系復(fù)雜的退化性質(zhì)，集中論述了退化發(fā)生頻率、故障出現(xiàn)時(shí)間的數(shù)學(xué)關(guān)聯(lián)，引入航空系統(tǒng)的故障和預(yù)測(cè)模型，而且針對(duì)退化數(shù)據(jù)中的分布參數(shù)的小樣本、非線性的特征，使用最小支持向量機(jī)預(yù)測(cè)的算法針對(duì)故障進(jìn)行了預(yù)測(cè)，基于退化量與突發(fā)故障的相關(guān)性，預(yù)測(cè)了突發(fā)故障預(yù)測(cè)模型相關(guān)參數(shù)。最后實(shí)例分析結(jié)果表明：本文設(shè)計(jì)的方法預(yù)測(cè)精度更高，更加合理有效。

［1］吳翊，李永樂(lè)，胡慶軍.應(yīng)用數(shù)理統(tǒng)計(jì)［M］.北京：國(guó)防科技大學(xué)出版社，2008.

［2］張廣明，袁宇浩，龔松建.基于改進(jìn)最小二乘支持向量機(jī)方法的短期風(fēng)速預(yù)測(cè)［J］.上海交通大學(xué)學(xué)報(bào)，2011，45（8）：1125-1129.

［3］叢林虎，徐廷學(xué)，楊繼坤，董琪.導(dǎo)彈退化故障預(yù)測(cè)方法研究［J］.電光與控制，2014，21（5）：78-82.

［4］洪杰，韓磊，苗學(xué)問(wèn)，馬艷紅.基于支持向量機(jī)的滾動(dòng)軸承狀態(tài)壽命評(píng)估［J］.北京航空航天大學(xué)學(xué)報(bào)，2010，36（8）：896-899.

［5］唐杰明，劉俊勇，楊可等.基于灰色模型和最小二乘支持向量機(jī)的電力短期負(fù)荷組合預(yù)測(cè)［J］.電網(wǎng)技術(shù)，2009，33（3）：63-65.

［6］Lawless J F.壽命數(shù)據(jù)中的統(tǒng)計(jì)模型與方法［M］.北京：中國(guó)統(tǒng)計(jì)出版社，1998.

［7］蘇春，張恒.基于性能退化數(shù)據(jù)和競(jìng)爭(zhēng)失效分析的可靠性評(píng)估［J］.機(jī)械強(qiáng)度，2011，33（2）：196-200.

［8］孫紹輝.面向預(yù)知維修的航空發(fā)動(dòng)機(jī)可靠性評(píng)估與預(yù)測(cè)研究［D］.南京：南京航空航天大學(xué)，2013.

［9］GOLDSTEIN M，WOOFF D.Bayes linear statistics，theo?ry and methods［M］.New York：Wiley Publisher，2007.

［10］喻天翔，孫玉秋，張祖明.多模式失效的機(jī)械零件可靠度計(jì)算新理論［J］.機(jī)械工程學(xué)報(bào)，2003（03）：34-35.

［11］蘇春，張恒.基于性能退化數(shù)據(jù)和競(jìng)爭(zhēng)失效分析的可靠性評(píng)估［J］.機(jī)械強(qiáng)度，2011（02）：22-23.

［12］李方方，趙英凱，顏昕.基于Matlab的最小二乘支持向量機(jī)的工具箱及其應(yīng)用［J］.計(jì)算機(jī)應(yīng)用，2006（26）：358-360.