缺失數據場合下Frechet分布參數的逆矩估計

何 亮,徐曉嶺,吳生榮

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

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

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

1 Frechet分布及參數的逆矩估計

1.1 Frechet分布的數學特征

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

其中α>0為形狀參數，β>0為尺度參數。

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

2) 當α>2時，

方差為：

變異系數為：

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分布變異系數、偏度及峰度函數關于α單調減少，如圖1所示。

引理1：記非負隨機變量T的密度函數、分布函數及失效率函數分別為f(t)，F(t)，λ(t)，密度函數在(0,+∞)上連續并二階可導。記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的失效率函數為

由引理1，有：

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

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

證明：T的平均失效率函數為

則有：

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

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

1.2 缺失數據壽命分布參數的逆矩估計

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

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

1.2.1 逆矩估計I

證明：首先有

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

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

Q1(α)=E1

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

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

表1 Q1，Q2的分位數

1.2.2 逆矩估計II

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

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

證明：對Q3(α)求導：

所以對于1

Q3(α)=E3

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

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

表2 Q3，Q4的分位數

2 Monte-Carlo模擬

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

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

表3 Monte-Carlo模擬對比結果

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

算例：取真值為α=β=1.0，樣本容量為n=20，ri取2,5,6,8,13,15,16,17,18,20。通過Monte-Carlo方法產生10個服從Frechet的隨機數：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 結論

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