何朝兵,劉華文
(1.安陽師范學(xué)院數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,河南安陽455000;2.山東大學(xué)數(shù)學(xué)學(xué)院,山東濟(jì)南250100)
左截?cái)嘤覄h失數(shù)據(jù)下幾何分布參數(shù)的點(diǎn)估計(jì)
何朝兵1,劉華文2
(1.安陽師范學(xué)院數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,河南安陽455000;2.山東大學(xué)數(shù)學(xué)學(xué)院,山東濟(jì)南250100)
證明了左截?cái)嘤覄h失數(shù)據(jù)下幾何分布參數(shù)極大似然估計(jì)的存在唯一性,給出了由EM算法得到的參數(shù)的迭代公式,進(jìn)行了隨機(jī)模擬試驗(yàn),結(jié)果表明參數(shù)的MLE和EM估計(jì)的精度都較高.
左截?cái)嘤覄h失;幾何分布;極大似然估計(jì);EM算法;隨機(jī)模擬
幾何分布是一種很重要的離散型壽命分布,它是指數(shù)分布的離散化,是特殊的離散型威布爾分布,在排隊(duì)論和可靠性數(shù)學(xué)等分支中有著很廣泛的應(yīng)用.文獻(xiàn)[1-6]對一般截尾壽命試驗(yàn)下幾何分布的參數(shù)進(jìn)行了統(tǒng)計(jì)分析.近些年來,對左截?cái)嘤覄h失數(shù)據(jù)的研究比較多,文獻(xiàn)[7-13]深入研究了被觀察變量是連續(xù)型分布的情況,但還沒有文獻(xiàn)涉及幾何分布的情況.本文首先證明了左截?cái)嘤覄h失數(shù)據(jù)下幾何分布參數(shù)極大似然估計(jì)的存在唯一性,然后給出了由EM算法得到的參數(shù)的迭代公式,最后進(jìn)行了隨機(jī)模擬試驗(yàn),結(jié)果表明參數(shù)的MLE和EM估計(jì)的精度都較高.
設(shè)(X,Y,T)是一離散型隨機(jī)變量,X的分布函數(shù)為F(x,p)=P(X≤x),分布律為f(x,p).這里:p是參數(shù);Y是一右刪失隨機(jī)變量,分布函數(shù)為G(y),分布律為g(y);T是一左截?cái)嚯S機(jī)變量,分布函數(shù)為H(t),分布律為h(t),且Y,T的分布與參數(shù)p無關(guān).假定X,Y,T是相互獨(dú)立取正整數(shù)的隨機(jī)變量,X是我們感興趣的隨機(jī)變量.左截?cái)嘤覄h失數(shù)據(jù)的試驗(yàn)?zāi)P褪牵簝H在Zi≥Ti時得到觀察數(shù)據(jù)(Zi,Ti,δi),而在Zi<Ti下無法得到任何觀察值,其中:下面求樣本的似然函數(shù):
其中
為了敘述與書寫方便,我們假定前n1個樣本有觀察值,剩下的n2個樣本沒有觀察值(n1+n2=n).則基于數(shù)據(jù){(Zi,Ti,δi),1≤i≤n1}的似然函數(shù)為
其中
且A與參數(shù)p無關(guān).
當(dāng)Xi服從幾何分布Geo(p)時,我們介紹參數(shù)p點(diǎn)估計(jì)的兩種方法.
2.1 極大似然估計(jì)法
基于數(shù)據(jù){(Zi,Ti,δi),1≤i≤n1}的似然函數(shù)為
其中
由于
所以φ′(p)<0.則φ(p)是p的嚴(yán)格減函數(shù),并且
由于u(p)與刪失和截?cái)嘧兞康姆植加嘘P(guān),從而導(dǎo)致u(p)的表達(dá)式一般會很繁瑣,所以^p的顯式表達(dá)式通常很難得到,這時可以通過數(shù)值方法(如二分法)獲得,這個^p就是p的MLE.實(shí)際操作中可利用R軟件中的uniroot( )函數(shù)得到^p.
注1 如果Y服從Geo(p1),T服從Geo(p2),則
把u(p)和u′(p)帶入(1)式,可得到對數(shù)似然方程
2.2 EM算法
由于左截?cái)嘤覄h失數(shù)據(jù)下幾何分布的似然函數(shù)比較復(fù)雜,下面添加部分缺損的Xi的值以便獲得較簡單的似然函數(shù).
若第i個樣本有觀察值且δi=0,只知道Xi>zi,雖然Xi的值有缺損,但由于ˉF(zi,p)=qzi比較簡單,所以此Xi的值就不添加了.
若第i個樣本沒有觀察值,添加其觀察值為(Wi,αi),其中:
則
可得似然函數(shù)
對于p如果沒有任何先驗(yàn)信息可利用時,可取π(p)為(0,1)上的均勻分布;對于p如果有先驗(yàn)信息時,可取共軛先驗(yàn)分布貝塔分布Be(b,c),即π(p)∝pb-1qc-1,0<p<1,b>0,c>0.實(shí)際上(0,1)上的均勻分布是Be(1,1).
下面假設(shè)對于p沒有任何先驗(yàn)信息可利用,取π(p)為(0,1)上的均勻分布,則p的添加后驗(yàn)分布為
在第m+1次迭代中,假設(shè)有估計(jì)值p(m),則可通過E步和M步得到p的一個新的估計(jì).為了書寫方便,簡記(|p(m),δ,Z,T)為(|·).
E步
顯然,在p(m),δ,Z和T給定下,
故
在p(m),δ,Z和T給定下,(Wi,αi)的分布律為:
則(Wi,αi)關(guān)于Wi的邊緣分布律為:
則
M步
(3)式給出了由EM算法得到的參數(shù)p的迭代公式.
注2 由于L(p|δ,Z,T)=L(p)≤L(^p),其中^p是MLE,所以L(p|δ,Z,T)是有界函數(shù),由EM算法的性質(zhì)可知[14],上述的EM算法是收斂的.
由ψ(p)和ψ3(k,p)的表達(dá)式可以看出,選擇恰當(dāng)?shù)挠覄h失變量和左截?cái)嘧兞繒笶M迭代公式變得簡潔.若Y服從Geo(p1),T服從Geo(p2),下面求具體的迭代公式.
再根據(jù)(3)式,可得EM迭代公式
其中:
基于上面的討論,下面進(jìn)行隨機(jī)模擬試驗(yàn).
設(shè)Xi服從Geo(0.4),Yi服從Geo(0.6),Ti服從Geo(0.9),樣本容量分別取n=30,50,100,200,300,500,800.為了準(zhǔn)確性,對每一固定樣本容量產(chǎn)生1 000個隨機(jī)樣本,由每個樣本計(jì)算出一個估計(jì)值,然后取這1 000個估計(jì)值的算術(shù)平均值作為這一樣本容量下參數(shù)的估計(jì)值.運(yùn)用EM算法時從p(0)=0.5開始迭代.
隨機(jī)模擬的主要理論依據(jù)是(2)式和(4)式,編寫R程序時用到的函數(shù)主要有rgeom( ),min( ),uniroot( ).
隨機(jī)模擬結(jié)果見表1.
表1 左截?cái)嘤覄h失數(shù)據(jù)下幾何分布參數(shù)估計(jì)的隨機(jī)模擬結(jié)果
由表1可以看出,p的MLE和EM估計(jì)的差別不大,與真值0.4的偏差都不超過1%;樣本容量對估計(jì)值的影響也不大.說明得到的估計(jì)值是比較穩(wěn)定的,并且精度也較高.在運(yùn)用EM算法進(jìn)行模擬時發(fā)現(xiàn),迭代6次左右即可收斂,收斂速度很快,并且實(shí)際操作很方便.相比之下,極大似然法在解對數(shù)似然方程時的數(shù)值計(jì)算速度就顯得稍慢些.
[1] BHOJ,DINESHS,ABSANULLAH M.Estimation of the generalized geometric distribution using ranked set sampling[J].Biometrics,1996(52):685-694.
[2] FERGUSON T S.A characterization of the geometric distribution[J].Amer Math Mothly,1972,27(2):256-260.
[3] 徐曉嶺,費(fèi)鶴良,王蓉華.幾何分布的兩個統(tǒng)計(jì)特征[J].應(yīng)用概率統(tǒng)計(jì),2006,22(1):10-20.
[4] 魏立力,張文修.幾何分布的一類貝葉斯停止判決法則[J].應(yīng)用數(shù)學(xué)學(xué)報,2003,26(3):181-185.
[5] 劉銀萍.截?cái)嗲樾蜗聨缀畏植嫉膮?shù)估計(jì)[J].東北師大學(xué)報:自然科學(xué)版,2009,41(3):14-16.
[6] 徐曉嶺,王蓉華,費(fèi)鶴良.幾何分布產(chǎn)品定數(shù)截尾場合下參數(shù)的點(diǎn)估計(jì)[J].強(qiáng)度與環(huán)境,2009,36(2):51-63.
[7] LAI T L,YING Z.Estimating a distribution function with thuncated and censored Data[J].The Annals of Statistics,1991,19:417-442.
[8] GU M G,LAI T L.Functional laws of the iterated logarithm for the product-limit estimatorof a distribution function under random censorship or truncated[J].The Annals of Probability,1990,18:160-189.
[9] ZHOU YONG,SUN LIU-QUAN.Sequential confidence bands for quantile densities under truncated and censored data[J].Acta Mathematicae Applicatae Sinica:English Series,2005,21(2):311-322.
[10] SUN LIU-QUAN.Fixed design nonparametric regression with truncated and censored Data[J].Acta Mathematicae Applicatae Sinica:English Series,2003,19(2):229-238.
[11] 周勇,吳國富.左刪失右截?cái)鄶?shù)據(jù)的分位數(shù)的固定寬度序貫置信區(qū)間估計(jì)[J].應(yīng)用數(shù)學(xué)學(xué)報,2002,25(2):204-216.
[12] 劉煥彬,孫六全.截?cái)嗯c刪失數(shù)據(jù)下的一個回歸方法[J].應(yīng)用數(shù)學(xué)學(xué)報,2005,28(1):1-10.
[13] 茍列紅.左截?cái)嘤覄h失數(shù)據(jù)下半?yún)?shù)模型風(fēng)險率函數(shù)估計(jì)[J].應(yīng)用數(shù)學(xué)學(xué)報,2005,28(4):675-688.
[14] MCLACHLAN G,KRISHNAN T.The EM algorithm and extensions[M].New York:John Wiley &Sons,2007:83.
Point estimation of the parameter of geometric distribution for truncated and censored data
HE Chao-bing1,LIU Hua-wen2
(1.School of Mathematics and Statistics,Anyang Normal University,Anyang 455000,China;2.School of Mathematics,Shandong University,Jinan 250100,China)
In this paper,we firstly prove the existence and uniqueness of MLE of the parameter of geometric distribution for truncated and censored data,then give the iterative formula of the parameter by EM algorithm.Finally random simulation tests are conducted,and the results show that the MLE and EM estimation of the paramete are both fairly accurate.
left truncation and right censorship;geometric distribution;maximum likelihood estimation;EM algorithm;random simulation
O 213.2 [學(xué)科代碼] 110·7120
A
(責(zé)任編輯:陶 理)
1000-1832(2014)02-0025-05
10.11672/dbsdzk2014-02-006
2013-04-17
國家自然科學(xué)基金資助項(xiàng)目(61174099);河南省教育廳自然科學(xué)基金資助項(xiàng)目(2011B110001).
何朝兵(1975—),男,碩士,講師,主要從事概率統(tǒng)計(jì)研究;劉華文(1964—),女,博士,教授,博士研究生導(dǎo)師,主要從事隨機(jī)數(shù)學(xué)研究.