缺失數(shù)據(jù)場合下Frechet分布參數(shù)的逆矩估計

何 亮,徐曉嶺,吳生榮

(1.上海對外經(jīng)貿(mào)大學 統(tǒng)計與信息學院， 上海 201620； 2.宜興出入境檢驗檢疫局， 江蘇 宜興 214206)

在極值分布理論中，一般有三種類型的分布函數(shù)，即Gumbel分布、Frechet分布以及Weibull分布[1]。1927年Frechet[2]中推導出在初值非負的前提下最大極值分布的漸近分布表達式。1984年Mann[3]給出了Frechet分布與Gumbel分布之間的關系及其參數(shù)估計。Frechet分布在諸多領域有著廣泛的應用，Longin[4]用Frechet分布擬合股市的極端價格變動走勢，Harlow在文獻[5]中說明Frechet分布經(jīng)常被用在工程應用材料屬性的統(tǒng)計模型中，Nadarajah[6]將Frechet分布用于社會學模型。Zaharim等[7]將對數(shù)正態(tài)分布以及Frechet分布用于風速數(shù)據(jù)的預測分析中。Mubarak[8]討論了Frechet分布在逐步增加的定數(shù)截尾數(shù)據(jù)并包含隨機移除場合下參數(shù)的最大似然估計以及區(qū)間估計。Abbas和湯銀才[9]研究了在定數(shù)截尾數(shù)據(jù)場合下參數(shù)的最大似然估計以及最小二乘估計。Sindhu等[10]分別研究了定時截尾以及定數(shù)截尾數(shù)據(jù)場合下，F(xiàn)rechet分布形狀參數(shù)已知情況下尺度參數(shù)的Bayes估計。曾林蕊等[11]給出了Frechet分布在定數(shù)截尾數(shù)據(jù)缺失場合下參數(shù)的近似極大似然估計。

本文首先研究了兩參數(shù)Frechet分布的失效率函數(shù)與平均失效率函數(shù)，說明其圖像均為“倒浴盆”曲線。給出了Frechet分布在缺失數(shù)據(jù)場合下兩種逆矩估計，利用Monte-Carlo方法與文獻[11]中的近似極大似然估計進行比較，說明了逆矩估計在樣本量較小以及缺失數(shù)據(jù)較多情況下在均方誤差意義上優(yōu)于近似極大似然估計。

1 Frechet分布及參數(shù)的逆矩估計

1.1 Frechet分布的數(shù)學特征

設產(chǎn)品壽命T服從兩參數(shù)Frechet分布，記為T～Fr(α,β)，其分布函數(shù)及密度函數(shù)分別為：

其中α>0為形狀參數(shù)，β>0為尺度參數(shù)。

1) 當α>1時，T的數(shù)學期望為：

2) 當α>2時，

方差為：

變異系數(shù)為：

3) 當α>3時，偏度為：

[Γ(1-3/α)- 3Γ(1-1/α)Γ(1-2/α)+2Γ3(1-1/α)]

4) 當α>4時，峰度為：

4Γ(1-1/α)Γ(1-3/α)+3Γ2(1-2/α)]-6

Frechet分布變異系數(shù)、偏度及峰度函數(shù)關于α單調減少，如圖1所示。

引理1：記非負隨機變量T的密度函數(shù)、分布函數(shù)及失效率函數(shù)分別為f(t)，F(xiàn)(t)，λ(t)，密度函數(shù)在(0,+∞)上連續(xù)并二階可導。記g(t)=[1-F(t)]/f(t)，η(t)=-f′(t)/f(t)，則有結論：若對所有t有η′(t)>0，則λ(t)單調遞增。若對所有t有η′(t)<0，則λ(t)單調遞減。若η(t)先減后增，即η(t)為“浴盆”曲線，則：

如果存在y0使得g′(y0)=0，則η(t)也為“浴盆”曲線，否則η(t)單調遞增。若η(t)先增后減，即η(t)為“倒浴盆”曲線，則：如果存在y0使得g′(y0)=0，則η(t)也為“倒浴盆”曲線，否則η(t)單調遞減。

證明：T的失效率函數(shù)為

由引理1，有：

故當00；當t>βα1/α時，η′(t)<0。由引理1可知失效率函數(shù)要么為“倒浴盆”曲線，要么單調遞減。又

因此存在y0∈(0,+∞)，使得λ(y0)>0，故T的失效率函數(shù)為“倒浴盆”曲線。

證明：T的平均失效率函數(shù)為

則有：

令u=β/t，則G(u)=-u/β·ln(1-e-uα)，記h(u)=uln(1-e-uα)，則：

Frechet分布失效率函數(shù)和平均失效率函數(shù)圖像如圖2所示,在不同的參數(shù)下均呈現(xiàn)“倒浴盆”形狀。

1.2 缺失數(shù)據(jù)壽命分布參數(shù)的逆矩估計

在產(chǎn)品壽命試驗中，由于試驗設備故障以及數(shù)據(jù)采集過程遺漏等原因，比如數(shù)據(jù)存儲的失敗,存儲器的損壞，數(shù)據(jù)錄入人員失誤漏錄等，均會造成試驗數(shù)據(jù)的缺失。較之全樣本數(shù)據(jù)場合下的統(tǒng)計分析，缺失數(shù)據(jù)只能得到其前后的壽命數(shù)據(jù)，其統(tǒng)計方法有所不同。缺失數(shù)據(jù)情形在現(xiàn)實中經(jīng)常發(fā)生，因此對于缺失數(shù)據(jù)的研究是非常必要的。

假設有n個產(chǎn)品進行壽命服從Fr(α,β)分布的壽命試驗，其次序失效數(shù)據(jù)為：t(1)≤t(2)≤…≤t(n)?，F(xiàn)由于某些原因造成部分數(shù)據(jù)缺失，剩余的k個數(shù)據(jù)記為0t(r0)

1.2.1 逆矩估計I

證明：首先有

要證明ln(Sk/Si)單調，則證明函數(shù)δ(α)單調即可。對δ(α)求導：

由于T(rk+1-j2)0，δ(α)關于α單調遞增。從而ln(Sk/Si)，i=1,2,…,k-1為關于α的單調函數(shù)。

Q1(α)=E1

由引理2可知上述方程有唯一正實根。

針對不同的樣本數(shù)n，在不同的r1,r2,…,rk下，進行10 000次的Monte-Carlo模擬，將Q1,Q2的分位數(shù)、期望及方差匯總至表1。

表1 Q1，Q2的分位數(shù)

1.2.2 逆矩估計II

引理3：若隨機變量T～Fr(α,β)，

T(1),T(2),…,T(n)為其前n個次序統(tǒng)計量。序列{ri}，i=1,2,…,k滿足0

證明：對Q3(α)求導：

所以對于1

Q3(α)=E3

根據(jù)引理3可知上述方程的解唯一。

針對不同的樣本數(shù)n，在不同的r1,r2,…,rk下，進行 10 000次的Monte-Carlo模擬，將Q3,Q4的分位數(shù)、期望及方差匯總至表2。

表2 Q3，Q4的分位數(shù)

2 Monte-Carlo模擬

為考察參數(shù)估計值的精確性，在不同的數(shù)據(jù)場合下進行模擬試驗。在真值為α=0.5,1.0,1.5，β=1，樣本數(shù)為n=10,20，不同的缺失數(shù)據(jù)情形下分別進行5 000次Monte-Carlo模擬試驗。對比逆矩估計I、逆矩估計II以及文獻[11]中的近似極大似然估計方法估計的結果,模擬結果見表3，括號中的數(shù)值表示估計值的均方誤差。從表3中可以看出，當樣本數(shù)為n=20,k=10時，AMLE的結果與逆矩估計I、逆矩估計II的精度近似；當n=20,k=6，即缺失數(shù)據(jù)較多時，AMLE的點估計誤差較大，而且逆矩估計I、逆矩估計II的結果比較滿意；當n=10,k=5，即樣本數(shù)較少、缺失數(shù)較多時，兩種逆矩估計的均方誤差均低于AMLE的均方誤差。

對于逆矩估計I及逆矩陣估計II，為考察其區(qū)間估計的精度，取真值α=β=1.0，樣本數(shù)為n=10,20，置信水平為0.90，分別在不同的缺失數(shù)據(jù)情形下進行10 000次Monte-Carlo模擬試驗，得到參數(shù)的區(qū)間估計，并記錄真值落入?yún)^(qū)間的次數(shù)以及區(qū)間的平均長度，結果如表4所示。從表4中結果可以看出：逆矩估計I的平均置信區(qū)間長度均小于逆矩估計II的平均置信區(qū)間長度，且真值落入?yún)^(qū)間的次數(shù)多于逆矩估計II。結合引理3對方程條件有一定限制，故實際應用中應考慮使用逆矩估計I。

表3 Monte-Carlo模擬對比結果

表4 Monte-Carlo模擬的區(qū)間估計(置信水平90%下)

算例：取真值為α=β=1.0，樣本容量為n=20，ri取2,5,6,8,13,15,16,17,18,20。通過Monte-Carlo方法產(chǎn)生10個服從Frechet的隨機數(shù)：t(2)=0.391 3，t(5)=0.494 0，t(6)=0.768 9，t(8)=1.143 5，t(13)=3.738 7，t(15)=5.351 6，t(16)=6.443 3，t(17)=8.144 9，t(18)=12.652 8，t(20)=65.408 0。

3 結論

對兩參數(shù)Frechet分布在缺失數(shù)據(jù)場合下提出運用逆矩估計的方法，提出兩種不同的逆矩估計方法，對參數(shù)進行了點估計以及區(qū)間估計，通過Monte-Carlo模擬試驗以及算例證明了逆矩估計在均方誤差意義上要優(yōu)于近似極大似然估計，并將兩種逆矩估計的方法進行對比，說明逆矩估計I更優(yōu)。