差異表達miRNA在胰腺癌預后判斷中的價值

高 杰，張 曉*，魏 超

(東南大學公共衛(wèi)生學院，江蘇 南京 210009)

胰腺癌是一種常見的消化系統(tǒng)高度惡性腫瘤，是中國致死率第7位的腫瘤，發(fā)病率呈現(xiàn)不斷上升的趨勢[1-2]。胰腺癌的診斷和治療手段雖然一直在進步，但效果仍然不盡如人意，胰腺癌早期難以診斷且患者經手術治療后的5年生存率低[3]。因此迫切需要確定可靠、敏感的生物標志物用于胰腺癌的早期診斷和治療[4-5]。微小RNA(microRNA，miRNA)是真核生物中存在的一類長度約為20個核苷酸的內源性非編碼RNA，可以與mRNA結合從而調控基因的表達[6]。目前已有的研究證實miRNA幾乎在所有生物途徑中都發(fā)揮著重要作用，與腫瘤的發(fā)生、發(fā)展密切相關，可以作為腫瘤診斷和治療的分子靶標[7-8]。目前miRNA已被廣泛應用于肺癌、胃癌、乳腺癌、肝癌等腫瘤研究中[9-11]，本研究應用基因芯片表達匯編數(shù)據(jù)庫(Gene Expression Omnibus，GEO)和腫瘤基因組圖譜數(shù)據(jù)庫(The Cancer Genome Atlas，TCGA)篩選胰腺癌預后相關的miRNA。

1 材料與方法

1.1 數(shù)據(jù)收集

從美國國立生物技術中心的GEO數(shù)據(jù)庫(https://www.ncbi.nlm.nih.gov/geo/)中下載胰腺癌GSE32678和GSE71533芯片數(shù)據(jù)集，這兩套數(shù)據(jù)集包括胰腺癌組織和正常組織的基因芯片表達數(shù)據(jù)[12]。GSE32678數(shù)據(jù)集包含胰腺癌組織樣本25例，正常組織樣本7例，GSE71533數(shù)據(jù)集包含胰腺癌組織樣本36例，正常組織樣本16例。從TCGA數(shù)據(jù)庫(https://cancergenome.nih.gov/)獲取178例胰腺癌樣本的miRNA表達及預后數(shù)據(jù)。

1.2 差異表達miRNA篩選

將GSE32678和GSE71533數(shù)據(jù)集中的數(shù)據(jù)進行標準化處理，使用R語言中的Limma軟件包(Limma package R 3.1.1)篩選差異表達miRNA。篩選條件為|log2(fold change)|>2且FDR(矯正后的P值)<0.05。結合GSE24279數(shù)據(jù)集[13]，獲得GSE32678、GSE71533和GSE24279數(shù)據(jù)集差異表達miRNA的交集。所篩選出的交集miRNA至少在2個及以上數(shù)據(jù)集中存在差異表達。對篩選出的差異表達miRNA繪制熱圖、火山圖和韋恩圖，熱圖的繪制使用R語言中的heatmap軟件包，火山圖的繪制使用R語言的ggplot2軟件包。

1.3 交集miRNA的生存分析

獲取TCGA數(shù)據(jù)庫中178例胰腺癌樣本的miRNA表達和預后數(shù)據(jù)，應用R語言中的Survival軟件包(R package survival)，采用Kaplan-Meier方法對篩選出的交集miRNA進行生存時間相關性分析，miRNA表達數(shù)據(jù)為自變量，生存時間為因變量，檢驗采用雙側log-Rank檢驗，P<0.05為差異有統(tǒng)計學意義。

1.4 miRNA潛在調控靶基因預測

應用靶基因預測數(shù)據(jù)庫miRanda(http://www.microrna.org/)、miRTarBase(http://mirtarbase.mbc.nctu.edu.tw/)、 miWalk(http://mirwalk.umm.uni-heidelberg.de/)和TargetScan(http://www.targetscan.org/)對交集miRNA進行潛在調控靶基因的預測，篩選出的靶基因需要同時包含于以上所有數(shù)據(jù)庫中，P<0.05為差異有統(tǒng)計學意義。

1.5 靶基因的功能分析

DAVID數(shù)據(jù)庫是一個生物信息數(shù)據(jù)庫，可以提供生物信息功能注釋信息[14]。使用DAVID(The database for annotation，visualization and integrated discovery)在線軟件(https://david.ncifcrf.gov/)對靶基因進行GO功能富集分析和KEGG通路分析，P<0.05為差異有統(tǒng)計學意義。

1.6 靶基因的蛋白互作網(wǎng)絡

STRING數(shù)據(jù)庫是一種預測蛋白質之間相互作用的軟件[15]。根據(jù)預測得到的靶基因，應用STRING在線軟件(https://string-db.org/)將靶基因所對應的蛋白質之間的相互作用關系可視化。

1.7 交集miRNA-mRNA調控網(wǎng)絡

根據(jù)胰腺癌3個芯片數(shù)據(jù)集的交集miRNA與靶基因之間的關系，應用Cytoscape 2.8.2軟件[16]繪制miRNA與mRNA之間的調控網(wǎng)絡圖。

1.8 miRNA-mRNA-pathway調控網(wǎng)絡構建

應用Cytoscape 2.8.2軟件將核心miRNA-mRNA對、mRNA-pathway對組合在一起，構建核心miRNA-mRNA-pathway調控網(wǎng)絡。

2 結果

2.1 胰腺癌中差異表達的miRNA

GSE32768數(shù)據(jù)集中共有100個miRNA在胰腺癌組織中表達呈顯著改變(|log2(fold change)|>2，P<0.05)，其中表達上調的miRNA有55個，表達下調的miRNA有45個，見圖1A與圖1C。GSE71533數(shù)據(jù)集中共有113個miRNA在胰腺癌組織中表達呈顯著改變[|log2(fold change)|>2，P<0.05]，其中表達上調的miRNA有34個，表達下調的miRNA有79個，見圖1B與圖1D。GSE24279數(shù)據(jù)集中共有59個miRNA在胰腺癌組織中表達呈顯著改變，與GSE32768、GSE71533數(shù)據(jù)集中差異表達的miRNA取交集，在兩個以上數(shù)據(jù)集中顯著表達的miRNA為交集miRNA，3個數(shù)據(jù)集中均存在差異表達的miRNA有4個，兩個數(shù)據(jù)集中均存在差異表達的miRNA有18個，共有22個交集miRNA，見圖1E。其中9個miRNA表達上調，13個miRNA表達下調。上調的miRNA為hsa-miR-10a-5p、hsa-miR-210-3p、hsa-miR-221-3p、hsa-miR-222-3p、hsamiR-125b-5p、 hsa-miR-199a-3p、 hsa-miR-23a-3p、 hsa-miR-31-5p、 hsa-miR-135b-5p； 下調的miRNA為hsa-miR-148a-3p、hsa-miR-216a-5p、hsamiR-216b-5p、hsa-miR-217、hsa-miR-29b-3p、hsamiR-29c-3p、hsa-miR-30a-5p、hsa-miR-30c-5p、hsa-miR-335-5p、 hsa-miR-375、 hsa-miR-551a、hsa-miR-1184、hsa-miR-552-3p。

圖1 差異miRNA表達圖(聚類圖、火山圖、韋恩圖)

2.2 交集miRNA的生存分析

結合TCGA數(shù)據(jù)庫178個胰腺癌樣本的miRNA表達及預后數(shù)據(jù)對22個交集miRNA進行Kaplan-Meier生存分析。結果顯示hsa-miR-30a-5p、hsa-miR-551a、hsa-miR-375與胰腺癌患者的生存預后密切相關。其中hsa-miR-30a-5p和hsa-miR-551a低表達組的生存時間均長于高表達組生存時間(P<0.05)，見圖2A與圖2C。hsa-miR-375高表達組生存時間長于低表達組生存時間(P<0.05)，見圖2B。

圖2 miRNA高表達組和低表達組的Kaplan-Meier生存分析

2.3 靶基因的功能分析

通過miRanda、miRTarBase、miWalk和TargetScan數(shù)據(jù)庫對22個交集miRNA進行靶基因預測，獲得645個靶基因，并進行GO與KEGG pathway分析。圖3是排名前20的靶基因顯著性GO功能和KEGG通路圖。靶基因主要富集于基因的轉錄和表達、細胞增殖、信號轉導、蛋白質磷酸化等生物學過程和功能(P<0.05)，見圖3A。靶基因主要參與癌癥通路、MicroRNAs通路、肝細胞癌通路、細胞衰老通路等通路(P<0.05)，見圖3B。

圖3 交集miRNA靶基因的功能分析

2.4 胰腺癌靶基因的蛋白互作網(wǎng)絡

將預測得到的靶基因數(shù)據(jù)上傳到STRING數(shù)據(jù)庫中，獲得靶基因對應的蛋白質之間的相互作用。蛋白質之間通常是相互作用的，兩個或兩個以上的蛋白質會通過非共價鍵形成蛋白復合體，蛋白質之間的相互作用越多，則該蛋白質越重要。CD44、STAT3、STK11、AKT2、FOXO3靶基因對應的蛋白質互作關系密切，這些靶基因編碼的蛋白質可能會影響胰腺癌的放療敏感性，進而對預后產生影響，見圖4。

2.5 胰腺癌交集miRNA-mRNA調控網(wǎng)絡

根據(jù)交集miRNA與其靶基因之間的靶向調控關系，應用Cytoscape 2.8.2軟件構建的交集miRNA-mRNA調控網(wǎng)絡如圖5所示。該圖將交集miRNA與靶基因之間的關系可視化。miRNA調控的下游靶基因數(shù)目越多說明該miRNA越重要。將miRNA-mRNA調控網(wǎng)絡中miRNA所調控的下游基因數(shù)目前5的miRNA篩選為核心miRNA[17]，分別是：hsa-miR-125b-5p、hsamiR-221-3p、hsa-miR-31-5p、hsa-miR-222-3p、hsa-miR-10a-5p。

2.6 miRNA-mRNA-pathway調控網(wǎng)絡

基于上述5個核心miRNA，繪制核心miRNA-mRNA-pathway調控網(wǎng)絡。其中節(jié)點基因較重要的有RAF1、MAPK10、FOXO3、TIAM1、BCL2等，且與癌癥通路、MicroRNAs通路、河馬信號通路、細胞衰老通路、趨化因子信號通路等顯著相關。

圖4 胰腺癌靶基因的蛋白互作網(wǎng)絡

3 討論

圖5 胰腺癌交集miRNA-mRNA調控網(wǎng)絡

圖6 胰腺癌miRNA-mRNA-pathway調控網(wǎng)絡

胰腺癌通常發(fā)病隱匿，難以診斷和治療，容易發(fā)生轉移，治療效果差，5年生存率僅1.2%~6.0%[18-19]?；?、放療、手術治療等傳統(tǒng)療法在提高胰腺癌患者生存率方面發(fā)揮的作用很小[20-21]，所以尋找胰腺癌診斷、治療、預后相關的分子靶標至關重要。許多研究表明miRNA可以作為腫瘤診斷和預后的生物標志物[22]。目前有一些研究表明特異性miRNA可用作胰腺癌預后和診斷的標志物，但是每個研究所包含的研究患者很少，可重復性很差[23]，也有一些基于TCGA數(shù)據(jù)庫的研究[24]，TCGA數(shù)據(jù)庫擁有豐富的胰腺癌患者基因測序數(shù)據(jù)和相應臨床數(shù)據(jù)，TCGA數(shù)據(jù)庫中胰腺癌組織樣本有178例，但正常組織樣本僅4例[25]，正常組織樣本過少，容易造成偏倚。本研究所選取的GEO數(shù)據(jù)庫中GSE32678和GSE71533數(shù)據(jù)集共含有61例胰腺癌樣本和23例正常組織樣本，樣本豐富。此外本研究構建了胰腺癌miRNA-mRNA-pathway調控網(wǎng)絡，既能看出miRNA對mRNA的調控作用，又可以看出其對通路的調控作用。

本研究應用生物信息學方法篩選出22個交集miRNA，其中hsa-miR-217可以抑制胰腺癌細胞的侵襲能力[26]，與患者的生存顯著相關。hsa-miR-30c-5p與胰腺癌預后的相關關系尚未見報道，但有研究表明其與乳腺癌的生存期密切相關[27]。

接著采用TCGA數(shù)據(jù)庫中178例胰腺癌樣本的miRNA表達及預后數(shù)據(jù)對交集miRNA進行生存分析。生存分析發(fā)現(xiàn)hsa-miR-375、hsa-miR-30a-5p、hsamiR-551a與胰腺癌患者的預后相關。hsa-miR-375在胰腺癌中表達上調，可以抑制腫瘤細胞的增殖、遷移，與患者生存顯著相關[28-29]。已有研究表明hsamiR-30a-5p在多種癌癥中異常表達[30-31]且可以作為胰腺癌預后的標志[24]。hsa-miR-551a在胰腺癌中的預后研究尚未見報道，但是研究[32]表明其可以顯著抑制胃癌細胞的侵襲和遷移。

進一步采用miRanda等數(shù)據(jù)庫對交集miRNA進行靶基因預測并進行功能分析。功能分析發(fā)現(xiàn)胰腺癌的差異表達miRNA可能與基因的轉錄和表達、細胞增殖等功能有關，與癌癥通路、MicroRNAs通路等有關，這些都與腫瘤的治療效果密切相關。更進一步構建靶基因對應蛋白質的蛋白互作網(wǎng)絡，將靶基因對應的蛋白質之間的互作關系可視化，CD44、STAT3、STK11、AKT2、FOXO3靶基因對應的蛋白質互作關系密切，相應蛋白質可能會影響胰腺癌放射敏感性，進而影響預后。

接著繪制交集miRNA-mRNA調控網(wǎng)絡，通過可視化的miRNA-mRNA調控網(wǎng)絡我們發(fā)現(xiàn)hsa-miR-125b-5p、hsa-miR-221-3p、hsa-miR-31-5p、hsamiR-222-3p、hsa-miR-10a-5p是miRNA-mRNA調控網(wǎng)絡中的核心miRNA。已有研究表明hsa-miR-125b-5p的表達與胰腺癌的侵襲性呈正相關[33]。hsa-miR-221-3p和hsa-miR-31-5p在胰腺癌預后中的作用尚未見報道，但研究表明hsa-miR-221-3p可以促進宮頸鱗狀細胞癌的轉移且可作為宮頸鱗狀細胞癌診斷和治療的生物標志物[34]，此外也有研究表明hsa-miR-31-5p在結直腸癌中高表達，在癌細胞增殖、遷移、侵襲和凋亡過程中發(fā)揮重要作用[35]。hsa-miR-222-3p與胰腺癌預后的關系尚未知曉，但有研究表明hsa-miR-222-3p可參與口腔鱗狀細胞癌的各種癌癥途徑[36]，可以作為口腔鱗狀細胞癌預后的生物標志。已有研究表明hsa-miR-10a-5p可以促進胰腺癌的轉移，與患者生存緊密相關[37]。

最后根據(jù)miRNA-mRNA調控網(wǎng)絡確定5個核心miRNA并構建miRNA-mRNA-pathway調控網(wǎng)絡。通過調控網(wǎng)絡，我們發(fā)現(xiàn)RAF1、MAPK10、FOXO3、TIAM1、BCL2為節(jié)點基因前5位。已有研究表明這5個基因可以作用于相應的信號轉導通路，發(fā)揮原癌基因或抑癌基因的作用[38-42]。

綜上所述，通過對GEO和TCGA數(shù)據(jù)庫的胰腺癌相關數(shù)據(jù)的分析，所篩選出的hsa-miR-30a-5p、hsamiR-551a、hsa-miR-375與胰腺癌患者的預后密切相關，可以作為胰腺癌預后的判斷標志。hsa-miR-125b-5p、hsa-miR-221-3p、hsa-miR-31-5p、hsamiR-222-3p、hsa-miR-10a-5p是miRNA-mRNA調控網(wǎng) 絡 中 的 核 心 miRNA。 RAF1、MAPK10、FOXO3、TIAM1、BCL2對胰腺癌患者預后的影響值得進一步關注。本研究的局限性在于所有分析均通過生物信息學方法完成，未進行相關實驗進行驗證，這些都有待在進一步的工作中繼續(xù)完成。