前列腺癌差異表達基因的分析及其miRNA和lncRNA的預測*

宋詠剛，余鵬 ，胡文慧，石玉玲，趙雪，胡祖權，王赟，曾柱,3**

(1.貴州醫科大學 基礎醫學院，貴州 貴陽 550025；2.貴州醫科大學 貴州省免疫細胞與抗體工程研究中心 生物與醫學工程重點實驗室，貴州 貴陽 550025；3.貴州醫科大學 生物與工程學院，貴州 貴陽 550025)

前列腺癌是全球第二大男性癌癥，同時也是全球第四大癌癥，有著極高的發病率和死亡率[1]。前列腺癌是一種無癥狀的腫瘤，在整個癌癥發生過程中沒有任何癥狀，通過臨床手段很難檢測。前列腺癌不斷增長時，會伴隨膀胱頸梗阻、癌細胞轉移等癥狀，這時通常不可能進行治愈性治療[2]。因此，關于前列腺癌的腫瘤標志物的研究至關重要要。miRNA和lncRNA是基因表達調控的重要組成成分，特別是在癌細胞中表達顯著，它能幫助理解疾病發生的分子機制[3]。前列腺癌作為一種遺傳性疾病，國內外對其基因水平的調控機制和新分子靶點的研究日益增加[4]，例如基于自身免疫療法的Sipuleucel-T前列腺癌疫苗，通過將VISTA或PARP等分子體外激活樹突狀細胞(dendritic cells，DCs)、將激活的細胞回輸到體內產生免疫應答而治療癌癥[5]。本研究從生物信息學角度出發，通過將基因表達數據庫(gene expreesion omnibus，GEO)與癌癥和腫瘤圖譜中心(the cancer genome atlas,TCGA)數據庫關于前列腺癌的樣本數據相結合，篩選前列腺癌發生過程中的關鍵基因及其上下游的調控分子，為臨床診斷和臨床靶向藥物的研發提供理論基礎。

1 資料與方法

1.1 基因芯片數據的獲取

登錄美國國家生物技術信息中心(national centerfor for biotechnology information,NCBI)的GEO數據庫(http://www.ncbi.nlm.nih.gov/，GEO database)，檢索Affymetrix芯片數據GSE55945數據包含13例前列腺癌組織和8例正常組織的樣本數據；登錄TCGA數據庫檢索前列腺癌的轉錄組數據并下載，包含504例前列腺癌組織和55例癌旁組織的樣本數據。

1.2 差異基因的分析

1.2.1差異基因的分析 使用RStudio軟件‘affy’包對GEO原始數據文件進行標準化處理，并檢驗歸一化處理的結果。使用‘limma’包[6]篩選表達有顯著差異的基因(differentially expressed genes,DEGenes)，以|logFC|≥2(fold change，FC)并且P<0.05 作為標準進行篩選[7]。通過RStudio軟件‘edgR’包對TCGA的前列腺癌數據進行分析[8]，同樣以|logFC|≥2并且P<0.05作為標準進行篩選。

1.2.2DEGenes的GO和KEGG富集分析和蛋白質相互作用網絡的構建 登錄KOBAS(http://kobas.cbi.pku.edu.cn/)：對篩選出的DEGenes的進行GO分析和京都基因和基因組百科全書(kyoto encyclopedia of genes and genomes ,KEGG)途徑分析[9]；登錄STRING在線工具(https://string-db.org/)[10]：對篩選出的DEGenes構建蛋白質相互作用網絡(protein-protein interaction network,PPI)[11]。將結果導入CYTOSCAP軟件，使用‘MCODE’插件[12]對該相互作用網絡進行K 核解析(k-core decomposition)，篩選出DEGenes的核心模塊。篩選條件為P<0.05并且基因數目>5[13]。

1.2.3核心基因相應miRNA的預測 將篩選出的核心DEGenes導入mirDIP數據庫，對核心DEGenes的miRNA進行預測[14]。

1.2.4miRNA的生存分析及lncRNA的預測 登錄STARBASE(http://starbase.sysu.edu.cn/)，使用Pan-Cancer工具對預測出的miRNA進行生存分析，選擇與前列腺癌患者生存率相關(P<0.05)的miRNA作為核心miRNA。對符合條件的miRNA進行下游lncRNA的預測，為提高預測結果的準確性，選擇至少在一種以上的腫瘤中被證實(P<0.05)的lncRNA作為結果。

1.2.5基因mRNA-pathwang網絡的構建 將符合標準的mRNA、miRNA、lncRNA以及重要的信號通路及之間的關系導入CYTOSCAPE軟件中進行可視化，建立mRNA-pathwang網絡圖。

2 結果

2.1 由 RStudio軟件獲得DEGenes

芯片數據(GSE55945)經分析獲得正常組織相對于前列腺癌組織的差異表達基因共431個(|logFC|≥2,P<0.05)。TCGA數據庫關于前列腺癌的樣本經分析，獲得癌旁組織相對于前列腺癌組織的差異表達基因共2 472個(|logFC|≥2,P<0.05)，部分差異表達基因及其對應信息見表1與表2。利用RStudio軟件‘venn’包對兩個數據庫中篩選出的差異基因使用Merge函數，得到兩個數據庫共同的137個DEGenes，并構建維恩圖(圖1)。

表1 GEO數據庫正常組織與前列腺癌組織的部分差異表達基因Tab.1 The partial DEGenes of prostate cancer tissue to normal tissue in GEO database

表2 TCGA數據庫中正常組織與前列腺癌組織的部分差異表達基因Tab.2 The partial DEGenes of prostate cancer tissue to normal tissue in TCGA database

圖1 共同差異表達基因的維恩圖Fig.1 Venn diagram of common DEGenes

2.2 GO和KEGG富集分析DEGenes

使用KOBAS網站對篩選出的137個DEGenes進行GO分析和KEGG分析，GO分析結果顯示有23個GO富集(表3)，KEGG富集分析有20個通路被富集(表4)。

2.3 構建DEGenes 的PPI和核心基因相應miRNA的預測

通過STRING在線分析工具繪制DEGenes的PPI(圖 2)，結果用CYTOSCAPE軟件的“MCODE”工具對該PPI進行K 核解析，以節點數>5、所得分數>0.4篩選出子網絡，對其中最穩定的子網絡作為核心模塊進行分析，核心模塊包含7個節點基因和 21條相互作用關系對(圖 3)。使用mirDIO數據庫對7個核心基因預測得到12個對應的miRNA，通過CYTOSCAPE軟件對結果進行可視化(圖3)。

2.4 miRNA的生存分析及lncRNA的預測

登錄STARBASE 網站對預測出的12個miRNA進行生存分析，其中has-let-7a-5p和has-miR-98-5是與前列腺癌患者總體存活率相關的miRNA(P<0.05，圖4、表5)；通過Pan-Cancer工具驗證其對應的基因在前列腺癌患者的表達量，AURKB、BUB1B和EZH2在前列腺癌患者體內高表達(圖5)。

2.5 構建lncRNA-mRNA pathway 網絡

綜合以上結果，使用CYTOSCAPE軟件對預測出的lncRNA、與前列腺癌患者生存率有關的miRNA、核心DEGene以及相關通路構建網絡圖(圖6)，該網絡涉及3個重要的信號通路、3個核心基因、2個miRNA和4個lncRNA。見圖6。

3 討論

近年來隨著分子生物學、測序和數據分析技術的發展，生物信息學被廣泛地應用在各種腫瘤標志物或治療靶點的挖掘中，它可以幫助克服癌癥的傳統療法中診斷延遲和主觀不可靠評估等缺點[15]，特別在前列腺癌的診斷中，其黃金標準仍會有20%以上的假陰性率[16]。本研究通過將TCGA數據庫和GEO數據庫中前列腺癌的樣本數據相結合，挖掘出前列腺癌高度相關的差異表達基因137個。對這些基因進行K-核解析后得到7個核心(KIF4A、AURKB、BUB1B、TPX2、NCAPG、HJURP和EZH2)基因，其中BUB1B與前列腺癌的潛在致病風險有關，BUB1B主要參與細胞周期的生物過程，通過與CCNB1、CDK1等基因協同作用，影響前列腺癌患者的疾病進展[17]。SHA等[18]發現PRKAR2B主要加速細胞周期生物學過程并調節多個細胞周期基因，通過過度表達促進趨勢抵抗性前列腺癌(castration-resistant prostate cancer,CRPC)增殖和侵襲，并抑制CRPC細胞的凋亡。組蛋白甲基轉移酶EZH2d的失調，會導致腫瘤抑制性miRNA的表觀遺傳學沉默而增強前列腺癌細胞增殖和轉移過程中的致癌作用[19]。

表3 差異表達基因的GO分析Tab.3 GO analysis of DEGenes

注：A為核心模塊，B為核心模塊相互作用網絡圖。圖2 差異表達蛋白的核心模塊與核心模塊相互作用網絡Fig.2 The central modules of PPI network of DEGenes constoucted by using CYTOSCAPE and core module

圖3 DEGenes 及其對應的miRNAFig.3 DEGenes and corresponding miRNA

注:A為has-miR-21-5p，B為has-miR-98-5p。圖4 miRNA的生存分析Fig.4 Survial analysis of miRNA

miRNAlncRNA癌癥類型Phsa-let-7a-5pXISTUrothelialbladdercancer(BLCA)0.034798hsa-let-7a-5pKCNQ1OT1Headandnecksquamouscellcarcinoma(HNSC)0.0002hsa-let-7a-5pNEAT1Clearcellkidneycarcinoma(KIRC)1.21E-04hsa-miR-98-5pTTTY15Breastcancer(BRCA)0.046063hsa-miR-98-5pKCNQ1OT1Headandnecksquamouscellcarcinoma(HNSC)0.0002hsa-miR-98-5pXISTLungadenocarcinoma(LUAD)0.019418hsa-miR-98-5pNEAT1Cutaneousmelanoma(SKCM)3.06E-02

注:A為BUB1B，B為EZH2，C為MYC。圖5 差異基因在前列腺癌中的表達Fig.5 DIfferential expression of DEGenes in prostatic cancer

圖6 lncRNA-mRNA pathway 網絡Fig.6 Network of lncRNA-mRNA pathway

通過對篩選137個DEGenes的GO功能分析和KEGG富集分析發現，這些DEGenes主要涉及細胞粘附，細胞周期和機體炎癥反應。CAV1是Src激酶的重要底物，能通過與CSD之間的相互作用抑制粘著斑的形成，引起癌細胞的運動[20]。CD40/CD40L主要參與血管生成，在前列腺腫瘤腺體組織中，CD40在激素難治性期分泌較高，能引起慢性炎癥并導致浸潤區組織損傷，對前列腺癌進展至關重要[21]。TPM1是一種具有促凋亡功能的肌動蛋白結合腫瘤抑制因子，但在前列腺癌組織中，TPM1發生特異性的可變剪接，通過改變TPM1的結構引起肌動蛋白的失調，導致腫瘤的侵襲和轉移，這種作用在腫瘤晚期尤其明顯[22]。

非編碼RNA在體內不能被翻譯為蛋白質，卻參與多種細胞過程，主要包括miRNA和lncRNA，由于具有時空或組織特異性，被廣泛應用于分子靶點或腫瘤標準物的研究中[23]。本研究通過mirDIP數據庫預測核心DEGenes的上游miRNA，并將篩選出的miRNA通過STARBASE數據庫中結合臨床樣本進行生存分析，結果發現hsa-let-7a-5p和hsa-miR-98-5p與前列腺癌患者的總體生存率相關，并且它所對應的3個下游基因AURKB、BUB1B和EZH2在前列腺癌患者體內高表達。其中，hsa-let-7a-5p已被作為I期肺癌的潛在生物標準物之一[24]。hsa-miR-98-5p主要通過下調NEAT1的水平引起非小細胞肺癌的發生發展[25]。對以上2種miRNA進行lncRNA預測，共發現4種lncRNA(XIST、KCNQ1OT1、NEAT1、TTTY15)，其中，XIST和KCNQ1OT1能通過促進細胞遷移和上皮-間質細胞轉化引起腫瘤的發展[26-27]。TTTY15來定位于Y染色體上，在前列腺癌細胞中高表達，臨床研究證實，對TTTY15和USP9Y的測定能預測前列腺癌活檢的結果[28]。

本研究從生物信息學角度出發，對前列腺癌的差異基因進行篩選和分析，篩選出前列腺癌發生過程中的關鍵基因及其上下游的調控分子，構建了一條包含lncRNA、miRNA以及mRNA的分子作用網絡，為臨床診斷和臨床靶向藥物的研發提供理論基礎。