基于TCGA和ICGC構建和驗證肝細胞癌腫瘤微環境相關基因的預后風險模型

李紹智 梁維仁 郭金友 鄭家平

肝細胞癌（hepatocellular carcinoma，HCC）是全球最常見的惡性腫瘤之一，其發病率占惡性腫瘤的第六位，病死率占第四位［1］。近年來有較多治療方式，包括手術治療和藥物治療等已經取得巨大進展，但HCC患者的預后仍不理想，其中主要原因是腫瘤的復發和轉移［2］。越來越多的研究表明腫瘤微環境（tumor microenvironment，TME）在腫瘤的發生、發展等過程中扮演著重要的角色［3-4］。腫瘤細胞可以通過分泌各種細胞因子、趨化因子等，改變TME并導致周圍細胞重新編程，使這些細胞能夠促進或抑制腫瘤細胞的生存和發展［3，5］。研究報道TME相關基因在預測不同癌種患者的生存、預后具有重要價值［6-8］。在肝癌上，盡管一些研究已經構建了TME相關基因的預后風險模型，但這些研究缺乏獨立數據集的驗證，并且預測的性能仍需提高［9-11］。生物信息學方法在探索和鑒定腫瘤預后相關生物標志物上具有重要的價值［12］。本研究基于TCGA和ICGC數據庫篩選TME相關的關鍵基因，構建和驗證預后風險評估模型，以期為HCC患者的預后評估提供參考。

1 材料與方法

1.1 數據下載和預處理 HCC轉錄組數據、臨床信息、隨訪資料和體細胞突變信息從癌癥基因圖譜（the cancer genome atlas，TCGA）數據庫（https：//portal.gdc.cancer.gov/）中獲取，其中對隨訪信息有缺失的患者予以刪除。其次，獨立數據集國際腫瘤基因組協作組（international cancer genome consortium，ICGC）中的HCC隊列（LIRIJP）患者的轉錄組數據也被下載用于驗證。根據先前已經發表的TME相關研究［13-19］，4,061個TME相關的基因被選擇用于后續分析。

1.2 非負矩陣分解算法驗證分子亞型 通過單因素Cox分析對TCGA中這4,061個TME基因表達水平進行分析，獲得與預后相關的TME基因（P＜0.01）。通過非負矩陣分解（non-negative matrix factorization，NMF）對樣本進行聚類［20］，聚類的數量k根據cophenetic和殘差平方和確定。最終選擇合適的聚類數量用于后續分析。為了比較不同聚類間微環境細胞群的差異，該方法［19］被用于計算每個樣本中各種免疫細胞的評分，包括CD8+T細胞、細胞毒性淋巴細胞、B細胞、中性粒細胞、單核細胞、骨髓樹突狀細胞等。不同免疫細胞在不同聚類之間的差異也被比較。

1.3 預測模型構建 對臨床數據整理過濾后，370例HCC患者（TCGA數據集）以3∶2的比例被隨機分為訓練組（Training cohort）和測試組（Testing cohort），并比較兩組間臨床信息以確定是否存在組內差異，進一步通過Lasso、Cox回歸尋找合適的預測相關TME基因［21］。Lasso回歸能夠消除過于相關的基因，防止過度擬合，并能找到誤差最小的基因。基于Lasso回歸的結果，多因素Cox回歸分析被用于進一步篩選重要的與HCC患者預后密切相關的TME樞紐基因。根據TME相關預后基因的表達及其多因素分析產生的回歸系數生成最終的預后特征：TME風險評分（TMEScore）=（βmRNA1×mRNA1的表達水平）+（βmRNA2×mRNA2的表達水平）+…+（βmRNAn×mRNAn的表達水平）。TMEScore的cutoff值是根據訓練集中的中位TMEScore決定。根據TMEScore的cutoff值，將患者分為高風險組或低風險組。

1.4 預測模型性能評估 基于建立的TMEScore，進一步比較訓練集、測試集、整個人群以及獨立數據集（ICGC-JP）中人群的預后差異以及預測模型對于1、3、5年預后評估的預測性能。單因素、多因素Cox回歸分析用于研究TMEScore在整個人群中的預后價值，其中在單因素分析模型中具有顯著差異的參數被納入多因素分析模型中進一步驗證。亞組分析用來明確TMEScore在不同的亞組人群中的預后價值。值得一提的是，該模型的預測性能也與已發表的預后模型相比較［9-11］。

1.5 列線圖構建和評估 列線圖是一種通過標尺的長度來衡量不同變量或者特征對結果的影響，是一種直觀的風險模型可視化方法。因此，應用TMEScore整合臨床變量構建預后預測的列線圖用于評估1年、3年和5年的總體生存率，主要臨床特征包括年齡、性別、分級、分期和TMEScore。校準曲線用于檢驗列線圖的擬合的準確性。

1.6 富集通路分析 為了明確高、低風險評分組的富集通路差異，從MSigDB數據庫中下載并獲得的“c2.cp.kegg.v7.4.symbols”基因集用于GSEA分析，并比較兩個風險組之間的富集通路評分差異，并用熱圖可視化。使用R包 “limma”“org.Hs.eg.db”“clusterProfiler”和“enrichplot”來確定高、低風險評分組之間的富集通路差異。

1.7 統計學方法 采用R（3.6.1版）和Strawberry Perl（5.32.0.1版）軟件。計量資料用t檢驗，計數資料用卡方檢驗或Fisher確切概率法。預后分析采用Kaplan-Meier曲線。P＜0.05為差異有統計學意義。

2 結果

2.1 鑒定HCC不同的分子亞型 在數據預處理后，370例和240例HCC樣本的轉錄組芯片分別從TCGA和ICGC數據庫中獲得。通過單因素Cox分析對TCGA隊列中的TME相關基因表達分析，共有528/4,061個基因被篩選為TME相關的預后基因。通過NMF分析將所有患者分成C1和C2兩個亞群（見圖1A），結果發現患者的OS和PFS在兩個亞群中有明顯的差異，C1亞群患者的OS和PFS高于C2亞群的患者（見圖1B和1C）。進一步比較兩個亞群之間的免疫細胞浸潤評分，發現CD8+T細胞、上皮細胞、單核細胞、骨髓樹突狀細胞、NK細胞和T細胞在兩個亞群之間有明顯的差異（見圖1D）。聚類分析也發現C1和C2兩個亞群分別匯聚與不同的免疫亞型（見圖1E）。

圖1 兩個聚類（C1和C2）生存、預后及免疫評分和亞型的比較。A. 通過非負矩陣分解算法進行聚類。B～C. 兩個聚類之間的總生存期（Overall survival，OS）和無進展生存期（Progression-free survival，PFS）。D. 腫瘤微環境中的不同免疫細胞評分在兩個聚類間的比較。E. 免疫亞型在兩個聚類間的比例

2.2 基于TME相關基因構建預后風險評分 為了構建TME相關的預后風險評分模型，隨機將TCGA隊列分成訓練集（224例）和測試集（146例）兩組，兩組間的臨床基線資料差異無統計學意義（見表1）。為了驗證與預后相關的TME基因，對訓練集的224例患者行單因素Cox分析，共得到488個預后相關的TME基因。為了進一步縮小TME基因的范圍，Lasso回歸和Akaike信息準則也被用于篩選TME相關基因。最終，2個關鍵的TME預后基因（CDCA8和CBX2）被確認，TMEScore的風險評分模型為CDCA8的表達水平×0.3533+CBX2的表達水平×0.5156。不同TMEScore對患者的預后影響以及該評分模型對預后的診斷性能被進一步評估。單因素和多因素分析均揭示TMEScore是預后相關的風險因素（見表2）。生存分析顯示基于TMEScore的高、低風險分組能夠很好區分訓練集（見圖2A）、測試集（見圖2B）、整個數據集（見圖2C）以及ICGC（見圖2D）中的人群預后，且高TMEScore患者比低TMEScore患者的預后更差。并且，TMEScore預測1、3和5年患者生存的AUC在訓練集中為0.818、0.721和0.699（見圖2E）；在測試集中分別為0.671、0.652和0.633（見圖2F）；在整個TCGA中分別為0.762、0.704和0.664（見圖2G）。在ICGC數據集中，TMEScore預測1、3和5年患者生存的AUC為0.770、0.783和0.390（見圖2H）。該預測評分在ICGC數據集中的5年AUC較低的主要是由于其中位隨訪時間相對較短，較少的事件數影響了模型預測的效果。為了進一步分析TMEScore在不同臨床特征的人群中的預測價值，亞組分析在TCGA的患者中開展。分層分析結果顯示TMEScore在年齡＞65歲、≤65歲、男性和女性、G1～2期、G3～4期、M0、N0、I～II期、III～IV期、T1～2期、T3～4期人群中均顯示出預后預測價值（見圖2I-2T）。基于TMEScore風險評分在T1期患者及其他T分期患者中差異有統計學意義，提示該預后評分對區分早期和晚期患者具有一定的價值（見圖2U）。

表1 訓練集和測試集患者的臨床基線特征［n（%）］

表2 單因素和多因素分析預后相關風險特征

圖2 基于TMEScore的風險評分的預測效能評估

2.3 整合臨床參數的列線圖預測預后 通過整合臨床參數和TMEScore，一個預測患者1、3和5年生存評估列線圖被建立（見圖3A）。擬合曲線提示該列線圖具有良好的擬合性（見圖3B）。進一步比較構建的列線圖與其他臨床參數的預測性能，結果提示該列線圖預測患者1年生存的AUC比其他臨床參數高（AUC=0.777）（見圖3C）。決策曲線分析也提示構建的列線圖的預測性能比其他臨床參數更好（見圖3D）。為進一步驗證構建的TMEScore在預測患者預后上的價值，本研究搜集了已發表的與TME相關的風險模型，并比較這些模型的預測性能。相比于其他3個已發表的TME相關預測模型，本研究的模型在1、3和5年的AUC上均高于這三個模型（見圖3E）。生存分析提示Deng等、Tian等及本研究構建的模型能夠良好區分高、低風險人群的預后（見圖3F）。以上結果也提示本研究構建的TMEScore以及列線圖模型能夠有效評估患者生存預后風險，對具有高死亡風險的患者進行有效的監督隨訪和必要的干預對提升這些患者的生存具有重要意義。

2.4 富集通路分析和相關性分析 進一步對高、低風險評分組內的樣本行基因富集分析以明確高、低風險評分組之間涉及的信號通路差異，結果提示高風險評分組與KEGG_CELL_CYCLE，KEGG_DNA_REPLICATION，KEGG_HOMOLOGOUS_RECOMBINATION，KEGG_OOCYTE_MEIOSIS和KEGG_PRIMARY_IMMUNODEFICIENCY通路相關（見圖4A）。低風險評分組與KEGG_COMPLEMENT_AND_COAGULATION_CASCADES，KEGG_DRUG_METABOLISM_CYTOCHROME_P450，KEGG_FATTY_ACID_METABOLISM，KEGG_PPAR_SIGNALING_PATHWAY和KEGG_RETINOL_METABOLISM通路相關（見圖4B）。進一步研究發現TMEScore與免疫檢查點（PDCD1、CD274、CTLA4）、DNA復 制（POLE2、FEN1、MCM6、POLD3）、錯配修復（MSH6、MSH2）和上皮-間質轉化（FAP、LOXL2）均呈正相關（見圖4C）。由于TMEScore與免疫檢查點成正相關，進一步分析其與免疫細胞之間的相關性。結果顯示TMEScore與T細胞、CD8+T細胞、細胞毒性淋巴細胞、單核細胞、髓系樹突狀細胞、纖維母細胞呈正相關，與中性粒細胞呈負相關（見圖4D）。

圖3 列線圖構建可視化風險評分模型。A. 構建整合TMEScore的風險值及臨床參數的列線圖。B. 患者1年、3年、5年的擬合曲線。C. 不同臨床參數之間預測患者1年生存的AUC比較。D. 決策曲線評估列線圖的臨床獲益。E. 基于TMEScore的風險評分模型預測性能與先前發表的TME相關預測模型的比較。F. 基于TMEScore的風險評分模型與先前已發表的TME相關預測模型對患者的預后評估比較

圖4 基于TMEScore的風險評分值在樣本中的分子機制。A-B. 高、低風險評分組中的樣本參與的重要通路。C. HCC樣本風險評分值與腫瘤學重要通路的靶點的代表性基因表達之間的相關性。D. HCC樣本風險評分值與免疫細胞評分之間的相關性

3 討論

大量文獻表明TME在HCC的病理生成及復發、轉移中起關鍵作用，這表明HCC的生物學行為不僅與腫瘤細胞有關，也與TME密切相關［22-23］。因此，為了明確HCC中與TME相關的預后生物標志物，本研究根據先前發表的研究整理和過濾了4,061個TME相關的基因。通過Lasso Cox回歸結果確定了2個關鍵的TME預后基因（CDCA8和CBX2）。其中，CDCA8是染色體乘客復合體（Chromosomal passenger complex，CPC）的一個重要組成部分，與有絲分裂過程中的細胞動態定位調節有關。JOEN等［24］報道，沉默CDCA8能夠恢復ATF3腫瘤抑制因子的功能和以及失活致癌信號通路AKT/β-catenin，從而抑制HCC的生長和干性，提示靶向CDCA8可能是治療原發性HCC復發、轉移的潛在靶點。CUI等［25］發現高表達的CDCA8與HCC患者更差的總生存率（P=0.0054）和無進展生存（P=0.0009）相關。體外實驗顯示抑制CDCA8能減緩細胞增殖，阻斷G0/G1階段的細胞周期。體內實驗表明抑制CDCA8能抑制腫瘤的生長。這些發現也提示CDCA8可能是HCC的一個有效治療靶點。另一個關鍵基因CBX2，是核心蛋白復合物1（polycomb repressive complex 1，PRC1）的一個重要組成部分。MAO等［26］發現高表達的CBX2與HCC患者的不良預后有關，且敲除CBX2能抑制HCC細胞的增值和促進HCC細胞的凋亡，提示CBX2是一個可能的治療靶標。

基于這兩個關鍵的預測子，本研究構建TMEScore并發現高TMEScore與患者更差的預后有關，且在獨立數據集和亞組人群中也驗證了其預后價值，提示該風險評分特征能作為預測患者預后的獨立因素。基于TMEScore構建的列線圖能很好評估患者的預后風險。盡管先前研究也報道TME相關基因在預后預測上的價值［9-11］，但通過比較后發現本研究的模型在預測性能上高于已發表的模型。這些結果也驗證了本研究的風險模型在患者預后預測上的價值。

TMEScore與免疫檢查點（PDCD1、CD274、CTLA4）、DNA復制（POLE2、FEN1、MCM6、POLD3）、錯配修復（MSH6、MSH2）和上皮-間質轉化（FAP、LOXL2）均呈正相關。這與GSEA富集分析的結果相一致。這些重要的機制也已經被報道與HCC的發生、發展及療效有關［27］。免疫相關分析也揭示TMEScore與多種免疫細胞相關，提示該風險評分可能與腫瘤浸潤免疫細胞（tumor-infiltrating lymphocytes，TILs）水平相關。TILs水平也被認為是重要的免疫治療生物標志物，能夠反映免疫治療療效［28］。因此，本研究構建的風險評分模型可能對預測免疫治療療效也有一定的價值，但仍需要后續的研究。

本研究通過生物信息學的方法構建和驗證了一個TME相關的雙基因預后風險評估模型（CDCA8和CBX2）。該模型能夠有效預測HCC患者預后風險水平，對于具有高死亡風險的人群，臨床醫師能夠采取必要的監督隨訪和合理的干預措施提升患者的生存。