結合自噬相關基因表達譜和臨床因素的肺腺癌預后預測模型

李麗希，黃 鋼

(1.上海理工大學 健康科學與工程學院，上海 200090；2.上海市分子影像重點實驗室，上海 201200)

肺癌是全球死亡率最高的癌癥之一，非小細胞肺癌(Non-small cell lung cancer，NSCLC)是肺癌中最常見的類型，約占所有肺癌病例的80%[1]。肺腺癌(Lung adenocarcinoma，LUAD)是非小細胞肺癌的主要亞型之一，對全球不吸煙者而言是致死率最高的疾病[2-3]。由于LUAD在早期容易轉移復發，LUAD患者的預后效果很差，平均5年生存率低于20%[4]。在臨床實踐中，腫瘤分期系統已廣泛應用于癌癥患者的指導治療和預后評估。然而，預后的判斷通常只基于固有的解剖學信息，由于肺腺癌的異質性，很難預測疾病的發展。因此，迫切需要尋找有效的預后生物標志物來幫助臨床醫生做出準確的肺腺癌診斷，預測臨床結果，為個體化醫學提供參考。

過去幾年里，研究發現自噬在腫瘤的發生過程中發揮了重要的作用[6-9]。自噬是一個復雜的生理病理過程，自噬的溶酶體降解功能在細胞生理學中起著至關重要的作用，如適應代謝應激、清除危險物質(如蛋白質聚集體、受損或老化的細胞器、細胞內病原體)、細胞分化和發育過程中的更新等[10]。在癌癥中，自噬具有雙重作用，它既能夠抑制良性腫瘤的生長，也可以促進晚期癌癥的發展[11]。目前，許多研究小組已經確定把自噬作為癌癥治療的潛在靶點。

本項研究構建了一個結合多個自噬相關基因和臨床參數的模型來預測LUAD患者的預后。從TCGA數據庫的LUAD數據中篩選出表達具有顯著差異的自噬相關基因，對差異自噬相關基因進行單因素Cox回歸分析來確定與LUAD患者生存相關的候選基因，然后使用lasso回歸模型篩選出預后相關基因，對預后相關基因進行多因素Cox分析，構建風險評分模型，并對風險評分模型進行內部驗證和外部驗證。最后將風險評分與臨床參數結合，構建了預測患者生存概率的列線圖模型，使用一致性指數(Concordance index, C-index)、校準曲線和ROC曲線來評估模型的性能。

1 材料和方法

1.1 數據收集和數據預處理

在人類自噬專用數據庫HADb(http://www.autophagy.lu/)、ARN數據庫(http://autophagyregulation.org)、自噬數據庫(http://www.tanpaku.org/autophagy/index.html)上下載了自噬相關基因共1 417個。從TCGA(https://portal.gdc.cancer.gov/)下載LUAD的COUNT數據和FPKM數據各585例，包含526例肺腺癌樣本和59例癌旁正常肺組織。從UCSC Xena(https://xenabrowser.net/)下載TCGA-LUAD的生存信息641例和臨床數據706例。對于TCGA數據，過濾生存信息、腫瘤分期信息、年齡、性別和復發信息不完整的樣本，保留了TCGA的417例肺腺癌樣本和48例正常樣本。在TCGA數據中篩選出自噬相關基因的信息，并將癌癥樣本(n=417)隨機分配為訓練組和測試組，比例為1∶1。

同時，在GEO(https://www.ncbi.nlm.nih.gov/geo/)數據庫下載了GSE50081數據集用于外部驗證，該數據集包括了127例肺腺癌樣本和54例正常樣本。

1.2 預后相關基因的篩選

使用“limma”包對自噬相關基因進行差異分析，差異基因篩選標準為：|logFC|>1.5，P<0.05。

對表達具有顯著差異的自噬相關基因使用單因素Cox比例風險回歸分析篩選出候選基因，篩選閾值為：風險比HR≠1，p<0.05。

Lasso是一種高維預測回歸方法，并已被廣泛應用于高維數據生存分析的Cox比例風險回歸模型中[12]。為了進一步篩選出與LUAD生存顯著相關的基因，在訓練集(n= 209)中使用Lasso回歸模型對候選基因進行篩選，并進行十折交叉驗證，以確定最佳的預后相關基因。

1.3 多基因預后標志的建立和評估

對預后相關基因進行多因素Cox比例風險回歸分析，獲得預后相關基因的回歸系數。然后，采用predict函數將基因的表達水平和回歸系數進行組合算出每個患者的風險評分。

使用“survminer”包計算出最優cutoff值，以cutoff為臨界值，將訓練組分為高風險組和低風險組。為了確定風險評分在預測肺腺癌患者臨床預后中的作用，采用對數秩檢驗對訓練組進行了生存分析，比較高風險組和低風險組之間的生存差異。繪制了與時間相關的ROC曲線來進一步評估風險評分的預后性能，并計算了其3年和5年的AUC值。

此外，為了探討多基因預后標志在其他臨床參數中的診斷能力，進行了一項分層分析，以cutoff值為分界點進行分組，使用Kaplan-Meier曲線比較了stage亞組、年齡、性別亞組中高低風險組的生存差異。

1.4 多基因預后標志的驗證

使用內部驗證集(n=208)，外部驗證集GSE50081(n=127)，以及全集(n=417)來驗證風險評分的預測能力和適用性。在驗證集中，使用訓練集中獲得的回歸系數計算每個樣本的風險評分，然后根據cutoff值將患者分為高風險組和低風險組，采用對數秩檢驗進行生存分析，繪制與時間相關的ROC曲線。

1.5 列線圖的構建和驗證

對風險評分和一些臨床參數(stage、T期、N期、年齡、性別、復發)進行了單因素Cox回歸分析，以比較風險評分與臨床參數的預后能力。然后，使用多因素Cox回歸模型來確定風險評分是否具有臨床獨立性，其中，在單因素Cox回歸分析中具有顯著統計學差異(p<0.05)的臨床參數也被納入多因素Cox回歸模型中。

基于上述單因素和多因素Cox回歸分析，篩選出具有統計學差異的參數作為獨立預后參數，用于列線圖的構建，以預測患者3年、5年的生存概率。

為了評價模型的預測能力，計算出列線圖模型的C-index，并繪制其3年、5年的ROC曲線，同時繪制了3年時stage、風險評分和列線圖的ROC曲線，比較三者的預測能力。然后，使用校準曲線，通過500次重采樣，以3年、5年的觀察速率來可視化列線圖的性能，列線圖的預測結果和實際結果都能夠在校準曲線中進行比較，其中，45°線為最佳預測結果。在內部驗證集和全集中使用上述相同的辦法來驗證結果。

2 結 果

2.1 預后相關基因的篩選

在HADb數據庫、ARN數據庫和自噬數據庫中共下載了1 417個自噬相關基因，其中938個基因在TCGA數據中有表達。對938個自噬基因進行差異分析，獲得了38個上調基因和44個下調基因(見圖1a)，篩選條件為|logFC|>1.5，P<0.05。

在全集中，對差異基因進行單因素Cox回歸分析，發現有13個候選基因與肺腺癌生存相關(見圖1b)。為進一步確定與LUAD患者預后相關的基因，使用“glmnet”R包對候選基因進行了LASSO回歸分析以及十折交叉驗證，其結果顯示，當λmin=0.029時，模型性能達到最佳，此時篩選出了6個預后相關基因(見圖1c,1d)，即ARNTL2、NAPSA、ATG9B、CAPN12、MAP1LC3C、KRT81，這些基因中有4個(NAPSA、ATG9B、CAPN12、MAP1LC3C)的風險比小于1，表明它們的低表達與預后不良有關，而ARNTL2和KRT81的風險比大于1，表明它們的過度表達與低生存率有關。

圖1 回歸分析篩選與LUAD預后相關的自噬相關基因

2.2 多基因預后標志的建立和評估

對6個預后相關基因進行多因素Cox回歸分析(見圖2)，然后，使用predict函數結合多基因的回歸系數和表達量構建風險評分，通過“survminer”R包獲取風險評分的最優cutoff值，以cutoff值為分界點，將患者分為高風險組和低風險組，并展示了訓練集中患者的生存狀態和6個預后相關基因的熱圖(見圖3a)。對訓練組進行生存分析，結果顯示，與低風險組相比，高風險組的預后結果更差(見圖3b)。然后，我們構建了一個與時間相關的ROC曲線(見圖3c)，其3年、5年的AUC值分別為0.852、0.868，這表明這個多基因預后標志具有較好的預測能力。

圖2 預后相關基因的多因素Cox回歸分析

圖3 訓練集中多基因特征的預后分析

此外，對stage、年齡和性別進行了風險分層，以cutoff值為分界點，將訓練組的患者分為高風險組和低風險組，進行Kaplan-Meier生存分析(見圖4)。在stage Ⅰ/Ⅱ、stageⅢ/Ⅳ、男性、女性、年齡大于65歲和年齡小于65歲的亞組中，高風險組的生存率都顯著低于低風險組(p<0.05)。

圖4 風險評分在不同亞組中的生存分析

2.3 多基因預后標志的驗證

使用內部測試集(n=208)、外部測試集(n=127)和全集(n=417)來驗證風險評分的預測能力。與訓練集中的結果一致，測試集的生存分析曲線都顯示高風險組的預后結果比低風險組的差(見圖5a-5c)。ROC曲線顯示，內部測試集的3年、5年AUC值為0.863、0.938(見圖5d)，外部測試集的3年、5年AUC值為0.939、0.852(見圖5e)，全集的3年、5年AUC值為0.861、0.905(見圖5f)，以上結果都顯示風險評分在預測LUAD患者的預后方面表現良好。

圖5 風險評分的內部驗證和外部驗證

2.4 列線圖的構建與驗證

對風險評分和一些臨床參數(stage、T期、N期、年齡、性別、復發)進行了單因素和多因素Cox比例風險回歸分析，其結果顯示風險評分可以作為預測LUAD預后的獨立參數，而在傳統臨床參數中，stage和復發也可以作為獨立預后參數(見圖6a,6b)。我們將傳統臨床風險因素和風險評分相結合，構建一種能夠有效預測患者3年、5年生存率的列線圖(見圖6c)。列線圖的C-index指數為0.807，表明列線圖有較好的預測能力。校準曲線顯示，列線圖的預測結果與實際結果較為一致(見圖7a)。ROC曲線顯示，列線圖3年、5年的AUC值分別為0.898、0.88(見圖7d)。三年時，列線圖生存的AUC值遠高于風險評分模型和stage的AUC值(見圖7g)，這表明列線圖可能是預測LUAD預后生存的最佳方式。

圖6 臨床單、多因素Cox分析以及列線圖的構建

為了驗證列線圖的預測價值，使用內部測試集(n=208)和全集(n=417)來檢驗上述的發現。內部測試集和全集的列線圖的C-Index指數分別為0.8和0.792，校準曲線也顯示兩個測試集列線圖的3年、5年生存預測結果與實際結果有良好的一致性(見圖7b,7c)。列線圖的ROC曲線顯示，兩個測試集具有較好的預測準確度(見圖7e,7f)。同時，在3年期的生存預測中，列線圖無論在哪組都比風險評分和stage有更好的預測準確度(見圖7h,7i)。

圖7 列線圖預測LUAD生存率的性能以及列線圖、風險評分和Stage預測能力的比較

3 討 論

自噬是高度保守的代謝過程，在循環代謝能量以維持細胞內穩態方面起著關鍵作用[13]。有研究表明了多個自噬相關基因與肺癌的發生發展密切相關[14-16]，因此，決定把自噬相關基因作為肺腺癌治療的潛在靶點。通過對TCGA肺腺癌數據中的938個自噬相關基因進行差異分析，獲得了82個差異基因，然后對差異自噬基因進行單因素Cox回歸分析，篩選出了13個與LUAD生存相關的候選基因，然后使用lasso回歸進一步篩選出6個與LUAD預后相關的基因。通過多因素Cox回歸分析獲得每個預后相關基因的回歸系數，通過每個基因的表達量和回歸系數計算出每個患者的風險評分。在訓練集中，風險評分能夠很好地將高風險患者和低風險患者區分開，并且其預測性能也在內部、外部測試集中得到了驗證。同時，在分層分析中，風險評分在stage，年齡和性別亞組中的風險分層表現也很好，這意味著此風險評分模型可以根據亞組將LUAD患者分為高低風險組，幫助臨床醫生進行臨床決策。

用于構建風險評分的6個基因包括ARNTL2、NAPSA、ATG9B、CAPN12、MAP1LC3C和KRT81。ARNTL2屬于PAS超家族，在晝夜節律和缺氧過程中起著重要的作用，其在乳腺癌、腎細胞癌等人類惡性腫瘤中具有致癌作用[17-19]，目前已有研究報道ARNTL2的高表達與肺腺癌的低生存期相關，并且能夠影響肺腺癌的免疫浸潤水平[20-21]。NAPSA是天冬氨酸肽酶，其編譯的蛋白酶能夠參與肺表面活性物質蛋白B在肺中的蛋白水解過程，目前它已被證實是肺腺癌的生物標記物，并且已被用作識別原發性肺腺癌的免疫組化染色劑[22-24]。ATG9B是自噬相關基因，在自噬過程中起調節作用，與肝癌[25]、腎細胞癌[26]、胃癌[27]等多種癌癥的發生發展有關，但其在肺腺癌中的作用還尚未闡明。CAPN12是一種鈣蛋白酶，鈣蛋白酶能夠調節多種細胞生理過程，包括細胞增殖、細胞遷移、細胞侵襲、細胞自噬等，各種癌癥的發病機制也需要鈣蛋白酶系統，其可能起到促進癌癥發展的作用，最新研究也確定了CAPN12是新的結直腸癌易感基因[28-30]。MAP1LC3C是自噬蛋白ATG8的同源物，被用作自噬機制的生物標志物，有研究發現，MAP1LC3C介導了MET/HGF-RTK信號通道在癌癥中的作用，MAP1LC3C和MET復合物招募HGF并且激活MET-RTK信號通路從而進行自噬降解，進而影響腫瘤轉移[31]。在肺癌方面，有研究證明其與肺腺癌氧化磷酸化過程十分相關[32]。KRT81是一種角蛋白，相關研究發現其與肺腺癌腫瘤轉移相關[33]。上述結果都表明風險評分模型具有潛在的臨床應用價值。

最后，單因素和多因素Cox回歸分析表明，風險評分可以作為預后評估的獨立因素。為了提高風險評分的預測能力，結合臨床參數和風險評分構建了一個基于多基因預后標志的列線圖來預測患者生存率。通過比較，列線圖的預測性能高于風險評分和stage的預測性能，并在驗證集中得到同樣的結果。這表明，與單一的臨床參數相比，列線圖模型更能幫助臨床醫生預測LUAD患者的生存狀態，并為臨床醫生提供治療指導。然而，我們的研究還有一些不足之處，我們的數據只包含TCGA數據庫的mRNA數據，未來還可以從單核苷酸多態性、拷貝數變異數據、DNA甲基化等突變數據中進一步分析這6種新的生物標志物是否與上述突變相關。

4 結 論

基于6個基因的多基因預后標志來預測LUAD患者的生存風險，在訓練集和測試集中都表現出良好的準確率，并且獨立于其他臨床特征。然后，結合多基因預后標志和臨床特征構建了列線圖模型以預測LUAD患者的預后生存率，與單一臨床特征相比，列線圖模型具有更好的預測性能。因此，這6個基因很可能是LUAD的潛在生物標志物，基于多基因預后標志和臨床特征的列線圖模型很有可能用于評估LUAD患者的生存率，并幫助臨床醫生進行個體化治療的臨床決策。