區間刪失數據的3種統計模型分析及其SAS實現

張倩倩， 王純杰， 佟知真， 李純凈

(長春工業大學 基礎科學學院， 吉林 長春 130012)

?

區間刪失數據的3種統計模型分析及其SAS實現

張倩倩， 王純杰*， 佟知真， 李純凈

(長春工業大學 基礎科學學院， 吉林 長春 130012)

借助SAS9.4中PROC ICLIFETEST、PROC ICPHREG過程步編寫宏程序,同步實現了區間刪失數據的生存函數估計、廣義Log-Rank檢驗和PH類型回歸模型的統計推斷。結合回溯研究中368個樣本HIV-1感染時間的區間刪失數據給出實證分析。

區間刪失； ICLIFETEST； 廣義對數秩檢驗； ICPHREG； 宏程序

0 引 言

生存分析是對試驗或調查得到的人或生物的生存時間數據進行推斷,在醫學實踐中有著廣泛應用。一般稱給定事件的出現時間為生存時間[1]，分析生存時間數據通常意味著解決3個問題：估計生存函數，比較處理組或者生存函數，評估協變量的影響或者依靠生存時間的解釋變量。區間刪失數據是生存時間中越來越常見的一種數據，在過去幾十年里，出現了許多分析區間刪失數據的統計方法。Turnbull[2]找到了類似右刪失數據下的Kaplan-Meier估計的自相合算法來獲得生存函數估計；王弄升[3]2012年利用SAS軟件中宏程序%EMICM給出區間刪失數據生存函數的估計；Sun[4]等把Log-Rank檢驗推廣到區間刪失數據中，提出廣義對數秩檢驗；Finkelstein D M[5]給出區間刪失數據的COX回歸模型。但是基于SAS軟件還沒有完整的程序可以同步實現區間刪失數據3種統計分析任務。因此，文中借助SAS9.4中PROC ICLIFETEST[6]、PROC ICPHREG過程步編寫宏程序,實現了區間刪失數據的生存函數估計、廣義Log-Rank檢驗和PH比例風險類型的回歸模型統計推斷。

1 模型及方法介紹

1.1 區間刪失的數據類型

設T為非負的隨機變量，代表研究中個體的生存時間，對于區間刪失數據，只能知道T落在某個區間內，即T∈(L,R],在這里L≤R。區間刪失數據可分為Ⅰ型區間刪失與Ⅱ型區間刪失。Ⅰ型區間刪失數據可以表示{C,δ=I(T≤C)},C代表觀測時間，I(·)是示性函數。Ⅱ型區間刪失數據是包括有限區間的區間刪失數據(L,R)，假設每個個體觀測兩次，L,R是兩個隨機變量，L≤R以概率1成立。

1.2 區間刪失數據的非參數估計

1.2.1 EMICM算法[1]

SAS軟件計算非參數生存函數的方法是EMICM法。該算法是一個把自相合與ICM算法簡單結合的混合算法。考慮一個包含n個獨立個體帶有生存函數為S(t)的失效時間的研究。令Ti代表個體i(i=1,2，…,n)的生存時間。假設在Ti上的區間刪失數據被觀測如下：

O={(Li,Ri];i=1,2，…,n}

則區間刪失數據的似然函數為：

NPMLE算法就是要最大化此似然函數。在自相合算法中對數似然為：

Wellner,Zhan[7]指出，當NPMLE存在且唯一并且對數似然函數連續可微的時候,EMICM算法收斂到NPMLE。應用EMICM算法需要選擇一個收斂準則。這里選擇基于Robertson[8]1988年提出的Fenchel優化條件。這種準則下如果滿足

1.2.2 廣義對數秩檢驗[1]

j=1,2，…,m

類似的對于第l組，l=1,2，…,p+1，且j=1,2，…,m時，在sj處失效數和在風險數的估計為：

為了檢驗H0,給出統計量

此檢驗統計量是服從自由度為P的卡方分布。

1.2.3 PH比例風險模型[1]

則似然函數的對數形式為:

利用極大似然函數的得分方程估計β，基準生存函數S0可由連續的階梯函數來估計。

2 實證分析

2.1 區間刪失數據及變量情況

數據來源于美國和歐洲的16個研究中心的第5中心[9]，此中心采用回溯研究的方法,主要研究帶有血友病的病人感染HIV-1的風險。血友病人的治療血液制品來自于成千上萬個捐贈者的血漿制成的Ⅷ型凝血因子和Ⅳ型因子，所以這些病人都存在感染HIV-1的風險。文中對病人的HIV-1感染時間只得到了區間刪失數據，且病人依據他們每年獲得血液制品的平均劑量被分配到不同的組。來自第5中心的368個病人的觀測數據，不考慮進入研究的病人HIV-1抗體狀態，不接受或接受低劑量(1～20 000U)的Ⅷ型凝血因子(兩組病人的組別記為NF和LDF)見表1。

兩組病人的數量分別為236和13,這些數據以季度為單位，0代表研究開始時間(1978年1月1日)。

表1 HIV-1感染數據

續表1

2.2 SAS程序

對于上述血友病患者HIV-1感染的區間刪失數據，首先要分別建立兩組的SAS數據集，然后用ICLIFETEST過程步估計出兩個受試組HIV-1感染的時間生存函數，并檢驗兩組生存曲線的區別，最后用ICPHREG過程步建立PH比例風險回歸模型。區間刪失數據3種生存統計分析可由SAS宏程序%ICLIFETEST-PHREG同步實現。

2.2.1 創建NF組的數據集

data NF;

input lTime rTime @@;

Stage=0;

datalines;

55. 55. 56. 54. 53. 57.

56. 56. 54. 56. 54. 55.

...;

run;

2.2.2 創建LDF組的數據集

data LDF;

input lTime rTime @@;

Stage=1;

datalines;

720 920 025 57.

2326 821 2026 2527

...;

run;

2.2.3 估計兩組病人的生存概率及繪制生存曲線

data zq;

set LDF NF;

run;

proc iclifetest data=zq plot=s(cl) impute(seed=1234);/*ICLIFETEST過程步對應數據，畫出帶95%置信帶的生存曲線圖*/

strata stage;/*識別定義分組的變量*/

time (lTime,rTime);/*區間刪失的左右觀測時間*/

run;

2.2.4 兩組病人數據的廣義Log-Rank檢驗

proc iclifetest data=zq impute(seed=1234);

time (lTime,rTime);

test stage;/*檢驗的組別*/

run;

2.2.5 生存時間與協變量的PH比例風險回歸

proc icphreg data=zq;

class Stage / desc;

model (lTime,rTime) = Stage / basehaz=pch(intervals=(10));/*協變量Stage與區間刪失數據建立PH比例風險模型*/

hazardratio Stage;/*求風險比例*/

run;

上述5個程序是分別建立數據集，做區間刪失數據的3種統計模型分析的SAS程序，對于臨床試驗得到的這類區間刪失數據，如果利用SAS程序做分析，需要進行上述復雜的5個程序的運行，文中給出一個宏程序，同步實現上述模型分析。

2.2.6 區間刪失數據同步實現3種統計任務的SAS宏程序%ICLIFETEST-PHREG

%macro iclifetest-phreg(data,var1,var2);

proc iclifetest data=&data plot=s(cl) impute(seed=1234);

strata &var1;

time (lTime,rTime);

run;

proc iclifetest data=&data impute(seed=1234);

time(lTime,rTime);

test &var1;

run;

proc icphreg data=&data;

class &var2/ desc;

model (lTime,rTime) = &var2/ basehaz=pch(intervals=(10));

hazardratio &var1;

run;

%mend;

%iclifetest-phreg(zq,stage,stage)

宏程序%ICLIFETEST-PHREG可以對含二分變量和其他協變量的區間刪失數據同步實現生存函數的估計、廣義的Log-Rank檢驗及PH比例風險類型的回歸分析。其語法結構是：data：含二分變量和其他協變量的區間刪失數據集；Var1：某一個二分協變量；Var2：做回歸分析的所有的協變量集合。

2.3 結果分析

兩個不同組的非參數生存估計與生存函數曲線的宏程序結果分別如圖1和表2所示。

圖1 NF組與LDF組估計的生存函數

組別區間失效概率生存概率標準誤差NF4180.00550.99450.005419210.01950.98050.010922260.07420.92580.018727270.08560.91440.019928310.11500.88500.021132570.12340.87660.0214LDF0140.00001.00000.000015180.13900.86100.036819200.35350.64650.047121230.36480.63520.047224250.43730.56270.047426260.50620.49380.046027290.54560.45440.043930570.56060.43940.0432

圖1中相應的NF組的生存概率要大于LDF組的生存概率，即不接受Ⅷ型凝血因子的血友病實驗組的生存概率要高于接受低劑量Ⅷ型凝血因子的實驗組。表2中NF組的失效概率要小于LDF組的失效概率。

兩組生存概率是否相等的檢驗結果及回歸的顯著性檢驗結果分別見表3和表4。

表3 檢驗兩組生存概率是否相等

表4 回歸的顯著性檢驗

在顯著性水平為0.05的情況下，表3和表4的p值均小于0.000 1，即拒絕原假設，LDF組和NF組的生存概率有明顯區別。表4說明該回歸模型參數是顯著的。分組的風險比估計見表5。

表5 分組的風險比估計

表5結果顯示，患有血友病的病人接受低劑量的Ⅷ型凝血因子組感染HIV-1病毒的風險概率要高于不接受組病人感染HIV-1病毒的風險，且是后者的7.360倍。

3 結 語

基于區間刪失數據生存函數的估計算法除了EMICM算法還有自相合、ICM、Turnbull算法等。廣義對數秩檢驗也成為檢驗生存函數是否相等的有用工具。實例說明，NF與LDF治療組的生存函數有顯著差異。回歸模型檢驗出分組協變量對生存時間存在有效的影響。文中給出的SAS宏程序可以針對帶有協變量的區間刪失數據同步實現生存函數的估計、比較處理組及協變量對生存時間的影響這3種統計推斷，對臨床工作人員整理、分析實驗結果有很大幫助。

[1] Sun J. The statistical analysis of interval-censored failure time data[M]. [S.l.]： Springer Science & Business Media,Inc,2006.

[2] Turnbull B W. The empirical distribution function with arbitrarily grouped, censored and truncated data[J]. Journal of the Royal Statistical Society,1976,38(3):290-295.

[3] 王弄升,張晉聽,駱福添.含有區間刪失數據時的生存函數估計及其SAS實現[J].中國醫院統計,2012,19(1):1-4.

[4] Zhao Q, Sun J. Generalized log-rank test for mixed interval-censored failure time data [J]. Statistics in Medicine,2004,23(10):1621-1629.

[5] Finkelstein D M. A proportional hazards model for interval-censored failure time data [J]. Biometrics,1987,42(4):845-54.

[6] Guo C, Ying S, Gordon J. Analyzing interval-censored data with the ICLIFETEST procedure [J]. SAS Global Forum,2014,279:327-345.

[7] Jon A Wellner, Yihui Zhan. A hybrid algorithm for computation of the nonparametric maximum likelihood estimator from censored data [J]. Journal of the American Statistical Association,1997,439(92):945-959.

[8] Robertson T, Wright F T, Dykstra R L. Order restricted statistical inference [J]. Journal of the American Statistical Association,1990,85(410):111-112.

[9] Kroner B L, Rosenberg P S, Aledort L M, et al. HIV-1 infection incidence among persons with hemophilia in the United States and western Europe,1978-1990. Multicenter Hemophilia Cohort Study[J]. Journal of Acquired Immune Deficiency Syndromes,1994,7(3):279-86.

SAS-based study on three statistical models for interval censored data analysis

ZHANG Qianqian, WANG Chunjie*, TONG Zhizhen, LI Chunjing

(School of Basic Sciences, Changchun University of Technology, Changchun 130012， China)

Applying the procedures such as PROC ICLIFETEST and PROC ICPHREG in SAS9.4, we realize the statistic deduction for the survival function, Generalized log-rank test and PH regression for the interval censored data by means of macro programming. Interval censored data during the infection time from 368 samples are analyzed in the prospectively study.

interval censored; ICLIFETEST; generalized Log-rank test; ICPHREG; macro program.

2017-01-22

國家自然科學基金(青年基金)資助項目(11301037); 國家自然科學基金資助項目(11571051，11671054); 吉林省教育廳“十三五”規劃項目(2016317)

張倩倩(1991-)，女，漢族，山東德州人，長春工業大學碩士研究生，主要從事生存分析方向研究,E-mail:zhangqqxiao@163.com.*通訊作者：王純杰(1978-)，女，漢族，遼寧燈塔人，長春工業大學副教授，博士，主要從事生存分析方向研究,E-mail:cjwang2014@126.com.

10.15923/j.cnki.cn22-1382/t.2017.2.02

O 212.3

A

1674-1374(2017)02-0111-06