肝癌m6A RNA甲基化調(diào)控因子風(fēng)險預(yù)測模型的構(gòu)建及其臨床意義

王進(jìn)，毛德利，趙鑫，湯權(quán)，董曉強(qiáng)

1.蘇州市獨(dú)墅湖醫(yī)院(蘇州大學(xué)附屬獨(dú)墅湖醫(yī)院)普外科，江蘇 蘇州 215000；2.蘇州大學(xué)附屬第一醫(yī)院普外科，江蘇 蘇州 215006

肝癌是臨床上常見的惡性腫瘤之一，發(fā)病率居全球第5位，惡性程度高，預(yù)后差，全球每年新發(fā)肝癌病人超過74萬，每年因肝癌而死亡病人約70萬[1]。因此，尋找可以判斷肝癌的獨(dú)立預(yù)后因子，對于肝癌病人的后期治療尋找新的靶點(diǎn)及改善病人總體生存率有著至關(guān)重要的作用。RNA修飾是一種轉(zhuǎn)錄后水平的調(diào)控方式，其中RNA甲基化修飾是主要形式之一。mRNA是遺傳信息到蛋白質(zhì)的中間體，mRNA修飾直接導(dǎo)致蛋白質(zhì)改變從而影響細(xì)胞功能,并且研究顯示mRNA甲基化是可逆、動態(tài)的修飾過程。其中N6-甲基腺嘌呤(N6-methyladenosine，m6A)修飾是mRNA含量最豐富一種甲基化修飾形式[2]。m6A廣泛參與調(diào)控生命周期的各個階段，包括mRNA剪接、加工、翻譯和降解等，在惡性腫瘤的進(jìn)展中發(fā)揮重要作用[3-7]。

近年來，多項研究表明m6A調(diào)控因子突變與癌癥的發(fā)生發(fā)展密切相關(guān)，其中甲基化調(diào)控因子METTL3缺失增加了胰腺癌細(xì)胞對化療和放療敏感性，但是對胰腺癌細(xì)胞增殖無明顯影響[8]，Lin等[9]研究顯示METTL3在結(jié)直腸癌與肺腺癌中高表達(dá)，通過促進(jìn)表皮生長因子受體(EGFR)和河馬信號效應(yīng)器TAZ等mRNA翻譯，并促進(jìn)肺腺癌細(xì)胞生長、侵襲。甲基化修飾中沉默負(fù)責(zé)“擦除”的調(diào)控因子ALKBH5后可明顯抑制惡性膠質(zhì)瘤干細(xì)胞增殖。而ALKBH5通過使FOXM1 mRNA去甲基化來增強(qiáng)FOXM1表達(dá)[10]。Chen等[11]研究發(fā)現(xiàn)，在肝癌中METTL3被激活后甲基化修飾SOCS1 mRNA通過YTHDF2依賴性機(jī)制，并影響肝癌病人預(yù)后[11-12]。目前，尚無m6A調(diào)控因子與肝癌病人臨床病理及預(yù)后關(guān)系的綜合性報道。本研究從癌癥基因組圖譜(TCGA)數(shù)據(jù)庫(n=424)和國際癌癥基因組聯(lián)盟(ICGC)數(shù)據(jù)庫(n=445)中獲得基因表達(dá)譜數(shù)據(jù)及相應(yīng)的臨床病理數(shù)據(jù)進(jìn)行分析，證實了m6A調(diào)控因子在肝癌病人中存在差異性表達(dá)。利用單因素Cox比例風(fēng)險回歸分析方法篩選6個關(guān)鍵的m6A調(diào)控因子并建立風(fēng)險預(yù)測模型，為肝癌病人的治療及改善預(yù)后提供了重要的參考依據(jù)。

資料與方法

一、數(shù)據(jù)資料

50例正常組織和374例肝癌組織的基因表達(dá)譜數(shù)據(jù)和相應(yīng)的臨床資料獲取自TCGA數(shù)據(jù)庫(https://portal.gdc.cancer.gov/)，243例正常組織和202例肝癌病人基因表達(dá)譜數(shù)據(jù)和相應(yīng)的臨床資料獲取自ICGC數(shù)據(jù)庫(https://icgc.org/)。

二、m6A RNA甲基調(diào)控因子的挑選

修飾調(diào)控因子中主要負(fù)責(zé)“書寫”的甲基轉(zhuǎn)移酶是最主要的調(diào)控因子，包括KIAA1429[13]，METTL3[14-15]，METTL14[14-15]，RBM15B[16-17]，RBM15[16]，Wtap[13, 15]，ZC3H13[18]；負(fù)責(zé)“讀取”的調(diào)控因子包括HNRNPC、HNRNPA2B[19]，YTHDF、YTHDF2、YTHDF3[20]，YTHDC1[21]，YTHDC2[22]，LRPPRC[23-24]，F(xiàn)MR1[24]；負(fù)責(zé)“擦除”的調(diào)控因子主要由FTO[25]，ALKBH5[26]組成。

三、m6A RNA甲基調(diào)控因子差異性表達(dá)及相關(guān)作用的分析

利用edgeR函數(shù)包對從數(shù)據(jù)庫獲取的正常組織和肝癌組織基因表達(dá)譜數(shù)據(jù)進(jìn)行歸一化處理，提取18種m6A RNA甲基化調(diào)控因子表達(dá)量，分析m6A RNA甲基化調(diào)控因子在腫瘤組織及正常組織表達(dá)的差異性，其中l(wèi)og2fold change(log2FC)絕對值>1，P<0.05的基因為差異表達(dá)基因。Corrplot R包分析調(diào)節(jié)因子在肝腫瘤組織中表達(dá)相關(guān)性，并計算相關(guān)系數(shù)。使用基因集富集分析(Gene Set Enrichment Analysis，GSEA)來探索調(diào)節(jié)因子在肝癌中的相關(guān)功能。

四、風(fēng)險預(yù)后模型的建立及評估

為了確定m6A RNA甲基化調(diào)控因子在肝癌中的預(yù)后價值，使用單因素Cox回歸分析方法篩選出12個與總體生存率相關(guān)基因(P<0.05)。使用 LASSO(the least absolute shrinkage and selection operator)算法從12個基因中確定了6個基因來建立風(fēng)險預(yù)測模型。每位病人風(fēng)險值等于6個基因中每個基因表達(dá)量乘于系數(shù)之和(Coefi是風(fēng)險系數(shù)，Xi是基因表達(dá)量)。計算公式如下：

根據(jù)風(fēng)險值將TCGA和ICGC數(shù)據(jù)庫中病人分為高風(fēng)險組與低風(fēng)險組，分析不同風(fēng)險組病人臨床病理特征及總體生存率關(guān)系，并繪制Kaplan-Meier受試者工作特征(ROC)及校準(zhǔn)曲線。為了評估風(fēng)險值是否可以作為獨(dú)立預(yù)后指標(biāo)，我們對TCGA數(shù)據(jù)庫進(jìn)行單因素和多因素Cox回歸分析，并在ICGC數(shù)據(jù)庫進(jìn)行驗證。

五、統(tǒng)計學(xué)分析

Wilcoxon秩和檢驗比較正常組織與腫瘤組織基因表達(dá)差異性；皮爾遜(Pearson)檢驗用于計算兩組基因相關(guān)系數(shù)；使用LASSO Cox回歸算法來建立風(fēng)險預(yù)后模型；高風(fēng)險組和低風(fēng)險組的總體生存率由Kaplan-Meier方法完成； 由ROC曲線檢驗風(fēng)險評分預(yù)測的準(zhǔn)確性；使用單變量和多變量Cox回歸分析來評估風(fēng)險值與臨床病理特征在肝癌預(yù)后中的危險程度；Wilcoxon秩和檢驗用于比較不同風(fēng)險組及臨床分期之間風(fēng)險基因表達(dá)差異性。所有統(tǒng)計分析由SPSS(19.0版) and R v3.6.0軟件完成。P<0.05為差異有統(tǒng)計學(xué)意義。

結(jié) 果

一、m6A RNA甲基化調(diào)控因子在肝癌組織及正常組織中的表達(dá)差異性

TCGA數(shù)據(jù)庫的基因表達(dá)譜數(shù)據(jù)分析結(jié)果顯示：與正常組織相比，18個m6A RNA 甲基化調(diào)控因子在腫瘤組織中均高表達(dá)(log2FC>1)，除METTL14和ZC3H13之外，其余調(diào)控因子表達(dá)差異均有統(tǒng)計學(xué)意義。其中以下因子差異表達(dá)倍數(shù)較大(圖1A)：METTL3(log2FC=1.352，P=7.37E-24)；KIAA1429(log2FC=1.233，P=5.59E-24)；RBM15B(log2FC=1.013，P=6.26E-23)。ICGC數(shù)據(jù)庫的基因表達(dá)譜數(shù)據(jù)結(jié)果顯示：與正常組織相比，腫瘤組織ZC3H13調(diào)控因子顯著低表達(dá)(log2FC<-1)，其余調(diào)控因子均高表達(dá)，除METTL14之外其余表達(dá)差異均具有統(tǒng)計學(xué)意義。其中KIAA1429 (log2FC=1.153，P=4.14E-58)以及RBM15B (log2FC=1.01，P=3.96E-48)差異表達(dá)倍數(shù)較大(圖1B)。ICGC 數(shù)據(jù)集基因表達(dá)量結(jié)果與TCGA數(shù)據(jù)集基本一致。

圖1 癌癥基因組圖譜(TCGA)數(shù)據(jù)庫(A)和國際癌癥基因組聯(lián)盟(ICGC)數(shù)據(jù)庫(B)中m6A RNA甲基化調(diào)控因子在肝癌組織及正常組織的表達(dá)差異性

二、m6A RNA甲基化調(diào)控因子相互作用分析

通過Corrplot R包對TCGA數(shù)據(jù)庫中m6A RNA甲基化調(diào)控因子進(jìn)行相關(guān)性分析，結(jié)果顯示調(diào)節(jié)因子之間有著復(fù)雜的聯(lián)系，其中負(fù)責(zé)“擦除”的RBM15、METTL14、METTL3、KIAA1429和Wtap與其余調(diào)節(jié)因子之間彼此顯著相關(guān)(圖2A)。最后GSEA顯示HNRNPC、 HNRNPA2B1和METTL3與腫瘤發(fā)生及進(jìn)展有關(guān)，比如癌癥[標(biāo)準(zhǔn)化富集得分(NES)=1.908, 標(biāo)準(zhǔn)化P<0.001]，ERBB信號通路(NES=1.925，標(biāo)準(zhǔn)化P<0.001)，RIG-I-like-受體信號通路(NES=2.073，標(biāo)準(zhǔn)化P<0.001)，Wnt信號通路(NES=2.070，標(biāo)準(zhǔn)化P<0.001)(圖2B)。以上分析結(jié)果提示m6A RNA甲基化調(diào)控因子在肝癌進(jìn)程中和生存預(yù)后發(fā)揮著重要作用。

圖2 癌癥基因組圖譜(TCGA)數(shù)據(jù)庫中m6A RNA甲基化調(diào)控因子相互作用分析 A.Corrplot R包分析m6A RNA甲基化調(diào)控因子的相關(guān)性；B.腫瘤發(fā)生和進(jìn)展相關(guān)調(diào)控因子的基因集富集分析(GSEA)

三、建立預(yù)后風(fēng)險模型

利用單因素Cox回歸分析m6A甲基化調(diào)控因子對病人生存預(yù)后的影響。結(jié)果顯示YTHDF2、LRPPRC、YTHDF1、METTL3、RBM15B、KIA-A1429、HNRNPA2B1、HNRNPC、WTAP、YTHDC1、RBM15、ZC3H13顯著影響病人的總體生存率(圖3A)，差異均有統(tǒng)計學(xué)意義(P<0.05)。其中YTHDF2、LRPPRC、YTHDF1、METTL3、RBM15-B、KIAA1429、HNRNPA2B1、HNRNPC、WTAP、YTHDC1、RBM15這些調(diào)控因子的風(fēng)險比(HR)>1，屬于高風(fēng)險基因；而ZC3H13調(diào)控因子的HR小于1，屬于低風(fēng)險基因。使用LASSO算法，從以上12個基因中篩選出6個基因(KIAA1429、ZC3H13、METTL3、YTHDF1、YTHDF2、LRPPRC)，建立風(fēng)險預(yù)測模型并計算風(fēng)險值。根據(jù)中位風(fēng)險值將肝癌病人分為高風(fēng)險組與低風(fēng)險組(圖3B～D)。GSEA顯示風(fēng)險值與WNT信號通路及泛素介導(dǎo)的蛋白質(zhì)水解密切相關(guān)(圖3E、F)。

圖3 在癌癥基因組圖譜(TCGA)數(shù)據(jù)庫中篩選6個肝癌預(yù)后相關(guān)基因建立風(fēng)險預(yù)測模型 A.單因素Cox回歸分析m6A甲基化調(diào)控因子對病人的生存預(yù)后的影響；B、C.LASSO Cox回歸分析篩選6個預(yù)后相關(guān)基因建立預(yù)后模型；D.被用來建立預(yù)后模型的6個基因系數(shù)；E、F.基因集富集分析(GSEA)顯示風(fēng)險值與Wnt信號通路及泛素介導(dǎo)的蛋白質(zhì)水解密切相關(guān)

四、高、低風(fēng)險組肝癌病人預(yù)后生存分析以及ROC曲線的繪制

使用Kaplan-Meier曲線對高風(fēng)險組與低風(fēng)險組預(yù)后生存分析顯示，在TCGA與ICGC數(shù)據(jù)庫中高風(fēng)險組與低風(fēng)險組總體生存時間差異有統(tǒng)計學(xué)意義，高風(fēng)險組病人的5年生存率明顯低于低風(fēng)險組(P<0.001)(圖4A)。同時使用ROC曲線檢測風(fēng)險值預(yù)測肝癌病人5年生存率的準(zhǔn)確性，結(jié)果顯示：曲線下面積(AUC)分別為61.9%(TCGA)與76.5%(ICGC)，提示風(fēng)險預(yù)測模型可以有效預(yù)測肝癌病人5年生存率。兩組數(shù)據(jù)AUC值相差較大可能是不同區(qū)域來源的病人存在的較大異質(zhì)性導(dǎo)致。TCGA肝癌數(shù)據(jù)主要來源于美國，ICGC肝癌數(shù)據(jù)來自于日本，統(tǒng)計發(fā)現(xiàn)TCGA數(shù)據(jù)肝癌病人5年死亡比例占39%，ICGC數(shù)據(jù)肝癌病人5年死亡比例占18%(圖4B)。校準(zhǔn)曲線顯示風(fēng)險模型對肝癌病人5年總體生存率預(yù)測效果較好(圖 4C)。

五、單因素與多因素Cox回歸分析

對TCGA和ICGC數(shù)據(jù)庫中肝癌病人臨床病理特征及風(fēng)險值進(jìn)行單因素與多因素Cox回歸分析。在TCGA數(shù)據(jù)庫中，單因素分析顯示分期(P<0.001)、T分期(P<0.001)、N分期(P<0.05)、M分期(P<0.05)、風(fēng)險值(P<0.001)與病人總體生存率顯著相關(guān)。將分期、T分期、N分期、M分期、風(fēng)險值納入多因素Cox回歸分析，結(jié)果顯示風(fēng)險值(P<0.001)與病人總體生存率顯著相關(guān)(P<0.001)(圖5A)。在ICGC數(shù)據(jù)庫中，單因素分析顯示分期(P<0.001)、性別(P<0.05)、風(fēng)險值(P<0.001)與病人總體生存率顯著相關(guān)。將分期、性別、風(fēng)險值納入多因素Cox回歸分析，顯示分期(P<0.001)、性別(P=0.006)，風(fēng)險值(P=0.001)與病人總體生存率顯著相關(guān)。以上結(jié)果提示風(fēng)險值可以作為肝癌的獨(dú)立預(yù)后因子(圖5B)。

注：AUC.曲線下面積；TCGA.癌癥基因組圖譜；ICGC.國際癌癥基因組聯(lián)盟。圖4 風(fēng)險值在肝癌在診斷及預(yù)后中的價值 A1、A2.不同風(fēng)險組肝癌病人總體生存率；B1、B2.受試者工作特征(ROC)曲線判斷風(fēng)險預(yù)測模型對肝癌病人5年生存率預(yù)測的準(zhǔn)確性；C1、C2.校準(zhǔn)曲線檢驗風(fēng)險預(yù)測模型對肝癌病人5年生存率預(yù)測的準(zhǔn)確性

圖5 在癌癥基因組圖譜(TCGA)數(shù)據(jù)庫和國際癌癥基因組聯(lián)盟(ICGC)數(shù)據(jù)庫中不同臨床病理特征及風(fēng)險值的單因素及多因素Cox回歸分析 A1.TCGA數(shù)據(jù)庫中單因素分析；A2.TCGA數(shù)據(jù)庫中多因素分析；B1.ICGC數(shù)據(jù)庫中單因素分析；B2.ICGC數(shù)據(jù)庫中多因素分析

六、不同風(fēng)險組肝癌病人臨床病理特征及預(yù)后的差異性

分析在不同風(fēng)險組中，臨床病理特征、生存預(yù)后以及6個調(diào)節(jié)因子表達(dá)的差異性。結(jié)果顯示YTHDF2、LRPPRC、YTHDF1、METTL3、KIAA1429在高風(fēng)險組高表達(dá)，而ZC3H13在高風(fēng)險組低表達(dá)。進(jìn)一步分析顯示：在TCGA數(shù)據(jù)庫中，不同風(fēng)險組之間分期(P<0.05)、生存時間(P<0.001)存在著顯著性差異(圖6A)。在ICGC數(shù)據(jù)庫中，不同風(fēng)險組之間分期(P<0.01)、分級(P<0.01)、年齡(P<0.05)、生存狀態(tài)(P<0.05)存在著顯著性差異(圖6B)。

圖6 在癌癥基因組圖譜(TCGA)數(shù)據(jù)庫(A)和國際癌癥基因組聯(lián)盟(ICGC)數(shù)據(jù)庫(B)中，不同風(fēng)險組臨床病理特征及預(yù)后的差異性

討 論

肝癌是全球第五大癌癥，具有高死亡率及侵襲性特點(diǎn)，占全球癌癥死亡率的9.1%，肝癌的發(fā)病率和死亡率在逐年增加[27]。肝癌病人往往確診時已晚期，同時術(shù)后易復(fù)發(fā)、轉(zhuǎn)移、耐藥，故死亡率較高。因此深入了解肝癌的發(fā)病機(jī)制，尋找新的診斷及預(yù)后標(biāo)志物對臨床診斷和治療工作有重要的指導(dǎo)意義。本次研究通過對肝癌中m6A異常甲基修飾調(diào)控因子進(jìn)行分析，證實了其在肝癌發(fā)生發(fā)展中的重要性。從中篩選出具有影響肝癌病人預(yù)后的6個調(diào)控因子，并建立風(fēng)險預(yù)測模型，其在臨床診療中具有很大程度的參考意義。

既往已有大量研究報道m(xù)6A RNA甲基化調(diào)控因子在多種癌癥中發(fā)揮重要作用[8-9, 25]。在乳腺癌中，METTL3可通過m6A修飾促進(jìn)HBXIP的表達(dá)，從而促進(jìn)腫瘤增殖[28]。對METTL3進(jìn)行的m6A甲基化修飾是急性髓細(xì)胞白血病發(fā)生所必須的方式[29]。METTL3、YTHDF1和YTHDF2甲基化調(diào)控因子已被證明影響肝癌病人預(yù)后，并且可能作為肝癌新的診斷治療生物標(biāo)志[11-12]。Cheng等[30]發(fā)現(xiàn)，在肝癌中負(fù)責(zé)“書寫”的KIAA1429調(diào)控因子通過m6A修飾抑制癌基因ID2表達(dá)，從而導(dǎo)致肝癌侵襲力增加。在本次研究中，通過總結(jié)目前常見的18種m6A RNA甲基化調(diào)控因子，從TCGA和ICGC數(shù)據(jù)庫中獲取基因表達(dá)圖譜及臨床信息，并分析其在腫瘤組織及正常組織表達(dá)的差異性，發(fā)現(xiàn)其中負(fù)責(zé)“讀寫”的METTL3、KIAA1429和RBM15B調(diào)控因子表達(dá)顯著高于其他調(diào)控因子。

既往研究顯示，乙型肝炎病毒感染是肝癌的主要致病因素之一，乙型肝炎病毒感染通過下調(diào)肝癌病人中RIG-I-like受體信號從而抑制Ⅰ型干擾素的表達(dá)[31]。多種惡性腫瘤與ERBB信號失調(diào)相關(guān)，比如肺癌、卵巢癌、乳腺癌，ERBB信號失調(diào)往往導(dǎo)致更差的生存預(yù)后[32-33]。過去的研究認(rèn)為ERBB信號是肝臟急性損傷時的保護(hù)機(jī)制，隨著研究的深入，發(fā)現(xiàn)肝臟在慢性條件刺激下，ERBB信號參與肝癌的形成[34]。Wnt信號主要參與結(jié)直腸癌生物學(xué)過程，在肝癌中，β-連環(huán)蛋白N末端區(qū)域突變，導(dǎo)致β-連環(huán)蛋白穩(wěn)定性增加，在肝臟組織中累積，進(jìn)而促進(jìn)腫瘤進(jìn)展。在本次研究中，GESA結(jié)果顯示m6A RNA甲基化調(diào)控因子與上述介導(dǎo)肝癌發(fā)生發(fā)展的多種信號傳導(dǎo)通路具有很大程度上的相關(guān)性。以上的實驗結(jié)果提示m6A RNA異常甲基化調(diào)控因子可能參與肝癌發(fā)生、發(fā)展。在肝癌診斷和預(yù)后中起到至關(guān)重要的作用。

在TCGA數(shù)據(jù)庫中，利用Cox回歸分析篩選12個預(yù)后相關(guān)基因，進(jìn)一步的LASSO 算法確定了6個影響肝癌病人生存預(yù)后的關(guān)鍵調(diào)節(jié)因子，包括KIAA1429、ZC3H13、METTL3、YTHDF1、YTHDF2、LRPPRC。由此建立風(fēng)險預(yù)測模型，根據(jù)中位風(fēng)險值將肝癌病人分為低風(fēng)險組與高風(fēng)險組。分析發(fā)現(xiàn)高風(fēng)險組與低風(fēng)險組之間病人總體生存率差異有統(tǒng)計學(xué)意義，單因素與多因素Cox回歸分析結(jié)果均顯示風(fēng)險值可以作為獨(dú)立預(yù)后的指標(biāo)。然而，與既往的報道[34]不同的是，本研究的分析結(jié)果顯示，在TCGA數(shù)據(jù)庫中N分期未表現(xiàn)出獨(dú)立的危險因素，這可能與本研究納入的標(biāo)本量較小有關(guān)。為了驗證所建立模型的有效性，本研究使用ROC曲線和ICGC測試集數(shù)據(jù)進(jìn)行了檢驗，最終確定了6個m6A甲基化調(diào)控因子構(gòu)建的風(fēng)險預(yù)測模型對判斷肝癌預(yù)后的可靠性。

綜上所述，本研究通過生物信息學(xué)分析證實m6A RNA甲基化調(diào)控因子在肝癌腫瘤組織中與正常組織中表達(dá)具有顯著性差異，m6A RNA甲基化調(diào)控因子的異常表達(dá)與肝癌臨床病理及病人預(yù)后密切相關(guān)，基于篩選的6個關(guān)鍵m6A RNA甲基化調(diào)控因子建立的風(fēng)險預(yù)測模型對指導(dǎo)肝癌后期個性化治療及改善預(yù)后具有重要的參考價值。