基于生物信息學的肝細胞癌預后模型的構建

李青玲，李樹德，衛(wèi)小娟,2，黃映光，李思熳*

(1.昆明醫(yī)科大學 基礎醫(yī)學院 生物化學與分子生物學系，云南 昆明 650500；2.中國科學院 昆明動物研究所，云南 昆明 650201;3.云南省第一人民醫(yī)院 普外一科，云南 昆明 650034)

肝癌發(fā)病率高居世界第六，是癌癥相關死亡第二大原因.肝細胞癌(Hepatocellular carcinoma，HCC)占原發(fā)性肝癌的75%～85%[1].由于HCC病情進展迅速，并發(fā)癥多，大多數(shù)患者確診時已達中晚期[2]，且大量早期HCC患者在手術治療后還會復發(fā).目前對HCC的治療決策多依賴于病理及影像學證據(jù)，缺乏分子學指標.而基于HCC特異性分子構建的預后預測模型則可進一步輔助臨床醫(yī)師進行HCC診治[3].腫瘤突變負荷(tumor mutation burden，TMB)作為衡量免疫治療效果的可靠指標，在HCC預后評估中至關重要.TMB表示外顯子編碼區(qū)每兆堿基中替換、基因插入或突變的總數(shù)[4]，是預測免疫檢查點抑制劑治療多種腫瘤(如肺癌、子宮內(nèi)膜癌、乳腺癌和結腸直腸癌)療效的有效指標[5-9].多項研究[5-9]中已初步揭示TMB與免疫微環(huán)境之間的聯(lián)系及TMB與HCC中的免疫治療之間的聯(lián)系，但由于其測序的復雜性和高成本，不太可能廣泛應用于免疫治療療效的臨床評價.Zhang等[10]基于2020年公共數(shù)據(jù)進行了HCC中TMB相關分析，然而，近2年the Cancer Genome Atlas(TCGA)數(shù)據(jù)庫更新后鮮有相關研究報道.

本研究旨在基于(TCGA)數(shù)據(jù)庫探索HCC中腫瘤突變概況，然后通過基因差異分析及預后分析篩選出TMB相關的預后基因，并基于篩選的預后基因構建HCC預后預測模型.

1 材料與方法

1.1 數(shù)據(jù)來源

從TCGA數(shù)據(jù)庫(http：www.ncbi.nlm.nih. gov/geo)下載與HCC相關的基因組和臨床病理數(shù)據(jù)，基因組數(shù)據(jù)格式為FPKM.由于數(shù)據(jù)丟失，部分樣本被刪除，最終選擇357個腫瘤樣本進行后續(xù)分析.同樣，從TCGA數(shù)據(jù)庫中獲取HCC體細胞突變數(shù)據(jù)，并使用R軟件中的“maftool”包進行分析[11].

1.2 TMB計算和分組

運用Perl腳本計算TMB.并通過X-tile結合患者總體生存率(overall survival，OS)找到TMB的最佳截斷值，然后依據(jù)該值將患者分為TMB-H組和TMB-L組.

1.3 獲取差異基因

使用R包“l(fā)imma”處理數(shù)據(jù)，篩選出與TMB相關的差異表達基因(differentially expressed genes，DEGs).采用|log 2 FC|>1,false discovery rate(FDR)<0.05篩選出所有差異基因，并使用“pheatmap”包對差異基因進行分層聚類.

運用R包“org.Hs.eg.db”進行Gene Ontology(GO)富集分析和Kyoto Encyclopedia of Genes and Genomes(KEGG)通路分析，“clusterProfiler”和“ggplot2”包用于注釋和繪圖[12].采用P<0.05篩選出潛在通路.

1.4 TMB預后指數(shù)的構建

使用Gene Expression Profiling Interactive Analysis(GEPIA)數(shù)據(jù)庫(http://gepia.cancer-pku.cn/)[13]分析差異基因的OS率.然后根據(jù)logrankP<0.05篩選出關鍵基因.然后根據(jù)風險評分(risk score,RS)將所有患者分為高危組和低危組.采用K-M曲線表示RS水平與OS之間的關系.最后，繪制受試者生存工作特征曲線(receiver operating characteristic curve,ROC)，計算曲線下面積(caculate the area under the curve,AUC)，以評估預后模型的價值.

1.5 統(tǒng)計分析

使用R Studio(https://www.rstudio.com)進行統(tǒng)計分析[14].OS采用K-M和log-rank檢驗方法計算.2組連續(xù)變量比較采用t檢驗或單因素方差分析.P值<0.05，表示差異有統(tǒng)計學意義.

2 結果

2.1 TMB在HCC中的研究概況

首先，通過評估每個HCC樣本的突變情況發(fā)現(xiàn)：單核苷酸多態(tài)性(Single nucleotide polymorphism，SNP)、錯義突變和C>T突變是HCC中最常見的變異類型.其中，C>T變異最高，達 15 420.HCC中突變前10位的基因為：TP53、TTN、CTNNB1、MUC16、PCLO和ALB、RYR2、ABCA13、MUC4頻率均高于10(P<0.001)(圖1).

2.2 突變基因的瀑布效應和相互關系

圖2(a)描繪了HCC中前30個突變基因的瀑布圖.我們進一步分析前20個突變基因表達相關度，結果發(fā)現(xiàn):CTNNB1與TP53、AX1N1表達負相關(P<0.05)，CTNNB1與SPTA1、OBSCN；MUC16與SPTA1、ABCA13、RYR2；PCLO與CACNA1E、HMCN1；RYR2與XIRP2；APOB與LRP1B；FLG與CSMD3、OBSCN；LRP1B與CACNA1E;HMCN1與SPTA1；TTN與LRP1B表達正相關(P<0.05)(圖2(b)).

2.3 TMB與臨床相關性

利用X-tile確定TMB=4.61為最佳截止值，并據(jù)此將357個樣本分為TMB-H組和TMB-L組.TMB-H組和TMB-L組的平均值分別為8.19和2.79(圖3(a)).生存分析結果顯示TMB與HCC預后顯著負相關(P<0.05).TMB高的患者預后較差[15](P<0.05)(圖3(b)).圖3(c)～(e)反映了TMB與臨床病理特征的關系：在HCC患者中TMB與N分期(P<0.05)、年齡(P<0.05)、性別(P<0.05)相關.總的來說，老年(>65歲)無淋巴轉移的男性TMB較高，但TMB水平與T/M分期、臨床分期和分級無明顯相關性(圖3(f)～(i)).

(a)變異的分類 (b)變異類型 (c)單核苷酸多態(tài)性分類 (d)突變前10位的基因

(a)HCC中前30個突變基因的瀑布圖

(b)HCC中前20個突變基因的相互關系

(a) TMB在各組間的分布 (b)TMB與HCC預后相關性 (c)TMB與年齡的關系

(d)TMB與性別的關系 (e)TMB與N期的關系 (f)TMB與T期的關系

(g)TMB與M期的關系 (h) TMB與臨床分期的關系 (i) TMB與臨床分級的關系

2.4 差異基因的富集分析

首先通過基因差異分析，篩選出39個與TMB相關的DEGs.圖4(a)為DEGs的聚類熱圖.我們進一步通過GO功能富集分析研究DEGs相關的主要生物過程(biological process，BP)、分子功能(molecular function，MF)和細胞成分(cellular component，CC) (圖4(b)).KEGG通路分析顯示，DEGs主要富集于細胞外基質(zhì)組織、細胞外結構組織、硫化合物生物合成過程和硫化合物代謝過程(圖4(c)).Gene set enrichment analysis(GSEA)分析結果顯示，TMB-H組富集在RNA代謝過程，TMB-L組富集在其他途徑包括調(diào)控超分子纖維組織、調(diào)控細胞器組織、正調(diào)控細胞器組織、肌動蛋白絲組織、調(diào)控細胞骨架組織、肌動蛋白聚合或解聚、調(diào)控肌動蛋白絲基過程、調(diào)控肌動蛋白絲組織(圖4(d)).

2.5 COX回歸分析

基于生存分析，可進一步發(fā)現(xiàn)：SFRP4、IL7R、FBLN2、COLEC10和CHGA的表達可能是影響HCC患者預后的獨立因素(圖5).由此，可利用COX回歸分析建立預后預測模型，并利用ROC曲線驗證模型的準確性.通過R包“merge”下載了357例HCC病例的轉錄組數(shù)據(jù).根據(jù)多元COX回歸模型，計算TMBPI為:

PI=(0.088 945 38×SFRP4-0.009 912 19×COLEC10-0.000 690 51×CHGA-0.024 612 34×FBLN2-0.165 479 22×IL7R).

每個個體的RS通過TMBPI得到，然后可根據(jù)風險評分的中位數(shù)將樣本分為高風險組和低風險組(圖6(a)).該模型1,3,5年生存率的AUC分別為0.64、0.67、0.59(圖6(b)).結果顯示高危組患者預后較差(P<0.001)(圖6(c)).

(a)DEGs的聚類熱圖

(b)DEGs的GO功能富集分析 (c)DEGs的KEGG通路分析

(d)TMB組間的GSEA分析

(a)SFRP4的表達與HCC患者預后相關 (b) IL7R的表達與HCC患者預后相關 (c) FBLN2的表達與HCC患者預后相關

(d)CHGA的表達與HCC患者預后相關 (e) COLEC10的表達與HCC患者預后相關

3 討論與結論

HCC侵襲性強，且具有高度異質(zhì)性.多項研究[16-18]表明肝癌細胞內(nèi)信號通路與正常肝細胞迥異.例如，PI3K/AKT和IKK/NF-κB通路在HCC中被激活，從而促進細胞增殖并誘導上皮-間質(zhì)轉化[19,20].此外，免疫系統(tǒng)在腫瘤的發(fā)展過程中起著至關重要的作用.Marina等[21]觀察到MYC/CTNNB1通路的激活促進了免疫逃避和肝癌藥物耐藥性.TMB作為免疫系統(tǒng)識別和攻擊肝癌細胞的關鍵因素，在很大程度上影響患者的預后.此外，在各種癌癥中可觀察到TMB與客觀緩解率顯著相關[22].在本研究中，我們發(fā)現(xiàn)錯義突變是HCC中最常見的類型.研究表明TP53的錯義突變可將HCC分化為不同的亞型，并可能促進疾病的進展[23].在HCC中SNP變異占絕大多數(shù)，在我們的研究中，C>T變異數(shù)量達到 15 420，這與HCC的發(fā)生密切相關[24,25].

功能分析結果表明，DEGs主要與細胞外基質(zhì)組織和硫化合物代謝過程有關.細胞外基質(zhì)是HCC免疫微環(huán)境的組成部分，靶向給藥已取得一定成效[26].硫化合物代謝可能與HCC中的自噬密切相關[27].此外，我們發(fā)現(xiàn)RNA代謝過程的活性與TMB之間存在明顯的負相關.越來越多的證據(jù)表明，不同的非編碼RNA(non-coding RNAs，ncRNAs)參與多種肝臟疾病(包括乙肝、丙肝和非酒精性脂肪肝病)代謝狀態(tài)調(diào)控，并最終導致HCC[28-30].這可能是TMB-H組總體預后不佳的最主要原因之一.

SFRP4、IL7R、FBLN2、COLEC10和CHGA在HCC發(fā)生發(fā)展中發(fā)揮重要作用.SFRP4位于染色體7p14.1上，含有一個富含半胱氨酸的結構域.該結構域可通過與Wnt直接結合來調(diào)節(jié)Wnt信號通路，從而形成沉默復合體并抑制HCC[31].IL-7R在T細胞分化和淋巴細胞發(fā)育中能發(fā)揮作用[32].IL-7R上調(diào)可激活細胞內(nèi)通路，誘導相關分子表達，促進肝癌細胞增殖和遷移[33].COLEC10編碼肝臟膠原凝集素1，是C凝集素家族的一員[34].研究表明[35]，COLEC10受miR-452-5p的調(diào)控，促進肝癌細胞的增殖、侵襲和遷移.CHGA(或CGA)在HCC樣本中高度表達，既往研究[36]表明它可以作為HCC的輔助診斷分子.本研究顯示FBLN2高表達組與低表達組的預后差異有統(tǒng)計學意義.然而，F(xiàn)BLN2在肝癌中的作用及其應用價值尚未見報道，這可能是一個值得探索的靶點.

基于此，本文構建了多元COX回歸模型，實現(xiàn)了對患者個體生存概率的預測和分析.在本研究中，該模型1、3、5年生存率的AUC分別為0.64、0.67、0.59，具有一定的準確度.該模型將為探索HCC新的預后因素提供了新的思路.此外，該模型無須識別患者的體細胞突變，而是通過直接檢測TMB相關的5個基因，這使得基于關鍵基因的靶向測序技術更加常規(guī).基于5個基因的預后預測模型可有效評估患者預后，這可能有助于指導肝癌患者的臨床治療.在常規(guī)臨床實踐中，臨床醫(yī)師多通過病理分期決定患者治療手段.但HCC病理活檢損傷大，多數(shù)患者依賴術前影像學(如CT、MRI)或常規(guī)生化指標(如AFP,CEA,CA199)等評估病情.本研究創(chuàng)新性地將TMB、DEGs與HCC臨床數(shù)據(jù)相結合，通過COX回歸模型預測患者預后，其可操作性強，侵入性更小.在精準醫(yī)療時代，利用分子分析結合病理分期的預后預測將比傳統(tǒng)方法更準確.

該模型具有對HCC患者的實質(zhì)性預測能力和價值，然而其局限性也應當被意識到.首先，目前研究缺乏患者的免疫治療信息，這限制了我們對TMB介導的HCC免疫微環(huán)境變化的研究.其次，由于TCGA項目中只有一部分樣本具有臨床病理和TMB數(shù)據(jù)可用于分析，目前的樣本量存在限制.最后，由于缺乏肝癌患者的全外顯子測序信息，我們沒有通過本中心臨床數(shù)據(jù)來驗證結論.因此，預測模型的有效性值得在未來使用獨立隊列進行更多的外部驗證.

綜上，本研究建立了一個新的基于5個基因的HCC患者預后預測模型.該模型可能是HCC預后風險分層的有效工具，并為肝癌的臨床診斷結合病理和分子分析提供了新的思路.