部分狀態(tài)信息下的截尾數(shù)據(jù)分析

田 英,佘 陽(yáng)

(1.佛山科學(xué)技術(shù)學(xué)院 機(jī)械與電氣工程學(xué)院,廣東 佛山 528000;2.廣東南車(chē)軌道交通車(chē)輛有限公司，廣東 江門(mén) 529000)

?

部分狀態(tài)信息下的截尾數(shù)據(jù)分析

田 英1,佘 陽(yáng)2

(1.佛山科學(xué)技術(shù)學(xué)院 機(jī)械與電氣工程學(xué)院,廣東 佛山 528000;2.廣東南車(chē)軌道交通車(chē)輛有限公司，廣東 江門(mén) 529000)

針對(duì)暗含部分狀態(tài)信息的類(lèi)截尾數(shù)據(jù)在可靠性建模中的分析與處理問(wèn)題，提出等失效數(shù)據(jù)法。該方法通過(guò)給類(lèi)截尾數(shù)據(jù)增加一個(gè)代表的剩余壽命，使其轉(zhuǎn)換為等失效數(shù)據(jù)，代表的剩余壽命值可能是一個(gè)固定值或是與平均剩余壽命成一定比例的值，對(duì)這兩種情況分別引入中間變量表示暗含信息，由極大似然法得出最優(yōu)中間變量和分布參數(shù)，通過(guò)實(shí)例分析驗(yàn)證了該方法的有效性。

類(lèi)截尾數(shù)據(jù)；部分狀態(tài)信息；代表的剩余壽命；預(yù)防維修；中間變量

截尾數(shù)據(jù)在可靠性領(lǐng)域中普遍存在，其主要特征是不完整性。右截尾數(shù)據(jù)是典型的不完整數(shù)據(jù)，指在數(shù)據(jù)收集的時(shí)間窗口內(nèi)，只知道失效時(shí)間大于當(dāng)前觀(guān)察到的運(yùn)行時(shí)間的數(shù)據(jù)[1]，由于這種觀(guān)察是隨機(jī)的，因此這種數(shù)據(jù)除了知道系統(tǒng)的工作狀態(tài)外無(wú)任何其他狀態(tài)信息。然而，有些截尾數(shù)據(jù)卻隱含某些狀態(tài)信息，例如，預(yù)防維修活動(dòng)中，元件被檢測(cè)到即將失效而進(jìn)行預(yù)防替換時(shí)產(chǎn)生的截尾數(shù)據(jù)，它不是隨機(jī)右截尾數(shù)據(jù)，且不能表示元件的壽命，但卻隱含元件不久將失效的重要信息。為了把這類(lèi)數(shù)據(jù)與隨機(jī)右截尾數(shù)據(jù)區(qū)分開(kāi)，筆者稱(chēng)其為類(lèi)右截尾數(shù)據(jù)。

令t+表示截尾數(shù)據(jù)，如果t+為隨機(jī)右截尾數(shù)據(jù)時(shí)，將其試驗(yàn)到失效，元件壽命可表示為t′=t++m(t+)，其中m(t+)為平均剩余壽命(mean residual life,MRL)；如果t+為類(lèi)右截尾數(shù)據(jù)時(shí)，由于元件即將失效，則其實(shí)際失效時(shí)間可能遠(yuǎn)小于t′。因此，在可靠性分析與建模中，不能忽略這些重要的狀態(tài)信息而僅簡(jiǎn)單地把這類(lèi)數(shù)據(jù)等同于隨機(jī)截尾數(shù)據(jù)或失效數(shù)據(jù)。為有效合理地應(yīng)用包含部分狀態(tài)信息的截尾數(shù)據(jù)，筆者提出等失效數(shù)據(jù)法對(duì)類(lèi)截尾數(shù)據(jù)進(jìn)行分析，并通過(guò)實(shí)例分析來(lái)驗(yàn)證該方法的合理性與有效性。

1 等失效數(shù)據(jù)法

1.1 基本思想

等失效數(shù)據(jù)法的基本思想是在類(lèi)截尾數(shù)據(jù)上增加一個(gè)代表的剩余壽命(representative residual life,RRL)，把類(lèi)截尾數(shù)據(jù)轉(zhuǎn)換為等失效數(shù)據(jù)，確定RRL值包含兩種情況：①是一個(gè)固定值;②是與MRL成一定比例的值。對(duì)這兩種情況分別引入一個(gè)中間變量來(lái)表示隱含的部分狀態(tài)信息，應(yīng)用極大似然法和Excel的“規(guī)劃求解”工具可快速而精確地求得中間變量和分布參數(shù)。

(1)

式中,δj表示RRL，若δj被指定，則可根據(jù)極大似然估計(jì)法來(lái)求得各模型參數(shù)。

(2)

式中: f(t)為概率分布函數(shù)；θ為其參數(shù)集。根據(jù)δj的大小，可分為兩種特殊情況：

(1) δj=0，即截尾數(shù)據(jù)被視為失效數(shù)據(jù)，這時(shí)將低估產(chǎn)品壽命。

顯然，δj值應(yīng)居于兩種特殊情況之間，即：

(3)

1.2 等失效數(shù)據(jù)法1

假設(shè)維修工程師對(duì)某元件代表的剩余壽命有一個(gè)粗略參考值，元件運(yùn)行到某個(gè)時(shí)間后，如果估計(jì)的剩余壽命小于該參考值，則對(duì)元件執(zhí)行預(yù)防替換；否則，元件繼續(xù)運(yùn)行。設(shè)δ為一個(gè)比參考值稍小的固定值，則類(lèi)截尾數(shù)據(jù)可由式(4)轉(zhuǎn)換為等失效數(shù)據(jù)：

(4)

1.3 等失效數(shù)據(jù)法2

(5)

2 實(shí)例分析

2.1 數(shù)據(jù)背景

筆者所用數(shù)據(jù)引自文獻(xiàn)[4]，如表1所示，表示在3種不同氣壓z下壓力表的壽命時(shí)間t，包括元件失效替換時(shí)間和預(yù)防替換時(shí)間(帶“+”號(hào)值)，從表1可以看出，氣壓越高元件壽命越短。

表1 相應(yīng)氣壓下壓力表的壽命

在進(jìn)行可靠性分析建模時(shí)，文獻(xiàn)[4-5]均假定預(yù)防替換恰好是在元件即將要失效時(shí)進(jìn)行的，因此都是把預(yù)防替換的截尾數(shù)據(jù)當(dāng)成失效數(shù)據(jù)處理，顯然這將低估元件的壽命，筆者將使用等失效數(shù)據(jù)法重新分析這些數(shù)據(jù)。

2.2 分布模型

表1的數(shù)據(jù)中包含協(xié)變量z，令zi=2+i(i=1,2,3)表示3種氣壓，Ti為氣壓zi下的元件壽命。為了方便分析，假定T1,k1T2,k2T3為獨(dú)立同分布的隨機(jī)變量，其中k1,k2為待估“加速因子”，滿(mǎn)足1

(1)Ti均服從同一類(lèi)分布。

(2)Ti各分布函數(shù)中形狀參數(shù)一致，尺度參數(shù)不一致。

(3)對(duì)于等失效數(shù)據(jù)法1，令δ=δ(z1)，有：

δ(zi)=δ(z1)/ki-1,1≤i≤3,k0=1

(6)

為了確定元件合理的壽命分布，視類(lèi)截尾數(shù)據(jù)為失效數(shù)據(jù)，則壓力表失效數(shù)據(jù)的對(duì)數(shù)似然函數(shù)可表示為：

(7)

在元件可靠性分析中，威布爾分布、伽馬分布、正態(tài)分布和對(duì)數(shù)正態(tài)分布應(yīng)用十分廣泛，故選擇它們作為備選的壽命分布模型[6]。根據(jù)各分布函數(shù)的特性可知，威布爾分布和伽馬分布的形狀參數(shù)β獨(dú)立于z，但是尺度參數(shù)η的大小與z有關(guān)；正態(tài)分布中ρ=μ/σ獨(dú)立于z，但μ的大小與z有關(guān)；對(duì)數(shù)正態(tài)分布的σl獨(dú)立于z，但μl的大小與z有關(guān)。因此，對(duì)于4類(lèi)分布模型，式(7)包含了4個(gè)獨(dú)立參數(shù)。

通過(guò)Excel的“規(guī)劃求解”工具使ln(L)值最大，可得到各分布模型的參數(shù)值，如表2所示，根據(jù)似然函數(shù)值ln(L)的大小，可知對(duì)數(shù)正態(tài)分布模型是最佳的擬合模型。

表2 各模型的分布參數(shù)

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

設(shè)隨機(jī)變量t服從對(duì)數(shù)正態(tài)分布，其對(duì)數(shù)均值為μl，對(duì)數(shù)標(biāo)準(zhǔn)差為σl，則對(duì)數(shù)正態(tài)分布的累積分布函數(shù)可表示為：

F(t)=Φ[ln(t);μl,σl]

(8)

令：

(9)

由式(9)可以看出，1/σl類(lèi)似于威布爾分布的形狀參數(shù)β，eμl類(lèi)似于威布爾分布的尺度參數(shù)η。因此，假定1/σl和μl分別為對(duì)數(shù)正態(tài)分布的形狀參數(shù)和尺度參數(shù)。則對(duì)數(shù)正態(tài)分布的均值和標(biāo)準(zhǔn)偏差分別可表示為[7]：

(10)

其平均剩余壽命為：

(11)

式中：Φ(w;σl,1)為平均值為σl，標(biāo)準(zhǔn)差為1的正態(tài)分布函數(shù)；R(t)為可靠性函數(shù)。

由于氣壓zi下的平均壽命μi(1≤i≤3)與“加速因子”的關(guān)系為μ1/μ2=k1，μ1/μ3=k2，則對(duì)數(shù)正態(tài)分布暗含如下關(guān)系：

μl,2=μl-ln(k1),μl,3=μl-ln(k2)

(12)

2.4 兩種特殊情況的分析

δj包含兩種特殊情況，即把類(lèi)截尾數(shù)據(jù)當(dāng)成失效數(shù)據(jù)或當(dāng)成隨機(jī)截尾數(shù)據(jù)。兩種情況均可用極大似然法估計(jì)參數(shù)。情況(1)可根據(jù)式(7)直接由密度函數(shù)得到似然函數(shù)，而情況(2)中截尾數(shù)據(jù)的似然函數(shù)則是可靠性函數(shù)，即：

(13)

兩種情況下的參數(shù)特性值如表3所示，其中μi(1≤i≤3)為氣壓zi下壓力表的平均壽命，兩種結(jié)果對(duì)比可知兩者形狀參數(shù)一致，只是平均壽命不同，但相差不大，這說(shuō)明該實(shí)例樣本的截尾度不大。

表3 關(guān)于函數(shù)z的參數(shù)特性值

2.5 等失效數(shù)據(jù)法分析

根據(jù)兩種等失效數(shù)據(jù)法可求得參數(shù)特性值如表4所示。其中pi(1≤i≤3)為氣壓zi下RRL與MRL的比值，δ1為氣壓z1下的RRL。方法1中δ1是中間變量，方法2中p=pi是中間變量。

表4 等失效數(shù)據(jù)法的參數(shù)特性值

由表4結(jié)果可得出如下結(jié)論：

(1)兩種方法均滿(mǎn)足pi∈(0,1)，表明存在最優(yōu)中間變量，且其提供了關(guān)于RRL的信息。

(2)由對(duì)數(shù)似然值ln(L)的大小可知等失效數(shù)據(jù)方法2優(yōu)于方法1。

(3)兩種方法的結(jié)果相吻合。

2.6 兩種參數(shù)估計(jì)法的對(duì)比分析

當(dāng)截尾數(shù)據(jù)被視為隨機(jī)截尾數(shù)據(jù)時(shí)，可用兩種方法進(jìn)行參數(shù)估計(jì)，第一種是當(dāng)p=1時(shí)，求式(2)的極大對(duì)數(shù)似然函數(shù)，即EM算法，第二種是極大似然法(maximum likelihood method,ML算法)，此時(shí)截尾數(shù)據(jù)子集的似然函數(shù)是可靠性函數(shù)，見(jiàn)式(13)。表5為兩種方法所得參數(shù)，結(jié)果表明，形狀參數(shù)完全不同，且ML算法所得壽命估值大于EM算法，這可能是對(duì)數(shù)正態(tài)分布有重右尾性質(zhì)的緣故。

表5 兩種估計(jì)方法的參數(shù)特性值

圖1為3種情況下的密度函數(shù)曲線(xiàn)，曲線(xiàn)1為類(lèi)截尾數(shù)據(jù)被視為失效數(shù)據(jù);曲線(xiàn)2為類(lèi)截尾數(shù)據(jù)被當(dāng)成隨機(jī)截尾數(shù)據(jù)，且參數(shù)估計(jì)是用極大似然法;曲線(xiàn)3為根據(jù)等失效數(shù)據(jù)法2,當(dāng)p=0.446 5時(shí)所得。

圖1 3種情形下的密度函數(shù)曲線(xiàn)

由圖1可知，等失效數(shù)據(jù)法所得概率密度函數(shù)曲線(xiàn)居于兩種特殊情況之間。

3 結(jié)論

筆者定義了類(lèi)截尾數(shù)據(jù)，引入了代表的剩余壽命RRL來(lái)處理可靠性建模中的類(lèi)截尾數(shù)據(jù)，并提出了兩種確定RRL的等失效數(shù)據(jù)法，一種是以一個(gè)固定值作為RRL，另一種是以與平均剩余壽命MRL成一定比例的值作為RRL。筆者通過(guò)實(shí)例分析了這兩種方法，并由對(duì)數(shù)似然值可知方法2優(yōu)于方法1，兩種方法均引入一個(gè)中間變量表示數(shù)據(jù)隱含的狀態(tài)信息，實(shí)例證實(shí)了存在最優(yōu)中間變量，且通過(guò)分析發(fā)現(xiàn)：①由預(yù)防維修產(chǎn)生的截尾數(shù)據(jù)被視為失效數(shù)據(jù)或隨機(jī)截尾數(shù)據(jù)時(shí)，將低估或高估產(chǎn)品壽命。②由EM算法所得的參數(shù)估計(jì)結(jié)果可能與ML算法所得結(jié)果不同。

[1] 蔣仁言.威布爾模型族:特性、參數(shù)估計(jì)和應(yīng)用[M].北京：科學(xué)出版社，1998: 8-20.

[2] DEMPSTER A P, LAIRD N M, RUBIN D B. Maximum likelihood from incomplete data via the EM algorithm [J]. J R Stat Soc, 1977,39(B1):1-37.

[3] JIANG R, JARDINE A K S. Composite scale modeling in the presence of censored data[J]. Reliability Engineering and System Safety,2006(91):756-764.

[4] LOVE C E, GUO R. Using proportional hazard modeling in plant maintenance [J]. Quality & Reliability Engineering International, 1991(7):7-17.

[5] JIANG R, JARDINE A K S. Health state evaluation of an item: a general framework and graphical representation[J]. Reliability Engineering and System Safety, 2008(93):89-99.

[6] JIANG R. A trade off BXlife and its applications[J]. Reliability Engineering and System Safety, 2013(113):1-6.

[7] WANG W B. An overview of the recent advances in delay-time-based maintenance modeling [J]. Reliability Engineering & System Safety,2012(106):165-178.

TIAN Ying:Experimentalist; School of Mechanical & Electrical Engineering, Foshan University, Foshan 528000, China.

[編輯：王志全]

Analysis of Censored Data with Partial Condition Information

TIANYing,SHEYang

An equivalent failure data approach was presented to analyze the quasi-censored data which implicitly contain the condition information in reliability modeling. The approach transformed a censored observation into an equivalent failure observation by adding a representative residual life to the censored observation. The representative residual life can be a fixed quantity or proportional to the mean residual life. In each case, an intermediate variable which represents the condition information was introduced. The intermediate variable and distribution parameters were simultaneously estimated using the maximum likelihood method. The effectiveness of the approach was illustrated by a real-world example.

Quasi-censored data; partial condition information; representative residual life; preventive maintenance; intermediate variable

2015-03-10.

田英(1986-)，女，湖南婁底人，佛山科學(xué)技術(shù)學(xué)院機(jī)械與電氣工程學(xué)院實(shí)驗(yàn)師.

2095-3852(2015)05-0598-04

A

O211.9；TB114.3

10.3963/j.issn.2095-3852.2015.05.016