殷清燕 車露美 劉星宇
(1.西安建筑科技大學 西安 710055)(2.西安郵電大學 西安 710121)
生存分析(Survival Analysis)是研究生存現象和響應時間的一門統計學分支,已被廣泛運用于醫學領域,如臨床試驗、疾病診斷和預后分析等[1~4]。生存數據中經常出現刪失數據,導致生存時間信息不完全。傳統的統計學習方法無法很好地解決生存刪失數據的分析。設計適用于生存分析的機器學習算法,具有重要的研究意義和應用價值。
TCGA 是由美國國家癌癥研究所和人類基因組研究所合作建立的可公開獲取的癌癥組學數據庫[5]。TCGA 收錄了33 種癌癥數據,包含四種類型的 生 存 事 件:OS(Overall Survial),PFI(Progression-Free Interval),DFI(Disease-Free Interval),DSS(Disease-Specific Survival)。OS 指從觀察開始到死亡的時間;PFI 指從觀察開始到腫瘤進展的時間;DSS是指從診斷或治療開始到死亡的時間;DFI指治療后無病狀態到腫瘤進展的時間,進展事件可能是局部復發和遠處轉移。本研究首先選擇OS和PFI 事件進行初步Kaplan-Meier 估計,OS 事件進行深入的生存分析。

盡管Cox 比例風險模型在生存分析領域非常流行,但是比例風險假設在高維生存數據上難以成立。為了緩解這個問題,研究者們提出了具有較少的約束模型假設的非參數方法。隨機生存森林(Random Survival Forset,RSF)是Ishwaran提出的一種適用于生存數據的決策樹集成方法[9]。RSF 算法提出了CART 決策樹框架下的生存樹分裂規則,即將結點分裂成具有不相似生存時間的左右子結點,選取最大化結點間生存差異的特征作為最佳分裂特征。
RSF 算法通過Bootstrap 隨機采樣數據生成生存樹;分割節點時使用隨機選擇的特征子集;生存樹完全生長不進行剪枝;將多個生存樹的預測結果進行多數投票得到最終的預測。

表1 梯度提升生存樹的算法流程
成分梯度提升(CW-GBM)是2016 年提出的基于boosting和穩定性選擇的梯度提升方法[12]。可以提高計算可行性,在穩定性選擇中引入隨機排列來達到控制假陽性率。
原始的TCGA 數據集包含33 種癌癥類型的10915 個患者信息。本研究對原始數據進行數據清洗,得到預處理后的5384 個患者生存數據。生存信息包括OS,OS.time,PFI,PFI.time,臨床特征包括年齡、性別、腫瘤狀態和病理分期等。
Kaplan-Meier 估計通常是進行生存分析的第一步。根據離散的生存時間估計生存概率,繪制Kaplan-Meier 生存曲線,觀察和比較不同組患者的生存狀況。根據腫瘤狀態特征將患者分為有腫瘤和無腫瘤,分別繪制生存曲線(圖1)。上圖表示OS事件的生存曲線,無腫瘤的最低生存概率約為0.55,有腫瘤的生存曲線下降非常快;下圖表示PFI事件的生存曲線,無腫瘤的最低生存率約為0.72,有腫瘤的生存曲線下降得最快。有腫瘤的生存概率均在2000 天后低于0.1,有腫瘤和無腫瘤的生存曲線之間的差別非常大。有腫瘤的OS生存曲線和PFI 生存曲線具有相同的趨勢,PFI 生存曲線下降速度更快。

圖1 不同Tumor status下的K-M曲線
根據性別特征將患者分為男性(Male)和女性(Female),繪制生存曲線(圖2)。男性患者的OS和PFI 生存曲線均低于女性患者,男性患者的生存概率低于女性患者。男性和女性患者在PFI 生存曲線上的差別更大。

圖2 不同性別的K-M生存曲線
根據病理分期特征將患者分為六組。根據癌癥分期將患者分為六期,分別繪制曲線(圖3)。I/II/III 的最低生存概率分別約為0.51、0.28、0.21,癌癥四期(Stage IV)在t=2000 時的生存概率小于0.3。在t=3000 時不同病理分期的生存曲線之間有明顯的差別,病理分期等級越高的患者生存概率越低。

圖3 不同Stage下的K-M生存曲線
根據不同癌癥類型將患者進行分組,繪制生存曲線比較不同癌癥類型患者的生存狀況(圖4)。所有癌癥類型的患者在t=2000天(約5.5年)時達到最低生存概率。癌癥患者的生存概率與癌癥類型,性別,腫瘤狀況(Tumor status)和病理分期(Stage)相關。

圖4 不同癌癥類型的K-M生存曲線
本文使用OrdinalEncoder 對“cancer type”,“tumor_stage”,“tumor_status”進 行 編 碼;使 用One-HotEncoder 對“gender”特征進行編碼,使用StandardScaler 對“age”進行標準化。采用貝葉斯超參數優化方法對每個算法的超參數進行調整[13]。橫坐標表示參數的調參范圍,縱坐標表示一致性指數(C-index)值,繪制各個算法模型的超參數優化過程圖。
圖5 是CW-GBM 算法的調參結果圖。隨著“n_estimators”的增大,模型的C-index 保持增大。n_estimators 達到300 之后,C-index 指數趨于穩定,選擇最優n_estimators=493。圖6 和圖7 是GBM 和RSF算法的調參結果圖。“min_samples_leaf”的調參范圍是(1,25),min_samples_split 的調參范圍是(2,25),max_features 的調參范圍是(1,5),C-index 取值(0.80,0.82)。GBM 選取的最優超參數為min_samples_leaf=9,max_features=2,min_samples_split=11,n_estimators=188。RSF 算法選取的最優超 參 數 為min_samples_leaf=3,max_features=1,min_samples_split=16,n_estimators=21。

圖5 CW-GBM模型的超參數優化

圖6 GBM模型的超參數優化

圖7 RSF模型的超參數優化
C-index 在0.5~1 之間,0.5 為完全不一致,說明該模型沒有預測作用,1 為完全一致,說明該模型預測結果與實際完全一致。一般情況下C-index在0.50~0.70 為準確度較低:在0.71~0.90 之間為準確度中等;而高于0.90 則為高準確度。各算法在TCGA 測試集上的C-index 見圖8。RSF、GBM、COX、CW-GBM 在測試集上的一致性指數依次為0.87、0.85、0.79 和0.79。RSF 模型在預測測試集上的表現較好,具有較高的準確性。

圖8 各算法在TCGA測試集上的C-index
圖9 是各算法的time-dependent AUC 曲線圖,虛線表示AUC平均值[14]。Cox PH模型的AUC平均值為0.84,最低值為0.65。觀察初期的AUC值增長較塊,約1000 天后達到AUC 平均值。RSF 和GBM算法的AUC 曲線位于CoxPH 曲線上方,RSF 和GBM算法的性能優于CoxPH。因此,基于貝葉斯超參數優化的RSF和GBM具有較好的性能。

圖9 各算法在TCGA數據集上的AUC
選擇性能表現優異的隨機生存森林和梯度提升樹進行癌癥患者的生存風險分層。根據預測的風險評分值的四分位數,將患者分為四個風險組:低危組、中低危組、中高危組和高危組[15]。在隨機生存森林RSF中選取的風險評分分界點為0~1.84,1.84~8.45,8.45~22.53,22.53~100。在梯度提升樹GBM 算法中選取的風險評分分界點為0~15.96,15.96~29.33,29.33~43.54,43.54~100。
RSF 風險評分的四分位數分別為1.84,8.45,22.53,據此劃分的低危組1459 例、中低危組1458例、中高危組1167 例,高危組1750 例。GBM 風險評分的四分位數分別為15.96,29.33,43.54,據此劃分的低危組1461例、中低危組1458例、中高危組1165 例,高危組1750 例(表2)。RSF 和GBM 四個風險組的生存曲線如圖10 所示。可見兩個算法中相鄰風險組的生存曲線均有顯著性的差異。風險評分較高的風險組的生存曲線也較低,表明生存風險分層方法的有效性。

圖10 四個風險組的KM生存曲線

表2 TCGA患者的風險分層
結合表3 各風險組在性別、年齡、腫瘤狀態等方面的差異。在性別方面,我們可以看到各風險組中女性的數量普遍較低,尤其在RSF 低風險組中,女性人數僅為9 人,而男性人數為484 人。隨著風險程度的提高,男性人數逐漸增多,而在GBM 高風險組中,男性人數高達1078 人,女性人數也有667人。這表明,隨著年齡的增長,男性患腦腫瘤的風險較高。在年齡方面,我們可以看到各風險組的年齡分布有所不同。在RSF低風險組中,平均年齡為47.98歲,隨著年齡的增加,風險組平均年齡也逐漸上升。在GBM 高風險組中,平均年齡為64.7 歲,這說明隨著年齡的增長,患腦腫瘤的風險也在增加。

表3 各風險組的信息
針對不同風險組的特征信息進行分析,結果表明腫瘤狀態對患者生存風險的作用最顯著。隨著年齡的增長,男性患腦腫瘤的風險較高;患者的腫瘤狀態為“With tumor”時,患者的生存風險較高;其次年齡特征對患者生存風險的影響顯著。相比高風險組,低風險組患者的平均年齡較低。因此,針對不同風險組的人群,應采取有針對性的預防措施,以降低腦腫瘤的發病風險。
本研究針對TCGA 癌癥患者生存數據,比較CoxPH 模型、隨機生存森林RSF、梯度提升樹GBM和CW-GBM 算法的生存分析性能。對各個算法進行貝葉斯超參數優化,獲得C-index 和time-dependent AUC性能評價值。結果表明,隨機生存森林和梯度提升樹優于其他模型,在測試集上C-index 為0.87,time-dependent AUC 為0.90。其次,基于RSF和GBM 的生存風險評分,進行癌癥患者的生存風險分層,劃分為低危組、中低危組、中高危組和高危組。KM 生存曲線的實驗結果表明,隨機生存森林和梯度提升樹算法對TCGA 癌癥患者的生存風險具有很強的預測能力,在識別高危和低危患者群體方面有著顯著的判別作用。