信息區(qū)間刪失數據的參數估計及敏感性分析

李文靜，鄧文麗，章婷婷

(江西師范大學數學與信息科學學院，江西南昌330022)

0 引言

在臨床實驗或醫(yī)學研究中，由于客觀因素的限制，失效時間常常不能直接觀測到，而只能知道它在某一個區(qū)間內，這類數據在統(tǒng)計學上被稱為區(qū)間刪失數據(interval-censored data).如在一些傳染性疾病感染時間的研究中，實驗對象被放入傳染源后，只能知道在某個觀察點實驗對象是否已染上疾病，染上疾病的具體時間卻無法觀測到，所以只能推測出從接觸傳染源到染上傳染病所經歷的時間落在某個區(qū)間內.

區(qū)間刪失數據存在于許多應用領域中，因此，這引發(fā)了一些統(tǒng)計學者對相關問題的研究.Huang Jian等［1］對區(qū)間刪失數據的分類及對應的統(tǒng)計方法進行了較為詳細地描述.Sun Jian-guo［2］較為全面和系統(tǒng)地概括了區(qū)間刪失數據分析中涉及到的基本概念和方法.呂秋萍等［3］運用無偏轉換思想構造了區(qū)間刪失數據函數的均值估計，并在此基礎上對所構造的估計量方差進行了研究.在區(qū)間刪失數據的研究中，許多學者都是基于失效時間變量T和刪失時間變量C相互獨立的假定進行研究的，稱這種刪失情況為獨立刪失或非信息刪失(Independent Censoring，Noninformative Censoring).然而，在實際問題中，這個假定常常會遭到質疑.如在對某種疾病的治療中，由于病情惡化或者是已接受的治療方案不奏效，從而導致病人退出治療，這種情況通常預示著該病人的存活時間會比較短，即刪失的個體對應的生存時間更短.相反地，有些病人的退出可能因為病情好轉，不需要進一步治療，這種情況的刪失個體的生存時間可能會較長.和獨立刪失相反的是非獨立刪失(Dependent Censoring)，或稱為信息刪失(Informative Censoring).如果對信息刪失數據仍采用獨立刪失下的統(tǒng)計分析方法，則可能會得到有偏或者無效的結論.

在信息刪失數據的研究中，對失效時間和刪失時間相依性的假定是至關重要的.正確的假定可以提高估計的效率，得到更好的統(tǒng)計結論;不合適的假定可能會導致錯誤的結論.在實際應用中，由于造成信息刪失的應用背景和原因的不同，失效時間和刪失時間相依的形式和程度也變得非常復雜，很難準確估計.敏感性分析可以評價相依關系的假定對統(tǒng)計分析結果造成的影響.王純杰［4］基于Copula函數的一些性質，給出了非參數模型下的信息區(qū)間刪失數據分布函數的相合估計.F.Siannis等［5］對失效時間和刪失時間的相依關系進行了假定，引入了標示相關程度的參數和偏度函數，且對參數估計受相依程序的影響進行了敏感性分析.Y.Park等［6］在單個總體和2個總體的情形下，分別對獨立刪失和信息右刪失混合數據下的相關估計問題進行了敏感性分析.Zhang Zhi-gang等［7］在正態(tài)脆弱模型假定下對I型信息區(qū)間刪失數據的比例風險模型進行了敏感性分析.Huang Xue-lin等［8］基于連接函數(Copula)對信息右刪失數據下的比例風險模型的估計問題進行了敏感性分析.

本文擬基于連接函數對信息區(qū)間刪失數據下失效時間的生存函數進行估計，并在3種不同連接函數的情形下關于相依關系對參數估計所造成的影響進行敏感性分析.

1 方法

記Ti為第i個個體的失效時間，Ci為第i個個體的刪失時間.試驗中得到的獨立同分布觀測值為{(ci，δi)，i=1，2，…，n}，其中

假設Ti和Ci的邊際分布函數分別為F(·)和G(·)，g(·)是Ci的邊際密度函數，這里i=1，2，…，n.

基于觀測樣本，可以構造似然函數:

當Ti和Ci相互獨立時，基于觀測樣本的似然函數為

當Ti和Ci不相互獨立時，給定1個有參數α的連接函數H(u，v，α)，假設 Ti和 Ci的聯合分布函數為J(t，c)=P(Ti≤t，Ci≤c)=H(F(t)，G(c)，α)，聯合生存函數為

S(t，c)=P(Ti＞ t，Ci＞ c)=1-F(t)-

G(c)H(F(t)，G(c))，

則第i個個體被刪失的概率為

同理可得，第i個個體失效的概率為

由此可知，

綜上所述，當失效時間和刪失時間不相互獨立時，樣本的觀測似然函數可以表示為

在信息區(qū)間刪失數據中，刪失時間能夠完全觀測到，所以G(·)可以直接用它的經驗分布函數代替，其中

似然函數(1)可以轉化為

如果對失效時間的分布形式掌握的信息不多，則通常會考慮用非參數模型直接估計失效時間的分布函數.類似于文獻［2］給出的獨立區(qū)間刪失數據非參數極大似然估計的方法，可以在(2)式中利用鄧文麗等［9］提出的一類保序最優(yōu)化問題的迭代算法得到分布函數F的估計.當已知一些影響失效時間的協變量時，比例風險模型和加速失效模型是廣泛接受的2類半參數模型.如果假定協變量的影響滿足比例風險模型，則在似然函數的表達式中，邊際分布函數F(·)可以用含回歸系數和基準風險函數的分布函數表達式代替，然后在似然函數(2)中，通過迭代的方法得到相關的估計;如果假定協變量的影響滿足加速失效模型，則在似然函數的表達式中，邊際分布函數F(·)可以用含回歸系數和隨機誤差項的分布函數表達式代替，然后在似然函數(2)中，通過迭代的方法得到相關的估計.張連增等［10］基于極大似然法研究了Copula的參數和半參數方法的估計效果.如果失效時間的分布函數形式已知，而只待估其中包含的參數，則利用似然函數(2)就可以得出參數的極大似然估計.

2 數值模擬

在實際應用中由于失效時間和刪失時間相依的形式和程度非常復雜，很難準確估計，所以通過敏感性分析來評價相依關系的假定對統(tǒng)計分析結果造成的影響.

模擬計算中失效時間T采用威布爾分布隨機生成，因為它的危險率不是常數，所以，與指數分布相比，它有較廣闊的應用，將其用于調查深槽輪滾珠軸承的疲勞壽命，或將其用于描寫電子管的失效.威布爾的分布函數為 F(t)=1-e-(λt)γ，其中，γ 是分布曲線的形狀參數，λ是尺度參數.模擬計算中選取了γ=2，λ=0.5.刪失時間C的邊際分布選取的是(0，A)上的均勻分布，調整A的大小可以改變刪失的比例.

失效時間和刪失時間的相關性選用阿基米德連接函數來描述.這里選取了 Clayton、Gumbe-Hougard和Frank 3種連接函數.

D.G.Clayton［11］給出在 τ=1/(1+2α)下的Copula函數:

E.J.Gumbel等［12］給出 τ=1-1/α 下的Copula函數:

H(u，v，α)=exp{- ［(-log u)α+

(-log v)α］1/α}(α ≥ 1).

2.1.1 30 g/L甲基二磺隆對麥田禾本科雜草的防效經過藥前雜草基數調查可知，試驗藥劑處理區(qū)、對照藥劑處理區(qū)及空白處理區(qū)野燕麥發(fā)生基數為：17.92，18.83，20.42，20.83，20.25，19.83 株；雀麥的發(fā)生基數為：18.92，19.08，16.92，20.00，20.08，20.08 株。

M.J.Frank［13］給出的 Copula 函數:

H(u，v，α)=log{1+(αu-1)(αv-1)

(α -1)}(α ＞ 0，α ≠1)，其中τ=1+4γ-1［D1(γ)-1］，γ=-log α，D1(γ)=

R.B.Nelsen［14］對于連接函數的相關性質和特殊的連接函數進行了詳細介紹.

下面主要是通過數值計算分析連接函數的選取對參數γ和λ的估計產生的影響.這里的穩(wěn)健性分析包括參數敏感性分析和連接函數敏感性分析.

2.1 參數的敏感性分析

分別取 τ=0.8、τ=0.5、τ=0.2的 Frank Copula作為連接函數，T服從γ=2，λ=0.5的威布爾分布，C服從(0，37)上的均勻分布，生成容量為200的樣本，刪失比例P(T＜C)為0.5.在本文方法中，選取 Frank連接函數，τ=0.8.模擬次數為1000，得到上述情況下λ～和γ～的均值、標準差和偏差的估計值(見表1).

表1 總體的參數τ變化下參數γ和λ的估計

由表1可以看出:如果生成樣本的Frank連接函數的參數τ為0.2、0.5、0.8，采用獨立刪失的估計方法，得到的γ和λ估計量的偏差都較大，特別是γ估計值的偏差很大.而采用本文方法(選取Frank連接函數，參數τ=0.8)得到的估計量都比較理想，其估計值的偏差遠遠小于獨立刪失下估計值的偏差.由此可見，在失效時間和刪失時間不相互獨立的情況下，采用獨立刪失方法進行估計會得到不理想的估計結果，因此，應該采用帶相關性假定的模型進行分析.

2.2 改變連接函數下參數估計的敏感性分析

分別取Clayton Copula和Gumbel Copula作為連接函數，τ=0.8，T服從γ=2，λ=0.5的威布爾分布，C服從(0，37)上的均勻分布，生成容量為200的樣本，刪失比例P(T＜C)為0.5.在估計方法中采用Frank連接函數，τ=0.8.2種數據集下λ～和γ～的均值、標準差如表2所示.

表2 總體的連接函數形式變化時參數γ和λ的估計

由表2可以看出:當τ=0.8時，如果生成樣本的連接函數分別選取Clayton和Gumbel連接函數，采用獨立刪失的估計方法，得到的γ和λ估計量的偏差都很大.而采用本文方法(選取Frank連接函數，參數τ=0.8)得到的估計量都比較理想，估計量的偏差比較小.由此可見，在失效時間和刪失時間不相互獨立的情況下，采用獨立刪失方法進行估計可能會得到不理想的估計，所以應該采用帶相關性假定的模型進行分析.

其次，在連接函數的選取對估計量的影響方面，當連接函數的假定和總體不一致時本文方法能夠得到較穩(wěn)健的估計量.

2.3 不同刪失比例下參數估計的敏感性分析

選取C服從均勻U(0，3.3)和U(0，4.0)，刪失比例P(T＜C)分別為0.3、0.7.T服從γ=2，λ=0.5的威布爾分布，連接函數為Frank，τ=0.8，生成容量為200的數據集，估計不同數據集下λ～和γ～的均值、標準差，估計結果如表3所示.

表3 不同刪失比例下參數γ和λ的估計

由表3可以看出:如果連接函數及其參數τ的假定都是正確的，則在不同的刪失比例下，本文方法都能夠得到較好的估計，但采用獨立刪失方法得到的估計卻不理想.

綜合上述模擬計算的結果，可以得出:1)本文方法的參數估計效果比獨立刪失方法的參數估計效果更好;2)由連接函數不能準確識別或者參數τ不能正確識別所導致的偏差小于由獨立刪失錯誤假設所引起的偏差;3)當連接函數或參數τ的假定發(fā)生偏差時，本文方法依然能夠較穩(wěn)健.

3 討論

通過上面的敏感性分析，可以看出在帶信息的刪失數據下估計參數的標準差小于在獨立情況下估計參數的標準差.且對于不同連接函數、相關系數以及刪失比例，效果都較穩(wěn)健.另外，除了假設分布函數為威布爾分布外，還關于對數正態(tài)、指數分布等分布函數作了估計，效果也較好.因此，文本提供的方法具有一定的實用價值.

本文的工作還有許多地方可以進一步深入地研究，如在本文方法的框架下，繼續(xù)解決半參數模型的分布函數估計［15］;考慮有協變量影響下基于連接函數的帶信息刪失的非參數估計等.

［1］ Huang Jian，Wellner JA.Interval censored survival data:a review of recent progress［M］.New York:Springer-Verlag，1997:123-169.

［2］Sun Jianguo.The statistical analysis of interval-censored failure time data［M］.New York:Springer-Verlag，2006.

［3］呂秋萍，鄧文麗.區(qū)間刪失數據函數的均值估計［J］.江西師范大學學報:自然科學版，2011，35(1):96-100.

［4］王純杰.基于Copula函數的相依刪失數據的非參數統(tǒng)計推斷［D］.長春:吉林大學，2012.

［5］Siannis F，Copas J，Lu Baobing.Sensitivity analysis for informative censoring in parametric survivalmodels［J］.Biostatistics，2005，6(1):77-91.

［6］ Park Y，Lee Jenwei.One-and two-sample nonparametric inference procedures in the presence of amixture of independent and dependent censoring［J］.Biostatistics，2006，7(2):252-267.

［7］Zhang Zhigang，Sun Liuquan，Sun Jianguo，et al.Regression analysis of failure time data with informative interval censoring［J］.Statistics in Medicine，2007，26(12):2533-2546.

［8］ Huang Xuelin，Zhang Nan.Regression survival analysis with an assumed copula for dependent censoring:a sensitivity analysis approach ［J］.Biometrics，2008，64(4):1090-1099.

［9］鄧文麗，朱瑩瑩.一類保序最優(yōu)化問題的迭代算法［J］.統(tǒng)計與決策，2011(14):10-11.

［10］張連增，胡祥.Copula的參數與半參數估計方法的比較［J］.統(tǒng)計研究，2012，31(2):91-95.

［11］Clayton D G.Amodel for association in bivariate life tables and its application in epidemiological studies of familial tendency in chronic disease incidence ［J］.Biometrika，1978，65(1):141-151.

［12］ Hougaard P.A class ofmultivariate failure time distributions［J］.Biometrika，1986，73(3):671-678.

［13］Frank M J.On the sumultaneous association of F(x，y)and x+y-F(X，Y)［J］.Aequationes Mathematicae，1979，21(41):37-38.

［14］Nelsen R B.An introduction to copulas［M］.2nd ed.New York:Springer-Verlag，2006.

［15］楊金英，趙培信.缺失數據下ρ～混合誤差線性模型的參數估計［J］.西南大學學報:自然科學版，2012，34(9):35-37.