Aalen模型在醫(yī)學(xué)研究中的應(yīng)用

曹志強(qiáng)王 楊李 衛(wèi)△

Aalen模型在醫(yī)學(xué)研究中的應(yīng)用

曹志強(qiáng)1，2王 楊1李 衛(wèi)1△

Cox比例風(fēng)險(xiǎn)模型［1］是生存分析中最常用的模型，很多實(shí)際問(wèn)題中的協(xié)變量并不滿足比例風(fēng)險(xiǎn)，而且協(xié)變量的效應(yīng)可能隨時(shí)間變化。基于這些情況的考慮，Aalen提出了加法危險(xiǎn)率模型［2-3］，Aalen模型是Cox模型的補(bǔ)充。Aalen模型一個(gè)重要的特征就是其回歸系數(shù)是隨時(shí)間變化的函數(shù)，這種函數(shù)沒(méi)有特定的形式，也不依賴任何參數(shù)假定。相對(duì)于Cox模型的半?yún)?shù)本質(zhì)，Aalen模型是非參的，適合用于模型中含隨時(shí)間變化的協(xié)變量效應(yīng)的研究。

原理與方法

Aalen模型的基本形式如下［4］：

其中α0（t）是基線函數(shù)，Zj（t）是第j個(gè)協(xié)變量在t時(shí)刻的值。αj（t），j＝1，…，k是回歸參數(shù)，其作用等價(jià)于Cox模型中的回歸系數(shù)。在實(shí)際中，直接估計(jì)αj（t）是困難的，因而轉(zhuǎn)向估計(jì)與其等價(jià)的累積回歸系數(shù)，定義如下：

假設(shè)數(shù)據(jù)樣本的形式是［Tj，δj，Zj（t）］，j＝1，…，n，其中0＝T0＜T1＜T2＜…是排序好的時(shí)間，δj是判斷Tj是否刪失的示性變量是協(xié)變量向量。對(duì)于第j個(gè)個(gè)體，定義Yj（t）：如果個(gè)體j在t時(shí)刻是存活的，Yj（t）為1，否則為0。一般利用最小二乘法估計(jì)Aj（t），具體做法是先定義一個(gè)n×（k＋1）的設(shè)計(jì)矩陣X（t），其構(gòu)造如下：對(duì)于X（t）的第i行，令Xi（t）＝Y(jié)i（t）（1，Zj（t）），即如果第i個(gè)個(gè)體在t時(shí)刻是存活的，那么Xi（t）＝（1，Zj1（t），…，Zjk（t）），否則Xi（t）就是k＋1維的0向量。設(shè)I（t）是一個(gè)n× 1維的向量，如果第i個(gè)個(gè)體在t時(shí)刻死了，那么I（t）的第i個(gè)元素為1，否則為0。基于上面的構(gòu)造，累積回歸系數(shù)矩陣A（t）＝（A0（t），A1（t），…，Ak（t））T的最小二乘估計(jì)為：

在這個(gè)矩陣中，ID（t）是一個(gè)n×n的對(duì)角矩陣，其對(duì)角線元素等于I（t）。A（t）估計(jì)的最大時(shí)間Tmax是矩陣Xt（Ti）X（Ti）變?yōu)椴豢赡娴淖钚r(shí)間。

為了檢驗(yàn)協(xié)變量有無(wú)統(tǒng)計(jì)學(xué)意義，Aalen提出了下面的假設(shè)檢驗(yàn)，原假設(shè)為：

用向量U的第j個(gè)元素Uj檢驗(yàn)Hj：

其中，W（t）是一個(gè)對(duì)角矩陣的權(quán)重函數(shù)，對(duì)角線元素為Wj（t），j＝1，2，…，k＋1。這種非參數(shù)檢驗(yàn)方法只能檢驗(yàn)Xt（Ti）X（Ti）是滿秩的這段時(shí)間。Aalen考慮了兩種權(quán)重函數(shù)，一種是W（t）等于t時(shí)刻存活的人數(shù)，另一種為下面的（6）式。

Aalen從理論上證明了檢驗(yàn)統(tǒng)計(jì)量U服從漸近多元正態(tài)分布，用（6）式作為權(quán)重函數(shù)構(gòu)造出來(lái)的檢驗(yàn)統(tǒng)計(jì)量被稱之為T(mén)ST。對(duì)于模型的擬合優(yōu)度檢驗(yàn)，Aalen提出了廣義殘差法和Arjas plot法。第j個(gè)個(gè)體的廣義殘差定義如下：

其中Sj是第j個(gè)個(gè)體的確定或者刪失時(shí)間，Z0＝（1，Z01，…，Z0k）T是0時(shí)刻協(xié)變量的取值。如果模型擬合的好，那么可視為來(lái)自標(biāo)準(zhǔn)指數(shù)分布的樣本。Arjas圖的思想是比較累計(jì)度（cumulative intensity）和真實(shí)死亡數(shù)，假如模型正確，那么它們的值應(yīng)該差不多。具體做法是在Aalen模型中，對(duì)于每一確定時(shí)間Ti≤R，畫(huà)相對(duì)于i的圖像，如果模型正確，那么圖像近似為一條直線。

實(shí)例分析

一項(xiàng)以治療H1N1流感為目的的研究［5］，比較奧司他韋（達(dá)菲）和傳統(tǒng)中藥湯（麻杏石甘湯和銀翹散加減方）的治療效果，410例確診為輕癥H1N1流感的成年患者被隨機(jī)非盲分成4組：對(duì)照組、達(dá)菲組、中藥組、達(dá)菲加中藥組。目標(biāo)變量time為從入組治療到結(jié)束的時(shí)間；status指發(fā)熱是否消退；age是患者的年齡；g2、g3、g4是三個(gè)啞變量，分別指患者服用的是達(dá)菲、中藥湯、達(dá)菲加中藥；fb48h是發(fā)病至入組時(shí)間是否大于48小時(shí)的二值變量；s2、s3、s4是三個(gè)中心啞變量。為了解決結(jié)點(diǎn)問(wèn)題，將time每個(gè)值加上［0，1］之間的隨機(jī)數(shù)。

首先用age、g2、g3、g4、fb48h、s2、s3、s4作為協(xié)變量，進(jìn)行Cox回歸，結(jié)果見(jiàn)表1。

表1 Cox回歸結(jié)果

從表1可知，age、fb48h、s3無(wú)統(tǒng)計(jì)學(xué)意義，s2在0.1水平下有統(tǒng)計(jì)學(xué)意義，其他變量在0.05水平下均有統(tǒng)計(jì)學(xué)意義。本文研究的目標(biāo)變量是發(fā)熱持續(xù)時(shí)間，相對(duì)不吃藥，如果吃了某種藥后發(fā)熱時(shí)間能夠顯著降低，說(shuō)明該藥有效（此時(shí)HR＞1）。從結(jié)果來(lái)看，達(dá)菲、中藥湯、達(dá)菲加中藥都能有效治療H1N1流感。從表1還能得出，達(dá)菲的HR值比中藥湯的要高一點(diǎn)，達(dá)菲加中藥的HR值比單純達(dá)菲或中藥湯的都要高。然而，通過(guò)Wald檢驗(yàn)發(fā)現(xiàn)，達(dá)菲相對(duì)中藥湯的HR值不顯著，但達(dá)菲加中藥相對(duì)達(dá)菲的HR是顯著的，HR置信區(qū)間為［0.072，0.871］，P＝0.014；達(dá)菲加中藥相對(duì)中藥湯的HR值也是顯著的，HR值的置信區(qū)間為［0.098，0.919］，P＝0.009。

用log-log圖檢驗(yàn)Cox模型中的組別變量是否服從比例風(fēng)險(xiǎn)假定。可見(jiàn)中藥組的與其它三組有交叉，因此，模型的比例風(fēng)險(xiǎn)假定可能存在問(wèn)題。文獻(xiàn)［6］指出，當(dāng)Cox模型中的協(xié)變量不滿足比例風(fēng)險(xiǎn)時(shí)，可采用Aalen模型分析。

在具體運(yùn)用Aalen模型之前，根據(jù)模型原理可以推斷出累積回歸系數(shù)的一些特征。如果Aalen模型中某個(gè)協(xié)變量的回歸系數(shù)是常數(shù)，即α（t）＝a，那么其在t時(shí)刻的累積回歸系數(shù)應(yīng)為A（t）＝at這樣的一條直線。假設(shè)風(fēng)險(xiǎn)因子超過(guò)了t0時(shí)刻，比如t0＝20小時(shí)之后，對(duì)風(fēng)險(xiǎn)函數(shù)不再有影響，則其累積回歸系數(shù)在20小時(shí)后應(yīng)該等于常數(shù)。如果變量在模型中有統(tǒng)計(jì)學(xué)意義，那么其累積回歸系數(shù)的置信區(qū)間不應(yīng)該包含0。

表2 Aalen模型變量的檢驗(yàn)結(jié)果

表2是用TST統(tǒng)計(jì)量檢驗(yàn)Aalen模型中的變量有無(wú)統(tǒng)計(jì)學(xué)意義。從P值來(lái)看，age、fb48h、s3無(wú)統(tǒng)計(jì)學(xué)意義，其他變量均有統(tǒng)計(jì)學(xué)意義，這和Cox模型選擇變量的結(jié)果類(lèi)似。

三組藥的累積回歸系數(shù)圖展示了三組藥在各時(shí)段如何影響風(fēng)險(xiǎn)函數(shù)。達(dá)菲的累積回歸系數(shù)在前38小時(shí)持續(xù)遞增，但超過(guò)38小時(shí)后斜率趨平且有降低趨勢(shì)，表明達(dá)菲對(duì)治療H1N1流感的效果主要體現(xiàn)在前38小時(shí)。中藥湯的累積回歸系數(shù)在前19小時(shí)顯著遞增，之后出現(xiàn)類(lèi)似的轉(zhuǎn)平和減低趨勢(shì)，表明中藥湯的療效在前19個(gè)小時(shí)體現(xiàn)更為明確。至于達(dá)菲加中藥，其累積回歸系數(shù)在前38小時(shí)平穩(wěn)增加，之后增加的趨勢(shì)減緩。因此達(dá)菲加中藥在前38小時(shí)的療效顯著，38小時(shí)后療效趨弱。

達(dá)菲、中藥湯、達(dá)菲加中藥的累積回歸系數(shù)的斜率均為正數(shù)，說(shuō)明相對(duì)不吃藥，它們對(duì)治療H1N1流感都有效。有些學(xué)者［7］提出根據(jù)累積回歸系數(shù)的斜率來(lái)度量協(xié)變量的影響，目前這在文獻(xiàn)中不常見(jiàn)，在此我們不提倡根據(jù)斜率的大小就定量地確定達(dá)菲、中藥湯、達(dá)菲加中藥的療效到底有多好。

達(dá)菲加中藥的累積回歸系數(shù)圖的趨勢(shì)一直遞增，而達(dá)菲或中藥湯的在后面時(shí)段下降，以至于達(dá)菲在50小時(shí)時(shí)的置信下限包含0，中藥湯在37小時(shí)時(shí)的置信下限包含0。該信息是對(duì)各藥物特點(diǎn)（即持續(xù)有效時(shí)間）的體現(xiàn)，也說(shuō)明樣本量隨時(shí)間減少，在一定程度上影響估計(jì)的精度。達(dá)菲、達(dá)菲加中藥比較好地符合Cox比例風(fēng)險(xiǎn)假定，這一點(diǎn)與之前l(fā)og-log圖中所得的信息一致。

討 論

在醫(yī)學(xué)研究中，生存資料的多因素分析常采用Cox模型，其回歸系數(shù)度量的是相對(duì)風(fēng)險(xiǎn)。Aalen模型是從絕對(duì)風(fēng)險(xiǎn)的角度考慮生存時(shí)間和協(xié)變量之間的關(guān)系。在對(duì)兩種模型結(jié)果進(jìn)行比較時(shí)，如果基于非參數(shù)模型在同一變量上給出的P值，大于半?yún)?shù)和參數(shù)模型（或前者無(wú)、而后者有統(tǒng)計(jì)學(xué)意義），而下結(jié)論Aalen模型不如Cox模型好，這是不合理的，其原因是兩個(gè)模型對(duì)應(yīng)的檢驗(yàn)本身不同。

實(shí)例分析表明，當(dāng)Cox模型的比例風(fēng)險(xiǎn)假定不滿足時(shí)，我們可以用Aalen模型分析。在檢驗(yàn)協(xié)變量有無(wú)統(tǒng)計(jì)學(xué)意義時(shí)，兩個(gè)模型得出了相似的結(jié)果。如果協(xié)變量在不同的時(shí)段對(duì)風(fēng)險(xiǎn)函數(shù)有不同的影響，成比例風(fēng)險(xiǎn)或者不成比例風(fēng)險(xiǎn)，用Aalen模型分析累積回歸系數(shù)圖是有幫助的。Aalen提倡將累加回歸系數(shù)圖作為診斷工具，觀察各變量在各時(shí)段如何影響風(fēng)險(xiǎn)函數(shù)。在本篇的例子中，通過(guò)觀察累積回歸系數(shù)圖，可以判斷藥物間的療效是否大致符合比例風(fēng)險(xiǎn)，以及了解治療H1N1流感中各藥物療效發(fā)揮時(shí)段的特點(diǎn)。可以說(shuō)，累積回歸系數(shù)圖直觀展示了協(xié)變量影響風(fēng)險(xiǎn)函數(shù)的本質(zhì)。

實(shí)際應(yīng)用中，Aalen模型遠(yuǎn)沒(méi)有Cox模型應(yīng)用廣，一個(gè)原因是SAS和SPSS軟件中沒(méi)有現(xiàn)成Aalen模型的程序，而Cox模型的程序幾乎在每一款統(tǒng)計(jì)軟件中都早已存在。如今R軟件的兩個(gè)包survival和timereg提供了處理Aalen模型的程序，讓使用Aalen模型作為Cox模型的有效補(bǔ)充成為了可能。很多有關(guān)Aalen模型的文獻(xiàn)，包括Aalen的原創(chuàng)，都是提倡將Aalen模型視為Cox模型的補(bǔ)充，把兩個(gè)模型結(jié)合起來(lái)分析問(wèn)題。總之，這兩個(gè)模型提供了數(shù)據(jù)的不同信息，它們不應(yīng)該視為相互替代，而應(yīng)是相互補(bǔ)充。這樣，我們對(duì)數(shù)據(jù)和問(wèn)題就會(huì)有更全面和更深刻的認(rèn)識(shí)和理解。

1.Cox DR.Regression models and life-tables.Journal of the Royal Statistical Society，Series B（Methodological），1972：187-220.

2.Aalen OO.A linear regression model for the analysis of life times.Statistics in medicine，1989，8（8）：907-925.

3.Aalen OO.Further results on the non-parametric linear regression model in survival analysis.Statistics in medicine，1993，12（17）：1569-1588.

4.Klein J，Moeschberger M.Survival Analysis：Techniques for Censored and Truncated Data，Springer，1997.

5.Wang C，Cao B，Liu QQ，et al.Oseltamivir compared with the Chinese traditional therapy maxingshigan-yinqiaosan in the treatment of H1N1 influenza：a randomized trial.Annals of Internal Medicine，2011，155（4）：217-25.

6.Abadi A，et al.Comparison of Aalen′s additive and Cox proportional hazards models for breast cancer survival：analysis of population-based data from British Columbia，Canada.Asian Pacific Journal of Cancer Prevention，2011，12：3113-3116.

7.Torner A.Proportional hazards and additive regression analysis of survival for severe breast cancer.Stockholm University，2004.

附錄：本文的R程序

（責(zé)任編輯：郭海強(qiáng)）

1中國(guó)醫(yī)學(xué)科學(xué)院，北京協(xié)和醫(yī)學(xué)院，國(guó)家心血管病中心，阜外心血管病醫(yī)院，心血管疾病國(guó)家重點(diǎn)實(shí)驗(yàn)室（100037）

2北京師范大學(xué)數(shù)學(xué)科學(xué)學(xué)院

△通信作者：李衛(wèi)，E-mail：liwei＠m(xù) rbc-nccd.com