何朝兵
(安陽師范學(xué)院 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院 河南 安陽 455000)
右刪失左截?cái)鄶?shù)據(jù)下離散威布爾分布的參數(shù)估計(jì)
何朝兵
(安陽師范學(xué)院 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院 河南 安陽 455000)
研究了右刪失左截?cái)鄶?shù)據(jù)模型下離散威布爾分布參數(shù)的極大似然估計(jì)和漸近置信區(qū)間.介紹了參數(shù)估計(jì)的牛頓迭代方法和EM算法,給出了參數(shù)的漸近置信區(qū)間.隨機(jī)模擬的結(jié)果表明,牛頓迭代方法和EM算法得到的參數(shù)估計(jì)結(jié)果差別不大.
極大似然估計(jì); 牛頓迭代方法; EM算法; 缺損信息原則; 漸近置信區(qū)間
離散威布爾分布是一種很重要的離散型壽命分布[1-3],它是威布爾分布的離散化,在排隊(duì)論和可靠性數(shù)學(xué)等分支中有著很廣泛的應(yīng)用,并且?guī)缀畏植际翘厥獾碾x散威布爾分布.當(dāng)觀察壽命數(shù)據(jù)時(shí),右刪失與左截?cái)嗲樾谓?jīng)常同時(shí)發(fā)生.例如,對(duì)確診為艾滋病的患者進(jìn)行觀察研究,如果患者仍然存活或提前退出研究,這就是右刪失;如果患者在研究開始之前就已經(jīng)死亡,這就是左截?cái)?右刪失左截?cái)鄶?shù)據(jù)模型廣泛應(yīng)用于醫(yī)學(xué)、生物學(xué)、經(jīng)濟(jì)學(xué)等領(lǐng)域,近年來對(duì)此模型的研究比較熱.文獻(xiàn)[4-8]研究了此模型下壽命分布為連續(xù)型的情形,但對(duì)離散威布爾分布情形的研究尚未見報(bào)道.本文主要利用牛頓迭代(NR)方法和EM算法研究了右刪失左截?cái)鄶?shù)據(jù)模型下離散威布爾分布參數(shù)的極大似然估計(jì)和區(qū)間估計(jì),隨機(jī)模擬的結(jié)果表明,由這兩種方法得到的參數(shù)估計(jì)結(jié)果差別不大.
設(shè)(X,Y,T)是離散型隨機(jī)變量,X的分布函數(shù)為F(x,θ)=P(Xi≤x),分布律為f(x,θ),其中θ是向量參數(shù);Y是右刪失隨機(jī)變量,分布函數(shù)為G(y),分布律為g(y);T是左截?cái)嚯S機(jī)變量,分布函數(shù)為H(t),分布律為h(t),且Y,T的分布與參數(shù)θ無關(guān).假定X,Y,T是相互獨(dú)立且取正整數(shù)的隨機(jī)變量,X是所感興趣的隨機(jī)變量.右刪失左截?cái)鄶?shù)據(jù)的試驗(yàn)?zāi)P褪牵簝H在Zi≥Ti時(shí)得到觀察數(shù)據(jù)(Zi,Ti,δi),而在Zi 求以下樣本的似然函數(shù): P(無樣本觀察值)=P(Zi 為了敘述與書寫方便,假定前n1個(gè)樣本有觀察值,剩下的n2個(gè)樣本沒有觀察值(n1+n2=n),則基于數(shù)據(jù){(Zi,Ti,δi),1≤i≤n1}的似然函數(shù)為 若X的分布函數(shù)為F(x;α,β)=1-e-(x/β)α,x>0,且α>0,β>0,則稱X服從尺度參數(shù)為β,形狀參數(shù)為α的威布爾分布,記為X~Wei(α,β).若X的分布律為P(X=k)=q(k-1)α-qkα,k=1,2,…,且0 假設(shè)右刪失左截?cái)鄶?shù)據(jù)試驗(yàn)中,所感興趣的變量X服從形狀參數(shù)為α,尺度參數(shù)為q的離散威布爾分布,此時(shí)θ=(α,q),下面介紹參數(shù)α與q的估計(jì)方法. 2.1 牛頓迭代(NR)方法 基于數(shù)據(jù){(Zi,Ti,δi),1≤i≤n1}的似然函數(shù)為 其對(duì)數(shù)似然函數(shù)為 令其兩個(gè)偏導(dǎo)數(shù)為零,得到的兩個(gè)似然方程為超越方程,只能用數(shù)值方法求解,本文利用NR方法求解. NR方法是通過使似然函數(shù)最大化而獲得極大似然估計(jì)的一種直接方法,它需要計(jì)算對(duì)數(shù)似然函數(shù)關(guān)于參數(shù)的一階和二階偏導(dǎo)數(shù),所以十分煩瑣.這里,通過使用R軟件maxLik程序包中的maxNR函數(shù)來獲得q和α的極大似然估計(jì). 注 如果Y服從幾何分布Geo(p0),T服從取值1,2,…,s的離散均勻分布,則 2.2 EM算法 EM算法是一種迭代方法,最初由文獻(xiàn)[9]提出,主要用來求后驗(yàn)分布的眾數(shù).它的每一步迭代由兩步組成:E步(求期望)和M步(極大化).EM算法處理不完全數(shù)據(jù)非常方便[10-11],下面用EM算法來求α和q的極大似然估計(jì). 若第i個(gè)樣本沒有觀察值,添加其觀察值為(Wi,βi),其中:Wi=Xi∧Yi=min(Xi,Yi),βi=I(Xi≤Yi),i=n1+1,n1+2,…,n,則 可得似然函數(shù) 取(α,q)的先驗(yàn)分布為無信息先驗(yàn)分布π(α,q)=c,0 E步: 在給定θ,δ,Z和T下,(Wi,βi)的分布律為 則(Wi,αi)關(guān)于Wi的邊緣分布律為 P(Wi=k)=ψ1(k,θ)+ψ2(k,θ)ψ3(k,θ),k=1,2,…. 所以 注 如果Y服從幾何分布Geo(p0),T服從取值1,2,…,s的離散均勻分布,則 2.3 極大似然估計(jì)的漸近方差和協(xié)方差 利用NR方法求極大似然估計(jì)時(shí),I0可以直接得到.而利用EM算法時(shí),由于似然函數(shù)為基于完全數(shù)據(jù)的似然函數(shù),所以I0不能直接得到,但根據(jù)缺損信息原則可以獲得觀察信息陣I0.缺損信息原則為:觀察信息=完全信息-缺損信息.下面求基于EM算法的觀察信息,完全數(shù)據(jù)的似然函數(shù)為 為了方便書寫,記 φ2=n2{βln[q(w-1)α-wα-1]+wαlnq}, ψ(β,w)βψ1(w)+(1-β)ψ2(w),β=0,1, lnL2n2lnψ(β,w). 則完全信息為 (1) 損失信息為 (2) 根據(jù)缺損信息原則,可得觀察信息為 2.4 置信區(qū)間 式中:bq和σq分別為q的1 000個(gè)估計(jì)值的bootstrap偏差和方差;zβ/2為標(biāo)準(zhǔn)正態(tài)分布的上β/2分位點(diǎn).α的參數(shù)bootstrap置信區(qū)間的構(gòu)造與q的類似. 表1 NR方法和EM算法下參數(shù)估計(jì)的均值(M)、偏差(B)、均方誤差(MSE)和置信區(qū)間的覆蓋率(CP) 從表1可以看出,通過NR方法和EM算法得到的參數(shù)估計(jì)的均值、偏差、均方誤差以及置信水平為0.95的置信區(qū)間的覆蓋率都比較接近,說明這兩種方法差別不大.同時(shí)可以看出,樣本容量對(duì)估計(jì)值的影響也不大,說明得到的估計(jì)值是比較穩(wěn)定的,并且精度也較高. 表2 由3種方法構(gòu)造的置信區(qū)間 由表2可以看出,通過NR方法和EM算法構(gòu)造的置信區(qū)間是一樣的,雖然它們與參數(shù)bootstrap置信區(qū)間不一樣,但差別不是很大. [1] NEKOUKHOU V, BIDRAM H. The exponentiated discrete Weibull distribution[J]. SORT, 2015, 39(1): 127-146. [2] ALMALKI S J, NADARAJAH S. Modifications of the Weibull distribution: a review[J]. Reliab Eng Syst Safe, 2014, 124: 32-55. [3] ENGLEHARDT J D, LI R. The discrete Weibull distribution: an alternative for correlated counts with confirmation for microbial counts in water[J]. Risk Anal, 2011, 31(3): 370-381. [4] DEWAN I. Comments: EM-based likelihood inference for some lifetime distributions based on left truncated and right censored data and associated model discrimination[J]. S Afr Stat J, 2014, 48(2): 183-185. [5] SHEN P S. Aalen’s additive risk model for left-truncated and right-censored data[J]. Commun Stat Simul Comput, 2014, 43(5): 1006-1019. [6] BALAKRISHNAN N, MITRA D. Likelihood inference based on left truncated and right censored data from a gamma distribution[J]. IEEE Trans Rel, 2013, 62(3): 679-688. [7] SU Y R, WANG J L. Modeling left-truncated and right-censored survival data with longitudinal covariates[J]. Ann Stat, 2012, 40(3): 1465-1488. [8] AHMADI J, DOOSTPARAST M, PARSIAN A. Estimation with left-truncated and right censored data: a comparison study[J]. Stat Probab Lett, 2012, 82(7): 1391-1400. [9] 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(1): 1-38. [10] CHUNG Y, LINDSAY B G. Convergence of the EM algorithm for continuous mixing distributions[J]. Stat Probab Lett, 2015, 96(1): 190-195. [11]GRIGOROVA D, ENCHEVA E, GUEORGUIEVA R. EM algorithm for MLE of a probit model for multiple ordinal outcomes[J]. Serdica J Comput, 2013, 7(3): 227-244. (責(zé)任編輯:孔 薇) Parameter Estimations of Discrete Weibull Distribution with Left Truncated and Right Censored Data HE Chaobing (SchoolofMathematicsandStatistics,AnyangNormalUniversity,Anyang455000,China) The maximum likelihood estimation and asymptotic confidence intervals of the parameters of discrete Weibull distribution were mainly studied with left truncated and right censored data. Newton-Raphson method and EM algorithm of parameter estimation were introduced, and the asymptotic confidence intervals of the parameters were given. Random simulation test results showed that there was little difference on parameter estimation between Newton-Raphson method and EM algorithm. maximum likelihood estimation; Newton-Raphson method; EM algorithm; missing information principle; asymptotic confidence interval 2015-10-18 國家自然科學(xué)基金資助項(xiàng)目(61174099);河南省高等學(xué)校重點(diǎn)科研項(xiàng)目(16A110001). 何朝兵(1975—),男,河南周口人,講師,碩士,主要從事概率統(tǒng)計(jì)研究,E-mail:chaobing5@163.com. 何朝兵.右刪失左截?cái)鄶?shù)據(jù)下離散威布爾分布的參數(shù)估計(jì)[J]. 鄭州大學(xué)學(xué)報(bào)(理學(xué)版), 2016,48(2): 18-23. O213.2 A 1671-6841(2016)02-0018-06 10.13705/j.issn.1671-6841.20152202 離散威布爾分布的參數(shù)估計(jì)方法
0,則稱X服從尺度參數(shù)為q,形狀參數(shù)為α的離散威布爾分布,當(dāng)α=1時(shí),歸結(jié)為幾何分布Geo(1-q).容易證明:若X~Wei(α,β),且β=(-1/lnq)1/α,則Y=[X]+1服從尺度參數(shù)為q,形狀參數(shù)為α的離散威布爾分布.利用此結(jié)論可以產(chǎn)生離散威布爾分布的隨機(jī)數(shù).
0,c為常數(shù),則θ=(α,q)的添加后驗(yàn)分布為
3 隨機(jī)模擬結(jié)果