唐璐薇,呂 超,文秋月,韋程東

(1.廣西科技師范學(xué)院 數(shù)學(xué)與計算機(jī)科學(xué)學(xué)院,廣西 來賓 546100; 2.柳州市婦幼保健院,廣西 柳州 545000; 3.南寧師范大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院,廣西 南寧 530000)

0 引 言

生存分析是一種既考慮結(jié)果又考慮生存時間的統(tǒng)計方法,其可充分利用截尾數(shù)據(jù)所提供的不完全信息,對生存時間的分布特征進(jìn)行描述,以及對影響生存時間的主要因素進(jìn)行分析.比例風(fēng)險模型是研究生存數(shù)據(jù)最廣泛的統(tǒng)計半?yún)?shù)模型之一[1-2],但其對數(shù)據(jù)的擬合效果并不理想.因此,需要我們尋求一個比例風(fēng)險模型的替代模型,即加性風(fēng)險模型.與比例風(fēng)險模型不同的是,加性風(fēng)險模型是假設(shè)基底風(fēng)險函數(shù)與協(xié)變量之間的一個加性結(jié)構(gòu).在實際應(yīng)用中,加性風(fēng)險模型對數(shù)據(jù)的擬合效果更好,正是由于它這一特性,加性風(fēng)險模型中的回歸參數(shù)更容易解釋實際意義[3-4].

目前,運用不同模型對數(shù)據(jù)進(jìn)行擬合的研究有很多.顧劉金首先利用Cox模型對數(shù)據(jù)進(jìn)行擬合,但由于Cox模型具有比例風(fēng)險性,故基于SPSS采用Logistic回歸模型對數(shù)據(jù)進(jìn)行擬合得出結(jié)果[5].Jardim等利用微分方程模型對葡萄牙的COVID-19相關(guān)數(shù)據(jù)進(jìn)行擬合分析,同時指出該方法可針對不同的流行病學(xué)數(shù)據(jù)進(jìn)行擬合分析[6].此外,很多學(xué)者也運用B樣條擬合方法進(jìn)行研究.曾卓等提出了基于三次B樣條小波的曲線擬合方法,該方法行之有效,為曲線擬合提供了一種兼顧擬合精度、光順性與加工精度的方法[7].張永華等針對傳統(tǒng)幾何軌跡跟蹤算法切向角獲取依賴高精度慣導(dǎo)設(shè)備的問題,提出了基于三次B樣條曲線擬合的軌跡跟蹤算法,該算法有效解決了傳統(tǒng)算法的問題[8].Bi等基于全局?jǐn)M合法的缺點,提出一種通用、快速的B樣條誤差擬合方案,該方法與傳統(tǒng)的擬合方法相比,顯著地提高了工作效率[9].龐宇等提出一種結(jié)合幅度系數(shù)法和波形特征法對脈搏波采用三次B樣條曲線擬合的血壓檢測算法,該算法能夠提高特征點的準(zhǔn)確率[10].蘆穗豪等提出一種基于改進(jìn)麻雀搜索算法的B樣條曲線擬合方法,旨在利用最少控制點高效地達(dá)到曲線擬合的目標(biāo)精度,進(jìn)而提升傳統(tǒng)建模方法的精度和效率[11].徐超清等基于神經(jīng)纖維走向信息考量的問題,提出一種基于B樣條擬合與回歸模型的腦神經(jīng)纖維聚類方法,該方法在功能區(qū)層面的聚類可以更有效地分割出具有解剖學(xué)結(jié)構(gòu)的腦神經(jīng)纖維[12].

上述研究大多是運用傳統(tǒng)的模型對生存數(shù)據(jù)進(jìn)行分析,且運用三次B樣條平滑函數(shù)對生存數(shù)據(jù)進(jìn)行擬合,而利用加性風(fēng)險模型進(jìn)行相關(guān)研究的較少.本文運用三次B樣條擬合方法,研究經(jīng)藥物治療后白血病患者的生存周期與實際情況的擬合程度,從而判斷藥物治療的有效性,并運用單變量加性風(fēng)險模型和多變量加性風(fēng)險模型找出對生存時間影響顯著的協(xié)變量,作出生存時間的預(yù)測圖,從而為決策者提供一定的參考依據(jù).

1 模型介紹

1.1 加性風(fēng)險模型

假設(shè)研究隊列數(shù)據(jù)包含N個獨立的樣本,設(shè)T為失效時間,C為相應(yīng)的截尾時間,Z(t)(0≤t≤τ)為協(xié)變量過程的向量,其中τ<∞表示后續(xù)時間.由于數(shù)據(jù)右刪失,所以只能觀察X和δ,其中X=min(T,C′)和δ=I(T≤C).用Yi(t)=I(Ti≥t)表示歷險過程的示性函數(shù),Ni(t)=ΔiI(Ti≤t)表示計數(shù)過程,這里I(·)是示性函數(shù).假定T和C是條件獨立的,可考慮以下加性風(fēng)險模型.當(dāng)給定協(xié)變量Zi(t)時,失效時間T的風(fēng)險函數(shù)為:

(1)

其中,λ0(t)為一個未指定的基底風(fēng)險函數(shù),β=(β1,β2,…,βp)為未知回歸參數(shù)的p維向量.

1.2 三次B樣條插值函數(shù)

B樣條曲線是在Bezier曲線基礎(chǔ)上發(fā)展起來的,其克服了Bezier曲線整體控制性所帶來的不方便.它是通過逼近一組控制點生成的曲線,計算公式為:

(2)

其中,pk為輸入的一組數(shù)據(jù)中的第k個數(shù)據(jù)點;Nk,d為B樣條混合函數(shù),這里為加性風(fēng)險模型,k為第k個混合函數(shù),d為次數(shù).本文采用三次B樣條插值函數(shù)進(jìn)行平滑,即d=3.

B樣條曲線保留了Bezier曲線的優(yōu)勢,同時曲線在拼接時又比Bezier曲線方便,在修改時可做局部修改,而且在調(diào)整某一控制點時不會影響整條曲線的趨勢.其性質(zhì)主要有以下幾點:

(1)局部性.K階B樣條曲線上的一點至多與k個控制頂點有關(guān),與其他控制頂點無關(guān).

(2)幾何不變形.B樣條曲線的形狀和位置與坐標(biāo)系的選擇無關(guān),不管坐標(biāo)系如何變化,B樣條曲線的形狀仍保持原樣.

(3)凸包性.與Bezier曲線一樣,B樣條曲線落在Pi構(gòu)成的凸包中,其凸包區(qū)域小于或等于同一組控制頂點定義的Bezier曲線凸包區(qū)域.

2 白血病數(shù)據(jù)的三次B樣條平滑

常見的生存數(shù)據(jù)處理模型有比例風(fēng)險模型(Cox模型)、加性風(fēng)險模型等,二者在處理生存數(shù)據(jù)上都有一定的優(yōu)勢.下面給出Cox模型和加性風(fēng)險模型對白血病數(shù)據(jù)分析的結(jié)果.

(1)Cox模型數(shù)據(jù)擬合.根據(jù)已知白血病數(shù)據(jù)中的自變量,選取6個自變量作為協(xié)變量:age(年齡)、sex(性別)、 ph.ecog(ECOG評分)、ph.karno(醫(yī)師的Karnofsky評分)、meal.cal(用餐時消耗的卡路里)、wt.loss(最近6個月的體重減輕),并運用Cox模型對其進(jìn)行分析,結(jié)果見表1.

表1 Cox模型對6個協(xié)變量的分析

由表1可知,性別和ECOG評分與生存時間顯著相關(guān),且P值遠(yuǎn)小于0.05,模型擬合效果較好.

(2)加性風(fēng)險模型數(shù)據(jù)擬合.為與Cox模型擬合效果進(jìn)行對比,下面運用加性風(fēng)險模型對上述6個協(xié)變量進(jìn)行分析,查看加性風(fēng)險模型對數(shù)據(jù)的擬合效果及相關(guān)性,結(jié)果見表2.

表2 加性風(fēng)險模型對6個協(xié)變量的分析

由表2可知,加性風(fēng)險模型對數(shù)據(jù)的擬合效果與Cox模型一樣,性別和ECOG評分與生存時間顯著相關(guān),且P值遠(yuǎn)小于0.05,模型的擬合效果較好.但加性風(fēng)險模型還給出了對生存時間的整體方差解釋率,更進(jìn)一步說明了該模型對數(shù)據(jù)的擬合效果.因此,本文在三次B樣條平滑函數(shù)下,運用加性風(fēng)險模型來研究白血病的治療數(shù)據(jù),從而給出相關(guān)結(jié)論.

2.1 三次B樣條平滑

為減輕協(xié)變量對術(shù)后患者生存時間的影響,本文利用三次B樣條平滑法對受干擾的協(xié)變量進(jìn)行平滑,以消除協(xié)變量對響應(yīng)變量即生存時間的影響.下面首先討論利用三次B樣條平滑法對受干擾的加性風(fēng)險模型中的協(xié)變量進(jìn)行平滑,然后討論單變量加性風(fēng)險模型和多變量加性風(fēng)險模型中協(xié)變量與生存時間的關(guān)系.三次B樣條平滑過程分以下幾點討論:

(1)59歲以下的男性,其生存時間與用餐時消耗的卡路里、最近6個月體重減輕的關(guān)系見圖1.其中,第一行指定的平滑系數(shù)分別為k=3,k=5;第二行指定的平滑系數(shù)分別為k=5,k=7.由圖1可知,從整體看,白血病患者的生存時間隨著用餐時消耗卡路里的增加而增加,但不同平滑程度的選擇影響了局部趨勢的解釋,特別是不同年齡生存時間的變化趨勢是不同的,且最近6個月的體重減輕對生存時間的影響也較明顯.

圖1 59歲以下男性兩協(xié)變量與生存時間的關(guān)系

(2)59歲以上的男性,其余條件與(1)相同,得到的結(jié)果見圖2.由圖2可知,59歲以上的男性,隨著年齡的增加,其消耗的卡路里逐漸降低,而相對(1)來說,近6個月的體重減輕呈細(xì)微的增加趨勢.

(3)59歲以下的女性,其余條件與(1)相同,結(jié)果見圖3.由圖3可知,59歲以下的女性,隨著年齡的增加,其消耗的卡路里趨于平穩(wěn),而隨著近6個月體重減輕的越來越多,生存時間呈先下降后上增的趨勢.

圖3 59歲以下女性兩協(xié)變量與生存時間的關(guān)系

(4)59歲以上的女性,其余條件與(1)相同,結(jié)果見圖4.由圖4可知,59歲以上的女性,從整體看與圖3中59歲以下的女性分析結(jié)果相似,但在平滑系數(shù)k=7時,脂肪消耗量發(fā)生細(xì)微變化,隨著脂肪消耗的增多,生存時間呈細(xì)微上增的趨勢.

通過三次B樣條平滑對年齡(age)、性別(sex)、ECOG評分(ph.ecog,0=好,5=死)、醫(yī)師的Karnofsky評分(ph.karno,0=差,100=好)、用餐時消耗的卡路里(meal.cal)和最近6個月的體重減輕(wt.loss)等6個協(xié)變量中受干擾的協(xié)變量進(jìn)行平滑后,暫時不能確定哪些協(xié)變量對生存時間的影響顯著.下面研究生存時間的單變量加性風(fēng)險模型,并根據(jù)研究結(jié)果選擇協(xié)變量進(jìn)行多變量的加性風(fēng)險模型研究,以探索這些協(xié)變量對生存時間的影響.

2.2 單變量加性風(fēng)險模型研究

選擇協(xié)變量:年齡(age)、性別(sex)、ECOG評分(ph.ecog,0=好,5=死)、醫(yī)師的Karnofsky評分(ph.karno,0=差,100=好)、用餐時消耗的卡路里(meal.cal)和最近6個月的體重減輕(wt.loss),并用這些協(xié)變量擬合加性風(fēng)險模型,結(jié)果見表3.其中,s()為運用三次B樣條平滑后的符號表示.

表3 單變量加性風(fēng)險模型擬合結(jié)果

由表3可知,性別、年齡、ph.karno和ph.ecog變量具有較好的統(tǒng)計學(xué)意義(P<0.05),說明它們對生存時間的影響是顯著的.此外,年齡和ph.ecog具有正β系數(shù),而性別和ph.karno具有負(fù)β系數(shù).因此,年齡較大和ph.ecog較高與事件發(fā)生率呈正相關(guān),而女性(sex=2)和ph.karno則與事件發(fā)生率呈負(fù)相關(guān),即年齡和ph.ecog是死亡的危險因素,女性性別和醫(yī)師的Karnofsky評分是死亡的保護(hù)因素.

2.3 多變量加性風(fēng)險模型研究

為研究性別、年齡、ph.ecog和ph.karno如何共同影響生存時間,本文利用多變量加性風(fēng)險模型進(jìn)行分析,結(jié)果見表4.

表4 多變量加性風(fēng)險模型擬合結(jié)果

由表4可知,加性風(fēng)險模型對生存時間的整體方差解釋率達(dá)78.5%,擬合效果較好.在多元加性風(fēng)險模型分析中,協(xié)變量性別和ph.ecog保持著顯著性(P<0.05),但協(xié)變量年齡和ph.karno不顯著(年齡:P=0.171 226,ph.karno:P=0.186 82均大于0.05).

性別的P值為0.000 712,風(fēng)險比HR=exp(coef)=0.56,表明患者的性別與死亡風(fēng)險降低之間有很強的關(guān)系.在保持其他協(xié)變量不變的前提條件下,女性(sex=2)相比男性,其死亡風(fēng)險低44%.可見,性別為女性與良好的預(yù)后相關(guān).

同樣,ph.ecog的P值為0.000 323,風(fēng)險比HR=1.88,表明ph.ecog值與死亡風(fēng)險增加之間有很強的關(guān)系.相比之下,年齡的P值為0.171 226,風(fēng)險比HR=exp(coef)=1.01,95%置信區(qū)間為0.99～1.03.由于HR的置信區(qū)間為1,因此該結(jié)果表明,在調(diào)整ph.ecog值和患者性別后,年齡對HR差異的貢獻(xiàn)較小,且不顯著.

為研究白血病患者治療的生存數(shù)據(jù)與治療的關(guān)系,以及白血病治療的有效性,針對上述對協(xié)變量的研究,本文選取性別和ph.ecog作為協(xié)變量來擬合加性風(fēng)險模型.使用K-M估計方法,以存活天數(shù)為橫坐標(biāo),以生存率為縱坐標(biāo),給出估計的生存率和生存率的置信區(qū)間,見圖5.

圖5 存活天數(shù)與生存率的關(guān)系

如圖5所示,隨著存活天數(shù)的增多,術(shù)后存活率也逐漸下降,這與自然規(guī)律相符.此外,由上述對多變量加性風(fēng)險模型協(xié)變量的分析可知,性別和ph.ecog對模型具有顯著影響.性別和ph.ecog對生存時間和生存率的影響見圖6.

圖6 性別和ph.ecog對生存時間的影響

由圖6可知,女性的整體生存率高于男性,在ph.ecog評分表中對應(yīng)的評分為5,且生存率隨著生存時間的變化逐漸下降.但這種差異是否顯著并不確定,還需要進(jìn)行統(tǒng)計檢驗,檢驗結(jié)果見表5.

表5 性別和ph.ecog差異顯著性檢驗

由表5可知,將性別和ph.ecog綜合考慮,得到的P值遠(yuǎn)小于0.05,說明男女之間以及不同的ECOG評分,其生存率是有差異的.

與K-M生存曲線不同的是,加性風(fēng)險模型擬合曲線得到的生存率是在運用三次B樣條平滑其他協(xié)變量后所預(yù)測的生存率,并不是實際觀察到的生存率.加性風(fēng)險模型擬合生存率曲線見圖7.

圖7 加性風(fēng)險模型擬合生存率曲線

由圖7可知,預(yù)測患者的生存率隨著生存時間的增多而逐漸下降,與使用K-M估計方法對白血病患者的生存時間和生存率所作出的曲線趨勢一致.由于本文研究的是白血病治療數(shù)據(jù),涉及對醫(yī)療效果的評價,因此選擇精確度較高,且對解決高度非線性預(yù)測問題有突出能力的加性風(fēng)險模型來對數(shù)據(jù)進(jìn)行擬合是可行的.

3 結(jié) 論

本文研究的白血病數(shù)據(jù)涉及分析白血病藥物治療的效果,精度要求較高,需要一個能夠?qū)ι鏀?shù)據(jù)進(jìn)行分析且擬合效果較好的模型.本文利用加性風(fēng)險模型對其進(jìn)行擬合研究,并通過圖形展示患者生存時間的變化.為避免數(shù)據(jù)缺失給加性模型研究帶來誤差,利用三次B樣條平滑函數(shù)對其進(jìn)行平滑,并對刪失的數(shù)據(jù)進(jìn)行處理,以避免造成對已有數(shù)據(jù)的影響.通過對數(shù)據(jù)分析可知,不同平滑程度對擬合效果影響不大,不能直接表達(dá)哪些協(xié)變量對生存時間影響顯著.本文不足之處在于運用三次B樣條平滑函數(shù)對加性風(fēng)險模型進(jìn)行平滑后,只對白血病治療數(shù)據(jù)進(jìn)行研究,沒有考慮到其他疾病數(shù)據(jù)的擬合效果是否一樣.此外,針對白血病的治療數(shù)據(jù)也只運用了加性風(fēng)險模型對其進(jìn)行研究,沒有運用其他模型進(jìn)行擬合.今后的研究將針對這兩個問題進(jìn)行深入探討,從多方面進(jìn)行生存分析研究.