孫榮榮，楊宏偉，黃曉明

（濱州學(xué)院機(jī)電工程學(xué)院，山東 濱州 256600）

1 引言

航空產(chǎn)品的壽命統(tǒng)計(jì)是飛機(jī)可靠性分析工作的基礎(chǔ)，是優(yōu)化維修間隔，修訂維修方案的重要依據(jù)。由于飛機(jī)系統(tǒng)部件幾乎全部依賴進(jìn)口，缺少試驗(yàn)數(shù)據(jù)，其壽命統(tǒng)計(jì)的數(shù)據(jù)通常來自飛機(jī)運(yùn)營(yíng)過程中產(chǎn)生的外場(chǎng)數(shù)據(jù)。隨著飛機(jī)可靠性水平的提高，特別是一些具有安全性后果的航空產(chǎn)品，外場(chǎng)數(shù)據(jù)大部分為預(yù)防性維修數(shù)據(jù)，修復(fù)性維修只占很小的比例，由于這些預(yù)防性維修產(chǎn)品在維修時(shí)并未完全失效或已經(jīng)失效但未被發(fā)現(xiàn)，其壽命數(shù)據(jù)通常為右刪失數(shù)據(jù)或區(qū)間刪失數(shù)據(jù)。因此，如何充分有效地利用各種類型的外場(chǎng)數(shù)據(jù)，是航空產(chǎn)品可靠性分析工作的一個(gè)重要內(nèi)容[1]。

在醫(yī)學(xué)、經(jīng)濟(jì)學(xué)等領(lǐng)域，國(guó)內(nèi)外學(xué)者針對(duì)隨機(jī)右刪失[2]、區(qū)間刪失[3]等各種類型的刪失數(shù)據(jù)統(tǒng)計(jì)問題已經(jīng)進(jìn)行了大量研究。但是目前航空領(lǐng)域的壽命統(tǒng)計(jì)分析仍然以完全失效數(shù)據(jù)為主[4]，少數(shù)考慮到刪失數(shù)據(jù)[5]，尤其是對(duì)于隱蔽功能故障，往往忽略了其外場(chǎng)數(shù)據(jù)的區(qū)間刪失性[6]。另外，傳統(tǒng)的參數(shù)估計(jì)方法通常只能處理一種或兩種數(shù)據(jù)類型，而且對(duì)于存在大量刪失數(shù)據(jù)的樣本，且當(dāng)分布模型中包含多個(gè)分布參數(shù)時(shí)，似然函數(shù)的求解會(huì)變得非常復(fù)雜[7]。

這里指在建立一種混合數(shù)據(jù)類型下的壽命統(tǒng)計(jì)模型，既能處理完全失效或刪失數(shù)據(jù)中的某一數(shù)據(jù)類型，也能同時(shí)處理多種數(shù)據(jù)類型，并采用一種簡(jiǎn)單高效的粒子群優(yōu)化算法進(jìn)行模型求解，最后在獲得的多個(gè)壽命分布類型中選擇最佳分布類型，并利用計(jì)算機(jī)模擬進(jìn)行驗(yàn)證。

2 航空產(chǎn)品主要壽命分布類型

從統(tǒng)計(jì)理論的角度來看，產(chǎn)品的壽命是服從一定統(tǒng)計(jì)規(guī)律的隨機(jī)變量。國(guó)內(nèi)外研究表明，航空產(chǎn)品的壽命分布模型主要有以下幾種類型[8]。

2.1 指數(shù)分布

在民機(jī)系統(tǒng)可靠性分析中，指數(shù)分布是最基本的分布類型，其故障率為常數(shù)，特別適用于民機(jī)電子元器件和復(fù)雜項(xiàng)目的可靠性分析，其概率密度函數(shù)為：t

2.2 正態(tài)分布

正態(tài)分布主要用于描述那些由于磨損而發(fā)生故障的機(jī)件，其概率密度函數(shù)為：

2.3 對(duì)數(shù)正態(tài)分布

對(duì)數(shù)正態(tài)分布則常用于分析由于裂紋擴(kuò)展而引起故障的機(jī)件，若時(shí)間T服從對(duì)數(shù)正態(tài)分布，則X=lnT服從正態(tài)分布，因此，當(dāng)t＞0時(shí)，對(duì)數(shù)正態(tài)分布的概率密度函數(shù)為：

式中：μ—分布函數(shù)均值；σ—分布函數(shù)標(biāo)準(zhǔn)差。

2.4 威布爾分布

威布爾分布的適用范圍較其它幾種分布模型更為廣泛，它不僅適用于分析偶然故障，還適用于分析早期故障、耗損故障等數(shù)據(jù)類型。雙參數(shù)威布爾分布的概率密度函數(shù)為：

3 多數(shù)據(jù)類型下壽命統(tǒng)計(jì)方法

3.1 壽命統(tǒng)計(jì)模型的建立

3.1.1 完全失效數(shù)據(jù)

對(duì)于一些非計(jì)劃維修產(chǎn)品，是在產(chǎn)品完全失效后再采取維修工作，因此可以知道產(chǎn)品失效的準(zhǔn)確時(shí)間。

假設(shè)樣本i在某一時(shí)刻TFi完全失效，則樣本i的負(fù)對(duì)數(shù)似然函數(shù)可表示為：

3.1.2 右刪失數(shù)據(jù)

對(duì)于一些定時(shí)維修或存在潛在故障的產(chǎn)品，例如油濾的更換、輪胎的磨損等，往往在產(chǎn)品失效之前就采取了計(jì)劃維修。由于無法知道產(chǎn)品準(zhǔn)確的失效時(shí)間，但是可以確定產(chǎn)品在被刪失時(shí)刻是未失效的。假設(shè)樣本i在時(shí)刻TCi由于計(jì)劃維修被刪失，則該樣本i的對(duì)數(shù)似然函數(shù)可表示為：

3.1.3 區(qū)間刪失數(shù)據(jù)

部分航空產(chǎn)品的故障為隱蔽功能故障，例如一些備用設(shè)備等，正常情況下其功能故障對(duì)于工作人員是不明顯的，只有在下次預(yù)防性維修時(shí)才能被維修人員發(fā)現(xiàn)。假設(shè)樣本i在某一檢查間隔區(qū)間(TLi，THi)內(nèi)失效，并不知道其確切的失效時(shí)間，但是可以知道該樣本在時(shí)刻THi未失效且在時(shí)刻TLi失效，則該樣本i的負(fù)對(duì)數(shù)似然函數(shù)為：

如果樣本同時(shí)包含三種數(shù)據(jù)類型，則總樣本出現(xiàn)的概率為三種數(shù)據(jù)類型下的似然函數(shù)求和，因此，可得到總樣本的負(fù)對(duì)數(shù)似然函數(shù)為：

式中：n—樣本總個(gè)數(shù)。

將假設(shè)壽命分布模型的概率密度函數(shù)和累積分布函數(shù)代入式（8）～式（11）中，即可得到總樣本的負(fù)對(duì)數(shù)似然函數(shù)表達(dá)式。以威布爾分布為例，假設(shè)樣本含有N1個(gè)完全失效數(shù)據(jù)，N2個(gè)右刪失數(shù)據(jù)，N3組區(qū)間刪失數(shù)據(jù)，且每組刪失數(shù)據(jù)個(gè)數(shù)為k( )i，則三種數(shù)據(jù)類型下的負(fù)對(duì)數(shù)似然函數(shù)為：

3.2 基于PSO算法的模型求解

由于故障數(shù)據(jù)的復(fù)雜性和多樣性，混合數(shù)據(jù)下的負(fù)對(duì)數(shù)似然函數(shù)的求解無法用傳統(tǒng)的參數(shù)估計(jì)方法進(jìn)行求解。因此采用一種隨機(jī)搜索算法—粒子群優(yōu)化算法（PSO）進(jìn)行壽命統(tǒng)計(jì)模型的求解[10]。

為了獲得不同可靠性分布模型下最優(yōu)的參數(shù)估計(jì)，求解的目標(biāo)是搜索分布模型參數(shù)的值，以使負(fù)對(duì)數(shù)似然函數(shù)最小。因此，算法求解目的是尋找一定約束范圍下參數(shù)Z的值，使目標(biāo)函數(shù)式（11）達(dá)到最小。

基于PSO算法的參數(shù)估計(jì)是在構(gòu)建的負(fù)對(duì)數(shù)似然函數(shù)的基礎(chǔ)之上直接確定適用度函數(shù)，進(jìn)而以并行方式搜索全局極值和局部極值的方式來求其數(shù)字解。其求解的具體步驟為：

第一步：初始化種群的個(gè)體，設(shè)定分布模型各參數(shù)的取值范圍，并隨機(jī)初始化速度；

第二步：基于適應(yīng)度函數(shù)Fitness=-lnL( )θ，計(jì)算各粒子的適應(yīng)度，并初始化個(gè)體最優(yōu)位置Pi和全局最優(yōu)位置Pg；

第三步：依次迭代求解，不斷更新粒子，將每個(gè)粒子的適應(yīng)度與Pi和Pg進(jìn)行比較并對(duì)其進(jìn)行更新；

其中，第k代到k+1代的粒子更新速度為：

式中：α—速度加速因子，它從整體上對(duì)更新速度進(jìn)行加權(quán)。

第四步：如果滿足精度要求，則輸出目標(biāo)函數(shù)的解，否則返回第二步。

3.3 最佳壽命分布模型的選擇

由于航空產(chǎn)品故障機(jī)理較為復(fù)雜，不能簡(jiǎn)單地根據(jù)經(jīng)驗(yàn)假設(shè)某一特定的分布模型，這里將指數(shù)、威布爾、正態(tài)和對(duì)數(shù)正態(tài)四種分布模型分別代入壽命統(tǒng)計(jì)模型中進(jìn)行參數(shù)估計(jì)，然后基于BIC準(zhǔn)則選擇一個(gè)最佳的分布模型。BIC 準(zhǔn)則是日本統(tǒng)計(jì)學(xué)家Akaike 在1973 年提出的模型選擇的最小信息量準(zhǔn)則（Bayesian Information Criterion，BIC 準(zhǔn)則）［11］。BIC是一種評(píng)估模型最優(yōu)配置的指標(biāo)，它是擬合精度和未知參數(shù)個(gè)數(shù)的加權(quán)函數(shù)，使BIC值最小的分布模型被選擇為最佳分布類型，該準(zhǔn)則的計(jì)算公式如下：

通過以上分析，基于不同分布模型下參數(shù)估計(jì)的結(jié)果可計(jì)算出各備選分布模型的BIC值，并將其中BIC值最小的作為相對(duì)最優(yōu)的分布類型。根據(jù)航空產(chǎn)品壽命分布特點(diǎn)，在進(jìn)行參數(shù)估計(jì)時(shí)僅考慮了四種分布類型。在機(jī)械、醫(yī)學(xué)等其它領(lǐng)域應(yīng)用時(shí)，也可以考慮其它壽命分布模型進(jìn)行參數(shù)估計(jì)，然后基于BIC準(zhǔn)則選擇一個(gè)最佳的分布模型。

4 利用計(jì)算機(jī)模擬進(jìn)行模型驗(yàn)證

為了驗(yàn)證模型的合理性，以威布爾分布為例，取形狀參數(shù)β=2，尺度參數(shù)η=1000，利用MATLAB 產(chǎn)生100 個(gè)服從威布爾分布的失效數(shù)據(jù)，從這100個(gè)值中隨機(jī)選取10個(gè)值作為完全失效數(shù)據(jù)，并在剩下的90個(gè)數(shù)據(jù)中將＜1200的數(shù)據(jù)構(gòu)造區(qū)間刪失數(shù)據(jù)，區(qū)間長(zhǎng)度為100，將＞1200的數(shù)據(jù)構(gòu)造右刪失數(shù)據(jù)，獲得的樣本數(shù)據(jù)，如表1所示。

表1 計(jì)算機(jī)模擬樣本數(shù)據(jù)Tab.1 Computer Simulation Sample Data

下面采用本文提出的方法，對(duì)構(gòu)造的三種數(shù)據(jù)類型下的樣本進(jìn)行壽命統(tǒng)計(jì)模型的建立與求解。將表1中構(gòu)造的樣本數(shù)據(jù)代入式（12）～式（15），利用PSO算法和MATLAB工具進(jìn)行編程求解，選擇其適應(yīng)度函數(shù)Fitness=NLL，取學(xué)習(xí)因子c1=c2=1.4962，慣性權(quán)重ω=0.7298，最大迭代次數(shù)MaxDT=100，搜索空間維數(shù)（優(yōu)化參數(shù)個(gè)數(shù)）D=2，初始化種群個(gè)體數(shù)N=20。同樣的方法可進(jìn)行正態(tài)分布、對(duì)數(shù)正態(tài)分布和指數(shù)分布下的參數(shù)估計(jì)，從而得到四種分布類型下的參數(shù)估計(jì)結(jié)果和BIC值，如表2所示。

表2 參數(shù)估計(jì)結(jié)果和BIC值Tab.2 Parameter Estimation Results and BIC Values

四種分布函數(shù)模擬結(jié)果的迭代次數(shù)和最小適應(yīng)值關(guān)系以及擬合曲線與100個(gè)樣本數(shù)據(jù)的關(guān)系，如圖1、圖2所示。

圖1 迭代次數(shù)和最小適應(yīng)值關(guān)系Fig.1 Relations Between Iteration Number and Minimum Fitness Value

圖2 四種擬合分布函數(shù)與樣本數(shù)據(jù)關(guān)系Fig.2 Relations Between Four Fitting Distribution Functions and Sample Data

由參數(shù)估計(jì)和BIC準(zhǔn)則的分析結(jié)果可以看出，表1中數(shù)據(jù)的最佳擬合分布類型為威布爾分布，其形狀參數(shù)估計(jì)值β^=2.0253，尺度參數(shù)估計(jì)值η^=1031.2，近似等于實(shí)際值。

5 工程應(yīng)用實(shí)例

某航空公司E190型飛機(jī)的電傳系統(tǒng)備用電瓶故障為一隱蔽功能故障，針對(duì)該故障，該航空公司采取的維修策略是每1000飛行小時(shí)（FH）進(jìn)行一次操作測(cè)試。

因此，可以將備用電瓶在兩次檢查之間的故障數(shù)據(jù)視為區(qū)間刪失數(shù)據(jù)。

通過收集從2010年1月至2014年6月期間的外場(chǎng)維修數(shù)據(jù)并進(jìn)行篩選整理，共獲得28個(gè)故障數(shù)據(jù)，如表3所示。

表3 備用電瓶區(qū)間刪失數(shù)據(jù)Tab.3 Standby Battery Interval Censored Data

基于表3數(shù)據(jù)，利用這個(gè)方法進(jìn)行求解，可分別得到不同分布類型下的參數(shù)估計(jì)結(jié)果和BIC值，如表4所示。

表4 備用電瓶參數(shù)估計(jì)和BIC值Tab.4 Standby Battery Parameter Estimation and BIC Values

通過備用電瓶參數(shù)估計(jì)和BIC值的分析結(jié)果可以推斷出其壽命數(shù)據(jù)的最佳分布類型為正態(tài)分布，其概率密度函數(shù)為：

在獲得了產(chǎn)品的壽命分布模型之后，便可以進(jìn)行平均壽命、可用度等可靠性參數(shù)的求解，從而進(jìn)行維修決策。在本實(shí)例中，數(shù)據(jù)樣本只考慮了區(qū)間刪失數(shù)據(jù)，在后期運(yùn)營(yíng)過程中，如果獲得了產(chǎn)品更多的外場(chǎng)數(shù)據(jù)或通過試驗(yàn)等方法獲得了完全失效數(shù)據(jù)，還可以應(yīng)用本方法對(duì)結(jié)果進(jìn)行優(yōu)化處理。

6 結(jié)論

（1）提出了一種適用于多種數(shù)據(jù)類型的壽命統(tǒng)計(jì)方法，能同時(shí)處理完全失效、右刪失和區(qū)間刪失三種數(shù)據(jù)類型，解決了航空產(chǎn)品外場(chǎng)數(shù)據(jù)故障模式多樣化，存在大量刪失數(shù)據(jù)的問題。

（2）采用粒子群優(yōu)化算法進(jìn)行復(fù)雜似然函數(shù)的求解，實(shí)現(xiàn)了混合數(shù)據(jù)類型下的參數(shù)估計(jì)，并基于BIC準(zhǔn)則在多個(gè)備選壽命分布模型中選擇最佳分布模型。

（3）以威布爾分布為例，利用計(jì)算機(jī)模擬驗(yàn)證了模型的合理性，并以電傳系統(tǒng)備用電瓶為研究對(duì)象進(jìn)行了工程實(shí)例分析，說明了方法的合理性和實(shí)用性。