肝癌預后模型的構(gòu)建及其在免疫治療中的作用研究

李一鳴，陳文翔，王 禹，何 俊，周劍寅

（廈門大學附屬中山醫(yī)院肝膽胰外科，福建 廈門 361005）

肝細胞癌（hepatocellular carcinoma）惡性程度高，手術(shù)是其根治療法，但大部分患者發(fā)現(xiàn)時已處于晚期，失去接受根治性手術(shù)的機會[1]。因此，非手術(shù)治療方法成為近年來肝癌治療的研究熱點。免疫檢查點抑制劑（ICI）為肝癌的非手術(shù)治療帶來突破進展[2]，但目前相關(guān)藥物的臨床應用仍不足，藥物并發(fā)癥難以避免[3]，藥物療效不穩(wěn)定[4]。肝細胞癌具有多種免疫表型，其高度異質(zhì)性可能是ICI 效果不理想的根本原因[5，6]。腫瘤突變負荷（TMB）指編碼區(qū)內(nèi)的體細胞突變數(shù)目，通常以每百萬堿基中體細胞突變的總數(shù)來表示。TMB 可作為免疫療法的生物標志物，TMB 高的患者更有可能在免疫治療中獲益[7，8]。基因組不穩(wěn)定是癌癥的標志之一，對預后有重要影響[9-11]。有研究利用lncRNAs 建立風險模型來評價乳腺癌基因組不穩(wěn)定性，展現(xiàn)出分子特征測量基因組不穩(wěn)定性的潛在能力[12，13]。長非編碼RNA（lncRNAs）是一類長度超過200 個的核苷酸，且不具有蛋白質(zhì)編碼潛力的非編碼RNA[14]。有研究顯示，lncRNAs 在維持基因組不穩(wěn)定方面具有關(guān)鍵作用[15]。然而，基因組不穩(wěn)定性相關(guān)lncRNAs（genome instability-associated lncRNAs，GI-lncRNAs）與HCC的關(guān)系報道較少。本研究試圖建立基于GI-lncRNAs 的肝癌預后模型，以期為肝癌相關(guān)的研究及治療提供新的思路。

1 資料與方法

1.1 數(shù)據(jù)獲取 從腫瘤基因組圖譜數(shù)據(jù)庫（TCGA）下載TCGA-LIHC 項目數(shù)據(jù)（T：375、N：50），包含其RNA-seq 數(shù)據(jù)、臨床數(shù)據(jù)及體細胞突變數(shù)據(jù)（https://portal.gdc.cancer.gov/）；下載腫瘤免疫細胞浸潤性數(shù)據(jù)（http://timer.cistrome.org）和TCGA 所有腫瘤細胞免疫分型數(shù)據(jù)[16]（https://doi.org/10.1016/j.immuni.2018.03.023）。整合來自Ensemble 網(wǎng)站（http://asia.ensembl.org）的人類基因注釋文件，用于區(qū)分TCGA數(shù)據(jù)庫中獲得的基因，篩選出所需要長鏈非編碼基因lncRNAs。從GSEA 官網(wǎng)（https://www.r-project.org/）下載KEGG、GO 相關(guān)基因集文件。

1.2 GI-lncRNAs 的鑒定和驗證 處理突變數(shù)據(jù)得到TCGA 每個樣品中的累計體細胞突變頻率（n=364），根據(jù)累計體細胞突變頻率進行排序，前25%樣品定義為高突變組（HG：93 個），后25%樣品定義為低突變組（LG：90 個）。比較HG 組和LG 組之間lncRNAs表達差異，篩選差異lncRNAs 的標準為|logFc|≥1，F(xiàn)DR 值≤0.05，對于符合條件的lncRNAs 稱為GIlncRNAs。使用hclust 函數(shù)，采用最大距離法根據(jù)GI-lncRNAs 的表達水平對腫瘤樣品（374 例）進行層次聚類分析，聚類結(jié)果組間做累計體細胞突變頻率差異分析。

1.3 GI-lncRNAs 預后模型（GIlncsig）的構(gòu)建和驗證排除TCGA-LIHC 臨床數(shù)據(jù)中隨訪時間＜30 d 及生存數(shù)據(jù)缺失的患者，將剩余343 例患者（All 組）隨機分為Train 組和Test 組。Train 組用于模型構(gòu)建，Test 組和All 組用于驗證模型的預測能力。在Train組中，首先對GI-lncRNAs 進行單因素Cox 回歸分析，篩選與預后相關(guān)的GI-lncRNAs，再通過多因素Cox 回歸分析構(gòu)建評分模型。根據(jù)多因素Cox 分析的回歸系數(shù)和每個lncRNA 表達水平，計算每個患者的風險評分。風險評分計算公式如下：Risk score=，其中coef 為多變量Cox 回歸的回歸系數(shù)，Exp 為相對應的基因表達量。根據(jù)風險評分中位值區(qū)分風險亞組，分別評估風險亞組生存期、GIlncsig 預測能力及GIlncsig 與其他臨床因素（性別、年齡、組織學分級、病理分級等）的關(guān)系。

1.4 GIlncsig 亞組分子免疫特點、ICI 治療及化學治療的綜合分析 對GIlncsig 亞組行基因突變分析與GSEA 富集分析，評估其分子及相關(guān)功能差異。利用腫瘤免疫細胞浸潤性數(shù)據(jù)（Timer、CiberSort、XCell、QUANTISEQ、MCPCounter、EPIC 和CiberSort）和ICI相關(guān)生物標志物[17]進行差異分析，分析GIlncsig 與免疫細胞的關(guān)系，評估GIlncsig 對于免疫治療的預測潛力。使用“pRRophetic”包[18]，利用TCGA-LIHC項目中常用化療藥物的IC50值（半抑制濃度），評估GIlncsig 在肝癌治療中的應用價值。

1.5 統(tǒng)計學方法 本研究所有統(tǒng)計分析均采用R-4.1.1 進行。采用生存曲線（Kaplan-Meier）、對數(shù)秩檢驗（log-rank test）評估患者的生存差異，受試者工作特征（ROC）曲線分析評估模型預測能力，秩和檢驗評價相關(guān)性。單因素Cox 回歸分析、多因素Cox 回歸分析、分層分析、χ2檢驗評估GILncSig 與其他臨床因素的關(guān)系。

2 結(jié)果

2.1 GI-lncRNAs 的鑒定和驗證 共鑒定出88 個GIlncRNAs。在HG 組，32 個GI-lncRNAs 表達下調(diào)，56個GI-lncRNAs 表達上調(diào)。層次聚類分析結(jié)果顯示，所有腫瘤樣本可聚類分為LG-like 組及HG-like組，差異分析顯示：兩組間體細胞突變模式存在差異，HG-like 組體細胞累積突變中位值更高。

2.2 GIlncsig 的開發(fā) 將343 例肝癌患者按1∶1 隨機分為Test 組或Train 組，兩組間臨床特征比較，差異無統(tǒng)計學意義，見表1。Train 組單因素Cox 回歸分析顯示有17 個GI-lncRNAs 與患者預后相關(guān)（P＜0.05），將這些基因納入多因素Cox 回歸分析，共得到 4 個 基 因（MIR210HG、LUCAT1、AC245041.2、AC010643.1）可用于建立肝癌GI-lncRNAs 預后模型，各基因系數(shù)見表2。以Train 組風險評分中位數(shù)（0.809 185 8）為臨界值將Train 組分為高風險組（High-risk）和低風險組（Low-risk）。HR＞1 的基因為危險基因。根據(jù)表2，上述lncRNAs 均為危險基因，其高表達與預后不良有關(guān)。基因表達熱圖顯示隨著風險評分的增加，lncRNAs 表達增加（圖1A）。Kaplan-Meier 生存曲線顯示，低風險組患者生存時間明顯大于高風險組（圖1B）。對GIlncsig 的ROC 曲線分析得到1、3、5 年的曲線下面積（AUC）分別為0.741、0.687、0.659（圖1C）。進一步評估GIlncsig 的預后價值，將風險評分和臨床特征比較：單變量分析顯示，風險評分與總生存期（OS）密切相關(guān)；多變量分析顯示，風險評分是患者的獨立預后指標，見表3。

表3 train 組中GIlncsig 及臨床特征與OS 的相關(guān)性分析

圖1 train 組GIlncsig 的驗證評估

表1 肝細胞癌患者臨床信息統(tǒng)計表[n（%）]

表1（續(xù)）

表2 train 組中多因素Cox 回歸分析建立風險評分模型

2.3 GIlncsig 在Test 及All 的驗證 以Train 組風險評分中位數(shù)將Test 組的患者分為高、低風險組。在Test 組及All 組中，采用上述方法進行驗證，基因表達熱圖顯示隨著風險評分的增加，lncRNAs 表達增加；Kaplan-Meier 生存曲線顯示，低風險組生存時間大于高風險組；單變量分析顯示，GIlncsig 與OS 顯著相關(guān)；多變量分析顯示，GIlncsig 是患者的獨立預后因子。

2.4 進一步探究GIlncsig 的臨床價值 在All 組，多指標ROC 生存曲線顯示，風險評分AUC 值為0.760，大于其他臨床指標（以3 年為例，所有樣本平均生存時間2.8 年），見圖2A。除風險評分外，多因素Cox 分析顯示，病理分期（stage）也有獨立預后因子。因此，分層分析以確定風險評分是否具有獨立于病理分期的預后價值。根據(jù)患者病理分期，將All 組中病理階段為Ⅰ或Ⅱ的患者合并為早期組（n=238），病理階段為Ⅲ或Ⅳ的患者合并為晚期組（n=83），利用GIlncsig 區(qū)分不同病理分期組的高低風險評分，結(jié)果顯示，早期組中，風險亞組間總生存率比較，差異有統(tǒng)計學意義（P＜0.05）；而晚期組中，風險亞組間總生存率比較，差異無統(tǒng)計學意義（P＞0.05）。進一步分析肝細胞癌風險亞組與臨床特征的關(guān)系，風險亞組與臨床特征的相關(guān)性熱圖顯示，GIlncsig 與腫瘤病理分期（stage）、分級（grade）、大小（T）均密切相關(guān)，見圖2B，進一步探究GIlncsig 對于以上臨床特征的內(nèi)部區(qū)分能力，結(jié)果間圖2C、圖2D。

2.5 與其他lncRNAs 風險模型比較 將GIlncsig 的預測能力與最近發(fā)表的2 個lncRNAs 模型比較，來源于Liao L 等[19]的研究（以下簡稱liaolncsig）和Gu X 等[20]的研究（以下簡稱Gulncsig），研究均基于同一個TCGA 患者隊列。結(jié)果顯示，GIlncsig 的3 年AUC 為0.744，高于liaolncsig（AUC=0.635）和Gulncsig（AUC=0.694），見圖2E。

圖2 All 組中GIlncsig 的臨床價值

2.6 風險亞組的分子特征 深入研究高低亞組間體細胞突變差異，繪制風險亞組基因突變瀑布圖，結(jié)果顯示，高風險組突變數(shù)高于低風險組；探究基因突變對患者生存的影響，結(jié)果顯示有6 個基因突變（DOCK2、LRP1B、TP53、ARID1A、NPAP1、MUC5B）與患者OS 相關(guān)，結(jié)合風險亞組突變情況（低風險組中14%突變：高風險組中42%突變），TP53 基因格外特別。探究GIlncsig 臨床意義是否優(yōu)于TP53 突變狀態(tài)，根據(jù)GIlncsig 與TP53 突變狀態(tài)繪制生存曲線，結(jié)果顯示不同組間生存狀況存在差異。TP53 wild/low 組比TP53 mut/high 組的患者有更好的生存結(jié)局，TP53 wild/high 組生存結(jié)局最差，TP53 mut/low組（n=26）生存結(jié)局最好，見圖3。

圖3 風險亞組間基因突變分析

2.7 GSEA 富集分析探究風險亞組中富集的基因集以KEGG 基因集進行通路富集分析，高風險組樣本基因組富集于癌癥、細胞周期、DNA 復制、神經(jīng)活性配體受體相互作用及ECM 受體相互作用等相關(guān)途徑。低風險組樣本基因組富集于視黃醇的代謝、補體系統(tǒng)、脂肪酸代謝、氨基酸代謝等相關(guān)途徑；以GO基因集進行功能富集分析，高風險組樣本基因組富集于染色體分離、減數(shù)分裂細胞周期、微管細胞骨架組織及細胞器裂變相關(guān)功能，低風險組樣本基因組富集于一元羧酸分解過程、脂肪酸分解代謝過程、補體激活、氨基酸分解代謝等相關(guān)功能。

2.8 風險亞組間的免疫特征 GIlncsig 亞組能夠區(qū)分免疫亞型，進一步分析腫瘤浸潤性免疫細胞的差異。風險評分與T cell CD8+naive、T cell CD4+effector memory、T cell CD8+central memory 等呈負相關(guān)，與B cell、T cell CD4+memory、Common lymphoid progenitor 等呈正相關(guān)。探究GIlncsig 與ICI 相關(guān)生物標志物的關(guān)系：選擇CD274、CTLA4、HAVCR2、I DO1、LAG3 和PDCD1 作為免疫檢查點相關(guān)特征，選擇CD8A、CXCL10、CXCL9、GZMA、GZMB、IFNG、PRF1、TBX2 和TNF 作為免疫激活相關(guān)特征[17]。結(jié)果顯示，高風險評分與CD274、CTLA4、HAVCR2、IFNG、PDCD1、TNF 的高表達呈正相關(guān)，與其他基因無相關(guān)性，見圖4。

圖4 風險亞組間的免疫特征分析

2.9 GIlncsig 與化療藥物的相關(guān)性分析 分析常見化療藥物療效（包括阿霉素、絲裂霉素C、順鉑、長春花堿、索拉非尼）與GIlncsig 的相關(guān)性，結(jié)果顯示，高風險組中，化療藥物阿霉素、絲裂霉素C 的IC50值更低，見表4。

表4 GIlncsig 與化療藥物的相關(guān)性分析

3 討論

目前認為，具有高TMB 的腫瘤細胞具有較高的新抗原水平，理論上，TMB 越高，能被T 細胞識別的新抗原越多，接受ICI 的療效就越好[7，8]。

lncRNAs 是腫瘤生物學的重要組成部分，有望作為患者的預后標志物。因此，本研究利用TCGA 數(shù)據(jù)庫中肝細胞癌lncRNAs 表達和體細胞突變數(shù)據(jù)，建立GIlncsig（由MIR210HG、LUCAT1、AC245041.2、AC010643.1 構(gòu)成）。GIlncsig 將患者分為兩個存活率不同的風險組，在區(qū)別樣品體細胞突變頻率、臨床病理分級分期方面具有一定價值。進一步發(fā)現(xiàn)，GIlncsig 對患者預后的預測能力優(yōu)于其他臨床指標，并且是肝細胞癌患者的獨立預后因子。與其他有關(guān)lncRNA 預后模型比較，GIlncsig 也展示出更高的穩(wěn)定性。本研究中，肝細胞癌患者在病理分期晚期組、風險亞組間總生存率差異不明顯，可能原因是GIlncsig 對于晚期患者的區(qū)別能力不強，也不能排除因晚期組人數(shù)較少而引起的誤差。另外，本研究中TP53 mut/low 組患者的生存結(jié)局最好，推測因為人數(shù)較少而造成的誤差所致。可以推測，GIlncsig 比TP53 突變狀態(tài)具有更大的預后意義。

GI-lncRNAs 是基于累計體細胞突變頻率而篩選得到，故GIlncsig 與TMB 有關(guān)。有證據(jù)顯示，TMB對于免疫治療有預測價值，可能與腫瘤免疫微環(huán)境有關(guān)。在探究風險亞組間免疫學性質(zhì)中，基因突變結(jié)果顯示，TP53 突變在兩組患者間存在差異，其中高突變組患者預后更差，與生存結(jié)果一致。TP53 突變在腫瘤中具有普遍性，一般和較差的預后有關(guān)[21]。GSEA 富集分析顯示，高風險組基因集富集于細胞周期、DNA 復制及癌癥通路，以上結(jié)果指向風險亞組間具有免疫微環(huán)境差異。

腫瘤免疫細胞浸潤會影響ICI 的治療效果。有研究表明，CD8+T 細胞浸潤多的患者對培溴利珠單抗的治療反應優(yōu)于浸潤少的患者[22]。為探討GIlncsig與腫瘤免疫細胞浸潤的關(guān)系，本研究從TCGA 數(shù)據(jù)庫下載HCC 序列的免疫細胞浸潤數(shù)據(jù)，并進行風險亞組間比較，結(jié)果顯示，GIlncsig 和CD8+T 細胞浸潤結(jié)果呈負相關(guān)。而GIlncsig 高風險評分與CD274、CTLA4、HAVCR2、IFNG、PDCD1、TNF 的高表達呈正相關(guān)。該結(jié)果可解釋TMB 在臨床應用中的不確定性。因此，在明確TMB 臨床價值的研究中，需重視探究TMB 與免疫細胞的關(guān)系。

Hong W 等[18]研究報道，基于免疫基因組分析的免疫評分可以預測化療和免疫治療的療效，本研究探究了GIlncsig 是否具有相似的能力。結(jié)果顯示，高風險組中，化療藥物阿霉素、絲裂霉素C 的IC50值更低；同時，本次GSEA 富集分析顯示，高風險組基因集富集于細胞周期、DNA 復制等方面，兩者可相互印證，表明該模型具有預測化療敏感性的潛在能力。

研究顯示[23-25]，部分建模過程中被鑒定的lncRNAs，如MIR210HG、LUCAT1、AC245041.2 等在不同癌癥類型的惡性表型中均起著重要作用，其高表達可能促進了腫瘤的發(fā)生發(fā)展。而本研究提出的AC010643.1 則是首次被發(fā)現(xiàn)。因此，GIlncsig 可以識別新的生物標志物。可以推測，GIlncsig 不僅可以作為肝癌患者的獨立預后因子，用于預測患者預后及腫瘤分級分期，而且可能是評價肝癌患者基因組不穩(wěn)定性、分子特征及免疫特點的指標。

本研究的局限性：lncRNAs 數(shù)據(jù)不完整，缺少不同平臺數(shù)據(jù)庫的驗證，對基因組不穩(wěn)定與免疫治療的關(guān)系研需要進一步深入。GIlncsig 是基于體細胞突變的計算框架確定，需要實驗生物學家進一步進行功能研究，以了解其具體調(diào)控機制。

綜上所述，本研究利用體細胞突變計算框架識別肝細胞癌的基因組不穩(wěn)定性相關(guān)lncRNAs，建立了以基因組不穩(wěn)定lncRNAs 為基礎(chǔ)的風險評分模型，該模型可以較好的預測肝細胞癌患者的預后。