基于TCGA數(shù)據(jù)篩選潛在的宮頸癌miRNA預(yù)后標(biāo)志物

苑飛艷，張新民，陳明珠，楊鳴，張立會(huì)*

宮頸癌是嚴(yán)重威脅全球女性健康和生命的惡性腫瘤，是女性癌癥死亡的第四大原因，據(jù)估計(jì)，2020年全球新增病例為60.4萬(wàn)例，死亡人數(shù)為34.2萬(wàn)人[1]。我國(guó)宮頸癌發(fā)病數(shù)約占全球總發(fā)病數(shù)的30%[2]，近年發(fā)展中國(guó)家宮頸癌的發(fā)生率和死亡率仍居高不下，且逐漸趨于年輕化。

宮頸癌主要采用手術(shù)治療，但對(duì)于晚期宮頸癌治療需要輔以放化療[3]，目前癌癥發(fā)生發(fā)展的分子機(jī)制尚不明確，有研究證實(shí)miRNA參與癌癥的發(fā)生、進(jìn)展、轉(zhuǎn)移和耐藥等各個(gè)階段，將miRNA用作臨床生物標(biāo)記物用于診斷、預(yù)測(cè)和預(yù)后是可行的[4]。有研究表明，致癌miRNA上調(diào)或抑癌miRNA下調(diào)可促進(jìn)癌癥的發(fā)生，故在癌癥的靶向治療上可有兩種策略：即用miRNA模擬物用于替代、恢復(fù)下調(diào)的抑癌miRNA；或用拮抗劑減少、抑制上調(diào)的致癌miRNA[5-6]。目前已經(jīng)發(fā)現(xiàn)miRNA表達(dá)譜能夠根據(jù)組織來(lái)源對(duì)腫瘤進(jìn)行分類，且其表達(dá)變異與轉(zhuǎn)移過(guò)程有關(guān)，故miRNA可作為生物標(biāo)志物用于腫瘤的診斷及預(yù)測(cè)[7]。TCGA數(shù)據(jù)收集并整理了多種腫瘤的多組學(xué)及臨床數(shù)據(jù)，因此本文從TCGA數(shù)據(jù)庫(kù)中獲取宮頸癌樣本的miRNA表達(dá)，篩選出宮頸癌預(yù)后相關(guān)的生物標(biāo)記物，為臨床工作及研究提供理論依據(jù)。

1 材料和方法

1.1 數(shù)據(jù)的獲取及處理

從TCGA數(shù)據(jù)下載宮頸癌患者的miRNA表達(dá)數(shù)據(jù)、臨床數(shù)據(jù)，去除重復(fù)數(shù)據(jù)、缺失數(shù)據(jù)，選擇出腫瘤樣本300個(gè)，提取樣本中miRNA表達(dá)數(shù)據(jù)，將相同的成熟miRNA進(jìn)行合并，合并過(guò)程區(qū)分3p和5p異構(gòu)體，共有2 188種miRNA，將NA值轉(zhuǎn)換為0，計(jì)算miRNA在所有樣品種含量為0的樣品數(shù)，如果超過(guò)25%樣品為零，刪除該miRNA，miRNA數(shù)量變?yōu)?88個(gè)，使用miRNAmeConverter數(shù)據(jù)包將miRNA名稱由TCGA數(shù)據(jù)轉(zhuǎn)變成為通用名稱，最終成功轉(zhuǎn)化的585種miRNA。

1.2 數(shù)據(jù)分組

將篩選的300個(gè)腫瘤樣本使用R語(yǔ)言隨機(jī)分組函數(shù)(sample函數(shù))分為訓(xùn)練集和驗(yàn)證集，各150例，訓(xùn)練集用于構(gòu)建預(yù)測(cè)模型，驗(yàn)證集用于驗(yàn)證預(yù)測(cè)模型效能。

1.3 差異miRNA的篩選

對(duì)訓(xùn)練集使用R語(yǔ)言包Survival進(jìn)行單因素Cox生存回歸分析，計(jì)算每個(gè)miRNA與宮頸癌患者生存的風(fēng)險(xiǎn)比(hazard ratio,HR)和P值，從基因表達(dá)量角度分析，HR小于1表示對(duì)生存有利[8]。選出P值小于0.05的與宮頸癌預(yù)后相關(guān)的miRNA。

1.4 風(fēng)險(xiǎn)模型的構(gòu)建

對(duì)訓(xùn)練集將篩選得到的miRNA使用glmnet包進(jìn)行Lasso回歸分析，刪除回歸系數(shù)為零的miRNA，篩選出與生存相關(guān)的miRNA，將每個(gè)miRNA 的回歸系數(shù)作為風(fēng)險(xiǎn)系數(shù)(Coef)，建立風(fēng)險(xiǎn)值公式：Risk score=Coef1×miRNA1Exp+Coef2×miRNA2Exp+……+Coefn×miRNAnExp[9-10]，式中Coef為相應(yīng)miRNA的多因素回歸系數(shù)，Exp為相應(yīng)miRNA的表達(dá)量，根據(jù)Risk score數(shù)值的中位數(shù)，將CC患者分為高風(fēng)險(xiǎn)評(píng)分組和低風(fēng)險(xiǎn)評(píng)分組。使用R語(yǔ)言包Survival，Dplyr，Survival ROC，ggplot2，survivalsvm進(jìn)行數(shù)據(jù)整理,Kaplan-Meier曲線繪制，以及l(fā)ogRanK檢驗(yàn)。使用多因素Cox分析風(fēng)險(xiǎn)值在訓(xùn)練組和驗(yàn)證組中是否屬于獨(dú)立預(yù)后因素。最后繪制風(fēng)險(xiǎn)評(píng)分和基于miRNA臨床預(yù)后模型的 ROC曲線并計(jì)算AUC值以評(píng)估模型的預(yù)測(cè)能力。

2 結(jié)果

2.1 樣本統(tǒng)計(jì)

將TCGA的miRNA表達(dá)數(shù)據(jù)和臨床數(shù)據(jù)中的生存時(shí)間和生存狀態(tài)、種族、分期，是否經(jīng)過(guò)治療進(jìn)行合并后，數(shù)據(jù)特征如表1：

表1 篩選的300個(gè)樣本數(shù)據(jù)特征

2.2 樣本分組

兩組數(shù)據(jù)的種族、生存狀態(tài)、分期以及治療方法，采用χ2檢驗(yàn)，兩個(gè)分組中所有變量的P值都大于0.05，兩組數(shù)據(jù)中每個(gè)變量的組成,差異無(wú)統(tǒng)計(jì)學(xué)意義(見(jiàn)表2)，可以用于建模和驗(yàn)證分析。

表2 比較訓(xùn)練集和驗(yàn)證集差異的統(tǒng)計(jì)學(xué)分析

2.3 差異miRNA的篩選

使用Cox單因素方差分析對(duì)每種miRNA表達(dá)量與生存分析，發(fā)現(xiàn)70種miRNA的P值小于0.05，由于前30個(gè)P<0.01，因此選擇這30個(gè)miRNA作為初步候選miRNA，為了預(yù)后檢測(cè)更加高效易行，使用glmnet包初篩的miRNA進(jìn)行Lasso回歸分析，將回歸系數(shù)為零的miRNA逐漸剔除，最終篩選出的5種miRNA及風(fēng)險(xiǎn)系數(shù)(Coef)見(jiàn)表3。

表3 篩選的5種miRNA及對(duì)應(yīng)的風(fēng)險(xiǎn)系數(shù)

2.4 風(fēng)險(xiǎn)評(píng)分方程的建立及風(fēng)險(xiǎn)模型的建立

對(duì)上述5個(gè)miRNA進(jìn)行多因素Cox回歸分析，根據(jù)相應(yīng)的回歸系數(shù)，建立風(fēng)險(xiǎn)評(píng)分方程：Risk score=Coef1×hsa-miR-101-3pExp+Coef2×hsa-miR-532-5pExp+Coef3×hsa-miR-150-3pExp+Coef4×hsa-miR-378dExp+Coef5×hsa-miR-378a-5pExp，將每個(gè)病人的這5種miRNA風(fēng)險(xiǎn)值相加按照大小排序，以中位數(shù)為界限，將訓(xùn)練集和驗(yàn)證集分為高風(fēng)險(xiǎn)組和低風(fēng)險(xiǎn)組。采用Kaplan-Meier生存分析，繪制KM曲線對(duì)預(yù)后模型進(jìn)行驗(yàn)證，結(jié)果表明：在訓(xùn)練集和驗(yàn)證集中，高風(fēng)險(xiǎn)組的預(yù)后比低風(fēng)險(xiǎn)組差(生存時(shí)間單位為：天)，P值<0.05(見(jiàn)圖1)。

圖1 訓(xùn)練集(a)和驗(yàn)證集(b)中高低風(fēng)險(xiǎn)評(píng)分組患者的生存曲線

2.5 多因素Cox分析

使用多因素Cox分析種族、發(fā)病年齡、腫瘤分期、風(fēng)險(xiǎn)值這四個(gè)因素中，風(fēng)險(xiǎn)值在訓(xùn)練組和驗(yàn)證組中屬于獨(dú)立預(yù)后因素：訓(xùn)練集風(fēng)險(xiǎn)評(píng)分(HR=6.123，95%CI：2.714-13.814，P<0.001)，驗(yàn)證集風(fēng)險(xiǎn)評(píng)分(HR=3.230，95%CI：1.578-6.610，P=0.001)，見(jiàn)圖2，訓(xùn)練集與驗(yàn)證集中風(fēng)險(xiǎn)值的P值<0.05，說(shuō)明該模型屬于獨(dú)立預(yù)后因素。

圖2 訓(xùn)練集(a)和驗(yàn)證集(b)多因素Cox分析結(jié)果

2.6 風(fēng)險(xiǎn)模型的評(píng)價(jià)

對(duì)風(fēng)險(xiǎn)值分組進(jìn)行受試者工作特征曲線(receiver operating characteristic,ROC)檢驗(yàn)，繪制模型的ROC曲線，計(jì)算ROC曲線下面積數(shù)值A(chǔ)UC以評(píng)估模型的預(yù)測(cè)效能，在訓(xùn)練集中該模型預(yù)測(cè)1年、3年和5年生存率的AUC分別為0.9、0.796和0.863，在驗(yàn)證集中該模型預(yù)測(cè)1年、3年和5年生存率的AUC分別為0.601、0.686和0.761(見(jiàn)圖3，彩插2)，訓(xùn)練集、驗(yàn)證集的生存率曲線下面積都在60%以上，說(shuō)明能較準(zhǔn)確地預(yù)測(cè)預(yù)后。

2.7 對(duì)miRNA基因使用mitarbase數(shù)據(jù)庫(kù)進(jìn)行預(yù)測(cè)

通過(guò)對(duì)這些miRNA對(duì)應(yīng)的靶基因所在的基因通路、GO的生物過(guò)程ReactomePA通路進(jìn)行富集分析發(fā)現(xiàn)這些靶基因主要富集于調(diào)控細(xì)胞的分化、信號(hào)的轉(zhuǎn)導(dǎo)和凋亡等過(guò)程，例如：調(diào)節(jié)淋巴細(xì)胞和T細(xì)胞的分化，參與肽鏈蘇氨酸磷酸化等；在KEGG通路富集分析中，靶基因主要富集在乳腺癌、肝癌、MAPK信號(hào)通路和AMPK信號(hào)通路等癌癥相關(guān)通路中；另外，對(duì)疾病種類進(jìn)行富集分析，發(fā)現(xiàn)靶基因富集與腫瘤相關(guān)疾病，其中hsa-miR-150-3p主要富集于結(jié)直腸癌、乳腺癌及多發(fā)性骨髓瘤等常見(jiàn)的腫瘤中?；虮倔w論、通路、疾病富集分析結(jié)果都表明這5種miRNA的靶基因與腫瘤密切相關(guān)，因此這5種miRNA可能是潛在的腫瘤預(yù)后標(biāo)志物(見(jiàn)圖4，彩插2)。

3 討論

宮頸癌是世界范圍內(nèi)最常見(jiàn)的婦科惡性腫瘤之一，其發(fā)生率僅次于乳腺癌，居中國(guó)女性惡性腫瘤第二位[11]。宮頸癌被認(rèn)為是第一個(gè)病因明確的惡性腫瘤，幾乎所有宮頸癌都與HPV感染有關(guān)[12]。目前宮頸癌致癌機(jī)制仍未被闡明，研究證明miRNA在生長(zhǎng)發(fā)育、細(xì)胞增殖、凋亡等多種生物學(xué)過(guò)程中發(fā)揮重要調(diào)節(jié)作用，探索miRNA在宮頸癌發(fā)生、進(jìn)展、治療及預(yù)后等方面的作用，為尋找宮頸癌新的生物標(biāo)記物提供新的研究方向。

1993年，Lee RC等[13]在秀麗新小桿線蟲發(fā)育過(guò)程的研究中首次發(fā)現(xiàn)第1個(gè)miRNA，它是一類內(nèi)源性的短的非編碼miRNA，長(zhǎng)度約為18～25個(gè)核苷酸，可以調(diào)節(jié)大約60%的蛋白質(zhì)編碼基因[14]，參與了大多數(shù)關(guān)鍵細(xì)胞過(guò)程中基因表達(dá)的轉(zhuǎn)錄后調(diào)節(jié)，如細(xì)胞增殖、分化、凋亡等的調(diào)節(jié)[15]。miRNA的失調(diào)和許多人類疾病有關(guān)，特別是癌癥，一項(xiàng)研究表明，超過(guò)50%的miRNA基因位于癌癥相關(guān)的基因組區(qū)域或脆性位點(diǎn)[16]，成熟的miRNA與mRNA通過(guò)堿基互補(bǔ)配對(duì)，發(fā)揮其調(diào)節(jié)靶基因表達(dá)的作用，miRNA的異常表達(dá)通常會(huì)促使癌癥的進(jìn)展[17]。

近些年，隨著對(duì)miRNA研究的深入，其在癌癥的病因、進(jìn)展、轉(zhuǎn)移和預(yù)后的調(diào)節(jié)作用逐漸被人們重視，門婧睿等[18]通過(guò)構(gòu)建風(fēng)險(xiǎn)模型，最終確定hsa-miR132-3p、hsa-miR-139-5p、hsa-miR-3677-3p、hsa-miR-500A-3p可以作為預(yù)測(cè)肝癌預(yù)后的生物標(biāo)志物；miR-21已被證明與多個(gè)器官的惡性腫瘤有關(guān)，如在人類乳腺癌細(xì)胞系和組織中的表達(dá)增加，在發(fā)病的所有階段都起著關(guān)鍵作用，有助于預(yù)測(cè)宮頸癌的發(fā)生[19-20];在腎透明細(xì)胞癌中可以聯(lián)合miRNA-223評(píng)估淋巴結(jié)轉(zhuǎn)移、遠(yuǎn)處轉(zhuǎn)移及預(yù)測(cè)腎透明細(xì)胞癌的生存[9]。在乳腺癌的研究中還發(fā)現(xiàn)：miR-9、miR-10b和miR-17-5p可用于協(xié)助診斷；miR-148a和miR-335可評(píng)估乳腺癌的預(yù)后;miR-30c、miR-187和miR-339-5p則可以用來(lái)預(yù)測(cè)患者對(duì)特定治療的反應(yīng)[21]。除此之外，miRNA在胃癌、胰腺癌、結(jié)直腸癌等中均有大量研究[22-25]。

miRNA在宮頸癌中的研究也取得了一定的進(jìn)展，有研究證明miR127可作為宮頸癌淋巴結(jié)轉(zhuǎn)移的標(biāo)志物，miR-199a可能是宮頸癌治療的潛在靶點(diǎn)[26]；多種miRNA在宮頸癌中的表達(dá)異常可影響對(duì)化療和放療的敏感性[27]；本研究共篩選出5種miRNA：miR-101-3p、miR-532-5p、miR-150-3p、miR-378d、miR-378a-5p，這些miRNA 在宮頸癌中已有相關(guān)研究。miR-101-3p在宮頸癌中的表達(dá)水平與腫瘤分期及中晚期宮頸癌同步放化療抵抗相關(guān)，放療抵抗組表達(dá)水平較放療敏感組高[28]；miR-378已被證實(shí)在許多癌癥中都有異常表達(dá)[24,29]，在宮頸疾病的發(fā)生發(fā)展中，隨著宮頸病變的進(jìn)展，miR-378發(fā)揮致癌基因的作用，表達(dá)水平增加，可作為早期診斷的生物標(biāo)志物[29-32]；有研究證明miR-150能促進(jìn)HeLa和SiHA細(xì)胞的增殖、遷移和侵襲，在癌組織中表達(dá)水平增高[33]；miR-532在肝癌和宮頸癌中均有異常表達(dá)[34-35]。

本研究從TCGA數(shù)據(jù)庫(kù)下載 miRNA表達(dá)數(shù)據(jù)及患者臨床信息，通過(guò)統(tǒng)計(jì)學(xué)及生物信息學(xué)方法單因素和多因素Cox回歸分析、Lasso回歸分析篩選出5種與宮頸癌預(yù)后相關(guān)的miRNA，并構(gòu)建風(fēng)險(xiǎn)模型，該模型預(yù)測(cè)第1年、3年、5年生存率的AUC分別為：0.9、0.796、0.863，具有良好的預(yù)測(cè)性能。最后，本文作者對(duì)以上5種miRNA靶基因進(jìn)行了預(yù)測(cè)，發(fā)現(xiàn)它們所在基因通路、GO、Reactome通路，均與腫瘤密切相關(guān)。

綜上所述，本課題組構(gòu)建的宮頸癌患者生存的風(fēng)險(xiǎn)模型，具有較好的敏感度和特異性，但本研究是基于在線數(shù)據(jù)庫(kù)，預(yù)測(cè)模型還需要在實(shí)際的臨床操作中進(jìn)行驗(yàn)證，必要時(shí)需進(jìn)一步通過(guò)實(shí)驗(yàn)的方式證實(shí)其有效性。miRNAs作為基因表達(dá)的重要調(diào)節(jié)因子，有望作為新型分子標(biāo)志物應(yīng)用于宮頸疾病的臨床治療及科學(xué)研究中。