乳腺癌預后的相關基因生物信息學分析

左 　然 任曉虎 吳德生李 　萍 吳 　雯 謝 　妮袁建輝* 讓蔚清*

（1．南華大學公共衛生學院，湖南　衡陽　421001；2．深圳市疾病預防控制中心毒理研究所，廣東 　深 圳 　518055；3．深圳市第二人民醫院轉化醫學研究院，廣東 　深 圳 　518000）

乳腺癌是全世界女性最常見的侵襲性腫瘤之一，其發病率高居女性腫瘤榜首，對女性健康造成極大威脅[1-2]。乳腺癌的發病機制目前尚不明確，遺傳、雌激素水平、生長因子、細胞因子、激酶和非編碼RNA等均在腫瘤細胞生長中起重要作用[3-4]。microRNA是一類重要的內源性非編碼RNA，主要通過轉錄后水平調控基因表達[5]。microRNA和表觀遺傳調控過程之間存在復雜的反饋網絡并形成表觀遺傳學microRNA調控回路[6]。隨著對microRNA功能研究的不斷涌現，越來越多的研究表明microRNA與腫瘤之間存在密切關系。有報道發現，microRNA-214-3p通過抑制靶基因MELK表達從而抑制肝細胞癌進展[7]；而miR-137下調和肝癌淋巴結轉移被發現與不良預后密切相關[8]；另外microRNA-200可影響管腔祖細胞增殖進而促進乳腺癌生長和轉移[9]。盡管microRNA及其潛在靶向調控的基因已被廣泛應用于腫瘤研究，但類似的分子事件在乳腺癌預后中鮮有涉及。本研究通過挖掘腫瘤基因組圖譜計劃(the cancer genome atlas，TCGA)中的乳腺癌基因表達譜和microRNA測序數據，尋找乳腺癌預后相關的關鍵性分子事件。

本研究中，我們利用van Iterson等[10]提出的全局算法整合基因表達譜和microRNA測序數據，并預測給定臨床背景下的microRNA靶基因。該算法用到的統計學檢驗可以評估Pearson相關系數出現假陽性的情況，同時基于全局檢驗通過Lasso回歸預測數據庫中所有可能的目標mRNA[11-12]。為了進一步分析microRNA及相關靶基因在乳腺癌中起到的作用，我們利用通用基因集富集算法(general applicable gene set enrichment，GAGE)[13]分析microRNA及其靶基因涉及的異常變化通路。進一步證實所篩選的microRNA在乳腺癌組織和細胞中能負向調節其靶基因，并使用Cox回歸風險模型比例風險模型評估其作為預后標志物的潛力。本研究分別對乳腺癌相關microRNAs和基因的表達數據進行分析。

1　材料與方法

1.1　數據準備

從TCGA數據庫中(http://cancergenome.nih.gov/)下載經Lowess標準化和對數轉化的乳腺癌3級基因表達數據和microRNA測序數據。將所有的基因符號轉換為Entrez Gene ID后，整理成矩陣并導入結構調查語言(structure query language，SQL)數據庫中。

1.2　乳腺癌中的基因與microRNA表達的定量分析

將來自120對癌-癌旁組織的基因表達譜數據以及204對癌-癌旁組織的microRNA測序數據整理、優化后，利用R語言Limma軟件包(Limma package R 3.1.1)通過整合線性模型、經驗貝葉斯理論與傳統t檢驗對基因表達譜數據進行定量及統計分析。利用Holm-Bonferroni方法控制多次統計學檢驗可能產生的假陽性，將假陽性率(false discovery rate，FDR)控制在0.01。差異表達基因的篩選標準為：FDR值(即矯正后的P值)小于0.01，差異表達絕對值變化大于4。對于異常改變的microRNA，篩選標準為FDR值(即矯正后的P值)小于0.01，差異表達絕對值變化大于1.5[14-15]。

1.3　乳腺癌中microRNAs潛在調控靶基因的預測

收集388例乳腺癌和388例配對癌旁組織的基因表達譜和microRNA表達數據。從miRTarBase(http://mirtarbase.mbc.nctu.edu.tw/)，MSigDBC3 motif Gene set(http://www.broadinstitute.org/)，PITA(http://genie.weizmann.ac.il/pubs/mir07/mir07_data.html)，Microcosm(http://www.ebi.ac.uk/enright-srv/microcosm/htdocs/targets/v5/)和TargetScan(http://www.targetscan.org/) 數據庫獲得microRNA靶基因。將來自388名患者基因譜表達數據和microRNA測序數據整理為矩陣格式。利用van Iterson等提出的基于線性(或廣義線性)回歸模型的全局無效假設的microRNA靶基因預測算法分析microRNA的潛在靶基因，該方法旨在提供一種有效的方法來評估多種疾病之間的定量比較，并預測每個潛在的microRNA-mRNA的關聯性[10]。該算法在R軟件microRNA-mRNA包中實現(http://www.humgen.nl/Microarray Analysis Group.html)。選擇同時在2個或以上數據庫中預測到的microRNA靶基因作為基礎數據庫。將樣本隨機分為訓練集(所有樣品的2/3)和測試集(所有樣品的1/3)，把microRNA預測靶基因對引入線性回歸模型的全局無效假設方法應用于訓練集，找出關聯并在測試集中進行驗證。

1.4　基因集富集分析

傳統定量分析關注單個基因顯著變化卻忽視了多個基因之間的協同作用，因此通用基因集富集算法(GAGE)是對特定基因集合(如特定信號通路)的變化進行分析，而不考慮單個基因改變。將基因表達譜數據整理成矩陣(基因表達為行，樣本分組為列)導入R軟件，通過GAGE包進行分析。基因集來自MSigDB C2 Curated Gene set(http://www.broadinstitute.org/)。不同于傳統的基因集富集，GAGE算法允許基因同時向不同方向(上調或下調)變化，因此更適合于通路分析。Benjamini和Hochberg方法來控制多重測試的假陽性率，選擇校正后的P值(即Q值)＜0.05基因集進一步分析。基因集的大小(包括在基因集中的基因的數目)設定在50～500個之間，選擇差異變化大于1倍標準差的基因集。

1.5　Cox回歸風險模型篩選乳腺癌潛在評估預后的分子靶標基因

利用R survival包(R package survival)分析microRNA及潛在靶基因在乳腺癌預后中的作用。以microRNA測序數據及基因表達譜數據作為主要自變量，生存時間為因變量，相關臨床資料作為次要自變量，使用Cox回歸風險模型進行迭代分析，Kaplan-Meier進行統計學檢驗。

2　結果

2.1　乳腺癌中異常表達的基因和microRNA

對基因表達譜和microRNA測序數據的分析結果顯示，17 533個基因中只有344個基因表達呈顯著改變(矯正后p＜0.01，差異表達絕對值變化大于4)，見圖1A和表1。而microRNA中，135個在乳腺癌組織中表達發生異常改變化(矯正后p＜0.01，差異表達絕對值變化大于1.5)，見圖1B和表2。

圖1　mRNA和microRNA表達譜的火山圖

表1　microRNA及其潛在靶基因實驗驗證

2.2　乳腺癌中基于表達數據的microRNAs潛在靶基因預測

對338對乳腺癌-癌旁組織的microRNA及基因表達譜數據進行關聯性研究發現，31個microRNA與6 429個預測靶基因顯著相關。其中有部分microRNA對基因表達呈正相關，盡管研究報道發現microRNA可誘導上調靶mRNA對細胞周期阻滯[16-17]，亦可促進結合核糖體mRNA基因的5′UTR和增強翻譯[18]，但目前尚無證據表明在人群中microRNA上調能引起mRNA表達增加，因此正相關可能由于microRNA對基因的間接調控作用。然而為保證結果的可靠性，我們僅選擇與靶基因表達呈負相關的microRNA進行驗證。利用另外60對癌-癌旁的microRNA和基因表達譜數據，對14個microRNA和67個負相關的潛在靶基因進行驗證。結果顯示，其中6個microRNA表達上調，其潛在靶基因表達下調，而8個microRNA表達下調，其潛在調控靶基因上調。

2.3　microRNA及其潛在調控靶基因的通路富集

為了進一步探討microRNA及相關靶基因在乳腺癌中可能的作用機制，基因集富集算法對microRNA測序及基因表達譜數據進行富集分析，發現了8個乳腺癌相關microRNA及其調控的31個關鍵潛在調控靶基因參與139條乳腺癌富集通路(見圖2)。

表2　乳腺癌相關臨床資料的生存分析

2.4　乳腺癌潛在預后標志物

生存曲線顯示，乳腺癌患者從初始診斷起的6年，生存率明顯下降(2 000 d，圖3A)。為避免樣本變異和未記錄樣本的選擇偏倚，選擇有完整記錄的樣本(有臨床結局和治療方案)進行分析。由于某些類別的腫瘤分期(IV，Tis和X)樣本量非常少，因此未被納入生存分析。Cox回歸風險模型顯示，腫瘤邊緣狀態、更年期狀態和淋巴結數量等相關臨床資料對患者預后無明顯影響。而放射治療顯示腫瘤邊緣狀態呈保護作用(HR=0.432，P=0.058，圖3B)。進一步深入分析發現，分泌卷曲相關蛋白1(SFRP1)高表達呈保護作用(HR=0.9，P=0.015，圖3C)，has-miR-342-5p既沒有保護作用，也沒有危害效應(HR=0.99，P=0.144，圖3D)。然而，交互作用研究顯示，has-miR-342-5p抑制靶基因SFRP1時提示不良預后(見表3)。另外，hasmiR-342-5p對生長特異型抑制基因2(GASP2)的抑制作用則代表不良預后(Cox多因素分析，P=0.016，HR=18.88)。

圖3　乳腺癌患者的生存曲線

表3　乳腺癌潛在預后標志物

3　討論

在本研究中，我們將microRNA與其靶基因的關聯性研究與通路富集分析整合后篩選出有明顯改變的microRNAs及其潛在靶基因。此外，我們對所有預測的microRNA-靶基因中，發現9對microRNA-靶基因已在其他研究中被實驗驗證(表1)。其中7個microRNA及其30個潛在靶基因被發現在139個乳腺癌異常通路中發揮關鍵作用。Cox回歸風險模型分析揭示了來自這些microRNA及其靶基因參與的分子事件可作為乳腺癌重要的預后特征。使用Cox比例風險分析，SFRP1、CLDN19及KRT5在患者生存中被發現有保護作用。已知SFRP1是與乳腺癌細胞生長負相關的腫瘤抑制因子[19]。據報道，乳腺癌中SFRP1的失活與不良預后有關[20]。此外，SFRP1僅在2 400～3 500 d內表示出顯著的危害效應，而在3 500 d后危害效應降低。然而，由于miR-342-5p的調節，SFRP1的保護作用受到了影響。

目前已有一些研究表明microRNA與乳腺癌的相關性，在乳腺癌的發生過程中起到促癌或抑癌的作用。miR-342-5p是一種多功能的microRNA，參與細胞增殖、分化、信號轉導等許多重要的生物學過程。據報道miR-342-5p在由MCF-7細胞構成的三維組織模型中通過靶向DNA結合抑制劑4(ID4)和DNA甲基轉移酶1(DNMT1)來調節CEACAM1誘導的管腔形成[21]。此外，miR-342-5p還通過調節Notch信號傳導參與神經元干細胞分化[22]，并可作為阿爾茨海默癥的潛在血漿生物標志物[23]。其他研究表明，miR-342-5p可能通過下調阿爾茨海默癥轉基因模型小鼠中的錨蛋白G來促成軸突病變[24]。

SFRP1是富含半胱氨酸結構域的SFRP家族成員，可能被Wnt結合。同時SFRP1還介導視網膜中感光細胞的極性[20]。CLDN19是一個重要的膜蛋白，在細胞間緊密連接中起關鍵作用[25]，CLDN19的異常表達可促進癌前病變的細胞向惡性方向發展[26]。也有報道CLDN19 在乳腺癌中表達異常降低[27]，與我們的研究結果一致。

在乳腺癌中microRNA靶向調控特定基因作為預后關鍵分子事件的研究還相對較少。因此我們的研究不僅揭示了miR-342-5p、SFRP1、CLDN19以及KRT5在乳腺癌發生發展中具有的重要生物學功能，同時還表明miR-342-5p可能通過抑制乳腺癌中的SFRP1、KRT5以及CLDN19來提示不良預后，這為miR-342-5p相關的分子事件作為乳腺癌臨床監測的潛在指標提供了科學依據。

[1] SIEGEL R L，MILLER K D，JEMAL A. Cancer statistics，2017[J]. CA Cancer J Clin，2017，67(1)：7-30.

[2] CHEN W，ZHENG R，ZHANG S，et al. Cancer incidence and mortality in China，2013[J]. Cancer Lett，2017，401：63-71.

[3] JORDAN C T，GUZMAN M L. Mechanisms controlling pathogenesis and survival of leukemic stem cells[J]. Oncogene，2004，23(43)：7178-7187.

[4] TSATSANIS C，SPANDIDOS D A. The role of oncogenic kinases in human cancer (Review)[J]. Int J Mol Med，2000，5(6)：583-590.

[5] CHEN K，RAJEWSKY N. The evolution of gene regulation by transcription factors and microRNAs[J]. Nat Rev Genet，2007，8(2)：93-103.

[6] MOUTINHO C，ESTELLER M. MicroRNAs and epigenetics[J].Adv Cancer Res，2017，135：189-220.

[7] LI Y，LI Y，CHEN Y，et al. MicroRNA-214-3p inhibits proliferation and cell cycle progression by targeting MELK in hepatocellular carcinoma and correlates cancer prognosis[J].Cancer Cell Int，2017，17：102.

[8] CUI S， SUN Y， LIU Y， et al. MicroRNA137 has a suppressive role in liver cancer via targeting EZH2[J]. Mol Med Rep，2017，16(6)：9494-9502.

[9] SANCHEZ-CID L， PONS M， LOZANO J J， et al.MicroRNA-200， associated with metastatic breast cancer，promotes traits of mammary luminal progenitor cells[J].Oncotarget，2017，8(48)：83384-83406.

[10] VAN ITERSON M，BERVOETS S，DE MEIJER E J，et al.Integrated analysis of microRNA and mRNA expression：adding biological significance to microRNA target predictions[J].Nucleic Acids Res，2013，41(15)：e146.

[11] GOEMAN J J，MANSMANN U. Multiple testing on the directed acyclic graph of gene ontology[J]. Bioinformatics，2008，24(4)：537-544.

[12] LI X，GILL R，COOPER N G，et al. Modeling microRNA-mRNA interactions using PLS regression in human colon cancer[J]. BMC Med Genomics，2011，4：44.

[13] LUO W，FRIEDMAN M S，SHEDDEN K，et al. GAGE：generally applicable gene set enrichment for pathway analysis[J].BMC Bioinformatics，2009，10：161.

[14] TOMCZAK K，CZERWINSKA P，WIZNEROWICZ M. The Cancer Genome Atlas (TCGA)：an immeasurable source of knowledge[J]. Contemp Oncol (Pozn)，2015，19(1A)：68-77.

[15] LEE H，PALM J，GRIMES S M，et al. The cancer genome atlas clinical explorer：a web and mobile interface for identifying clinical-genomic driver associations[J]. Genome Med，2015，7：112.

[16] VASUDEVAN S，TONG Y，STEITZ J A. Switching from repression t o a ctivation：microRNAs c an u p-regulate t ranslation[J].Science，2007，318(5858)：1931-1934.

[17] ZHANG R，SU B. Small but influential：the role of microRNAs on gene regulatory network and 3'UTR evolution[J]. J Genet Genomics，2009，36(1)：1-6.

[18] OROM U A，NIELSEN F C，LUND A H. MicroRNA-10a binds the 5'UTR of ribosomal protein mRNAs and enhances their translation[J]. Mol Cell，2008，30(4)：460-471.

[19] KLOPOCKI E，KRISTIANSEN G，WILD P J，et al. Loss of SFRP1 is associated with breast cancer progression and poor prognosis in early stage tumors[J]. Int J Oncol，2004，25(3)：641-649.

[20] VEECK J， NIEDERACHER D， AN H， et al. Aberrant methylation of the Wnt antagonist SFRP1 in breast cancer is associated with unfavourable prognosis[J]. Oncogene，2006，25(24)：3479-3488.

[21] WENG C， NGUYEN T， SHIVELY J E. miRNA-342 Regulates CEACAM1-induced lumen formation in a threedimensional model of mammary gland morphogenesis[J]. J Biol Chem，2016，291(32)：16777-16786.

[22] GAO F，ZHANG Y F，ZHANG Z P，et al. miR-342-5p regulates neural stem cell proliferation and differentiation downstream to Notch signaling in mice[J]. Stem Cell Reports，2017，8(4)：1032-1045.

[23] LUGLI G，COHEN A M，BENNETT D A，et al. Plasma exosomal miRNAs in persons with and without Alzheimer disease：Altered expression and prospects for biomarkers[J].PLoS One，2015，10(10)：e0139233.

[24] SUN X，WU Y，GU M，et al. miR-342-5p decreases ankyrin G levels in Alzheimer's disease transgenic mouse models[J]. Cell Rep，2014，6(2)：264-270.

[25] TSUKITA S， FURUSE M. Pores in the wall： claudins constitute tight junction strands containing aqueous pores[J]. J Cell Biol，2000，149(1)：13-16.

[26] SINGH A B，SHARMA A，DHAWAN P. Claudin family of proteins and cancer：an overview[J]. J Oncol，2010，2010：541957.

[27] TOKES A M，SZASZ A M，JUHASZ E，et al. Expression of tight junction molecules in breast carcinomas analysed by array PCR and immunohistochemistry[J]. Pathol Oncol Res，2012，18(3)：593-606.