何文龍 隋杏玲

生物芯片高通量檢測卵巢癌轉(zhuǎn)錄組差異表達(dá)基因的生物標(biāo)記識別

何文龍 隋杏玲

目的探究卵巢癌中基因的表達(dá)情況并識別重要的代謝通路。方法從ArrayExpress數(shù)據(jù)庫下載卵巢癌芯片數(shù)據(jù)，使用rank product方法識別差異表達(dá)基因；利用STRING900、STRING、Bossi&Lehner和PSICQUIC四個數(shù)據(jù)庫對這些差異基因進(jìn)行KEGG通路富集分析，并運(yùn)用R語言運(yùn)算取這些通路的交集，以找出在卵巢癌發(fā)生中重要的代謝通路；最后構(gòu)建各通路中重疊基因的網(wǎng)絡(luò)節(jié)點(diǎn)圖，以識別通路中的重要基因。結(jié)果在卵巢癌中共檢測到1,732個差異表達(dá)基因和176條代謝通路。其中，p53信號通路、酪氨酸代謝、鞘脂類代謝是卵巢癌中最重要的代謝通路。這些通路中基因的表達(dá)均發(fā)生了顯著變化。結(jié)論所識別的重要代謝通路和這些通路中基因表達(dá)的變化可作為卵巢癌診斷的潛在生物標(biāo)記。

卵巢癌；差異表達(dá)基因；通路富集分析；網(wǎng)絡(luò)構(gòu)建

作者單位：300010 天津，天津醫(yī)科大學(xué)公共衛(wèi)生學(xué)院（何文龍）；250101 濟(jì)南，濟(jì)南基因云館數(shù)字圖書館（隋杏玲）（通訊作者：隋杏玲）

卵巢癌惡性腫瘤是危害女性健康與生命的惡性腫瘤之一，僅有20%左右的卵巢癌患者可在早期被發(fā)現(xiàn)，晚期患者的死亡率居婦科惡性腫瘤死亡率之首[1,2]。

目前，卵巢癌的早期診斷多采用血清CA125、癌胚抗原（carcino-embryonic antigen,CEA）、人附睪蛋白4（human epididymis protein-4,HE4）等作為重要標(biāo)記物。但這些標(biāo)記物普遍存在敏感性和特異性差、準(zhǔn)確率較低等問題[3]，遠(yuǎn)不能滿足臨床實際的應(yīng)用需要。代謝組學(xué)作為一種將圖像識別方法和生物信息學(xué)結(jié)合起來的分析技術(shù)，在卵巢癌的早期診斷、預(yù)后和治療評價方面得到了廣泛的應(yīng)用并已檢測出了可作為卵巢癌診斷的候選標(biāo)記物[4,5]。雖然這些研究為卵巢癌的生物學(xué)特性的異常改變及早期診斷提供了途徑，但對于卵巢癌發(fā)病的分子機(jī)理及各基因間的相互作用還不能清晰闡釋。已有的研究表明，卵巢癌的發(fā)生和發(fā)展是一個復(fù)雜的過程，具有復(fù)雜的生物學(xué)行為，涉及了多個基因的異常表達(dá)和相互作用。因此，尋找更加靈敏、特異、檢測方便的卵巢癌腫瘤標(biāo)記物，對于早期發(fā)現(xiàn)并提高治療效果具有重要意義。

高通量技術(shù)和生物信息學(xué)分析手段的聯(lián)合應(yīng)用，為從基因組水平研究疾病的分子基礎(chǔ)提供了途徑。本文運(yùn)用生物信息學(xué)方法，對從ArrayExpress數(shù)據(jù)庫選取的卵巢癌的微陣列芯片基因表達(dá)數(shù)據(jù)進(jìn)行挖掘，旨在通過對差異表達(dá)基因進(jìn)行通路富集分析，找出卵巢癌中重要的代謝通路，并對其中包含的基因進(jìn)行分析，以期為卵巢癌的診斷提供分子標(biāo)記。

1 材料與方法

1.1 數(shù)據(jù)收集與預(yù)處理 實驗所用數(shù)據(jù)芯片E-GEOD-10971[6]、E-GEOD-14001[7]、E-GEOD-18520[8]和E-GEOD-27651均選自ArrayExpress數(shù)據(jù)庫（http://www.ebi.ac.uk/arrayexpress/）。使用Affymetrix包讀取每個芯片數(shù)據(jù)，分別采用穩(wěn)健多陣列平均法（RMA算法）[9]和分位數(shù)法[10]對其進(jìn)行背景矯正與規(guī)范化，以消除非特異性雜交的影響；運(yùn)用mas5算法[11]對探針的匹配值與誤配值進(jìn)行修正，對于表達(dá)式的值則采用中位數(shù)方法進(jìn)行取值。使用GeneFilter包的featureFilter方法進(jìn)行數(shù)據(jù)篩選，去掉四分位數(shù)間距大于0.5的基因，若有多個探針對應(yīng)同一個基因，則保留其中間值。

1.2 差異表達(dá)基因的篩選 使用RankProd包中的rank product方法對4個芯片數(shù)據(jù)進(jìn)行合并，然后計算各基因的FoldChange（FC）值。|log2FC|＞2且P小于0.01的基因則視為差異表達(dá)基因，并區(qū)分出上調(diào)與下調(diào)的基因。

1.3 通路富集分析 利用在線分析軟件EnrichNet-Network-based enrichment analysis（http://www.enrichnet.org/）中的STRING900、STRING、Bossi&Lehner和PSICQUIC分子網(wǎng)絡(luò)數(shù)據(jù)庫對所得的|log2FC|≥3的差異表達(dá)基因進(jìn)行通路富集分析（Kyoto Encyclopedia of Genes and Genomes,KEGG），然后運(yùn)用R語言運(yùn)算取所有通路的交集，并計算每條通路的Xd值[12]的平均值，以找出在卵巢癌發(fā)生中比較重要的通路。

1.4 通路網(wǎng)絡(luò)構(gòu)建 從http://www.enrichnet.org/下載Xd值的平均值大于0.68的各通路的cytoscape-network數(shù)據(jù)，將其導(dǎo)入Cytoscape 3.1.0軟件中構(gòu)建它們的網(wǎng)絡(luò)節(jié)點(diǎn)圖，僅保留該通路基因集中的重疊基因，然后將四個數(shù)據(jù)庫中同一通路的網(wǎng)絡(luò)節(jié)點(diǎn)圖進(jìn)行合并，找出該通路的重要基因。

2 結(jié)果

2.1 差異表達(dá)基因的識別 通過對所得基因芯片數(shù)據(jù)進(jìn)行過濾篩選，共得到1,732個差異表達(dá)基因，其中上調(diào)基因580個，下調(diào)基因1,152個。前27個上調(diào)和前29個下調(diào)的差異表達(dá)基因的具體信息見表1。

2.2 通路富集分析 將|log2FC|≥3的656個差異表達(dá)基因輸入STRING900、STRING、Bossi&Lehner和PSICQUIC四個數(shù)據(jù)庫進(jìn)行在線KEGG生物通路注釋，分別得到了193條、196條、177條和192條通路。進(jìn)一步通過R語言運(yùn)算取交集后，共得到176條通路。比較重要的通路有p53信號通路、酪氨酸代謝、鞘脂類代謝、甘氨酸-絲氨酸-蘇氨酸代謝、細(xì)胞外基質(zhì)相互作用等。前10條通路的具體信息見表2。

2.3 基因的相互作用網(wǎng)絡(luò)分析 Cytoscape作為一個強(qiáng)大的分子相互作用網(wǎng)絡(luò)分析和可視化軟件，對生物信息學(xué)的研究具有重要幫助。連通度作為評價節(jié)點(diǎn)中心性的參數(shù)，其大小可以反映該節(jié)點(diǎn)在網(wǎng)絡(luò)中的重要性。通過Cytoscape軟件對卵巢癌中重要通路的差異表達(dá)基因進(jìn)行網(wǎng)絡(luò)圖譜分析可以發(fā)現(xiàn)，細(xì)胞周期蛋白依賴性激酶抑制劑2A（CDKN2A）在p53信號通路中與細(xì)胞周期蛋白（CCND2、CCNE1）和蛋白磷酸酶（PPM1D）密切相關(guān)，暗示該基因在p53信號通路改變中具有重要作用（圖1A）；單胺氧化酶（MAOA、MAOB）和乙醇脫氫酶（ADH1B）在酪氨酸代謝通路中起著中心性作用（圖1B）；在鞘脂類代謝通路中主要包含5個基因，各基因間緊密作用，共同調(diào)控疾病的發(fā)生（圖1C）。各通路中基因間的相互聯(lián)系及其具體參數(shù)見圖1和表3。

3 討論

作為女性腫瘤的高發(fā)疾病，卵巢癌的發(fā)生是多種體內(nèi)多種因素共同作用的結(jié)果。我們通過對卵巢癌基因表達(dá)芯片進(jìn)行分析，發(fā)現(xiàn)共有1,732個基因的表達(dá)發(fā)生了變化；進(jìn)一步利用不同數(shù)據(jù)庫對部分差異表達(dá)基因進(jìn)行通路富集分析表明，p53信號通路、酪氨酸代謝通路和鞘脂類代謝通路是卵巢癌發(fā)生中的重要通路，暗示這些通路的改變可能是卵巢癌發(fā)生、發(fā)展的重要原因之一。這與Fu等[13]的研究結(jié)果基本一致。

表1 卵巢癌中P值小于0.01的前27個上調(diào)和前29個下調(diào)的差異表達(dá)基因

表2 基于KEGG的通路富集分析

已知p53通路的異常改變是腫瘤細(xì)胞惡變轉(zhuǎn)化的主要原因之一。本研究通過對參與p53信號通路的基因的網(wǎng)絡(luò)分析發(fā)現(xiàn)，CDKN2A在p53信號通路起著重要作用。作為與細(xì)胞周期相關(guān)的抑癌基因，CDKN2A/p16在人類癌癥的發(fā)生中經(jīng)常發(fā)生突變[14]，但在早期卵巢癌和轉(zhuǎn)移性卵巢癌中發(fā)生突變的頻率卻很低[15]。Bian等[16]利用DNA重組技術(shù)構(gòu)建了人類卵巢癌細(xì)胞系UACC-1598的CDKN2A/p16-A148T突變體和CDKN2A/p16-野生型基因表達(dá)體系，發(fā)現(xiàn)在卵巢癌中CDKN2A/p16的第148位密碼子發(fā)生了單堿基突變。Kim等[17]在早期卵巢癌和復(fù)發(fā)性上皮性卵巢癌也發(fā)現(xiàn)了該基因的突變。在本研究中，CDKN2A（pfp=3.671164，P=0）是差異表達(dá)基因之一，其表達(dá)量上調(diào)可能是導(dǎo)致p53通路改變的重要原因。

表3 各通路中差異表達(dá)基因的信息

酪氨酸作為苯丙氨酸代謝的產(chǎn)物之一，其代謝的異?？梢鸺膊〉陌l(fā)生。Fong等[4]通過對卵巢癌患者的卵巢代謝物進(jìn)行分析發(fā)現(xiàn)，苯丙酮酸和苯乳酸的含量明顯升高；卵巢癌患者在接受治療后，血清中酪氨酸-賴氨酸-亮氨酸-40的含量明顯下降[18]。卵巢癌中酪氨酸代謝的異?？赡芘c苯丙氨酸代謝密切相關(guān)。我們的研究結(jié)果表明，酪氨酸代謝是卵巢癌中重要的代謝通路之一，與酪氨酸下游產(chǎn)物代謝相關(guān)的基因發(fā)生了明顯下調(diào)。研究表明，卵巢癌治療藥物藥效的發(fā)揮與酪氨酸激酶受體的下調(diào)有關(guān)[19]。目前已有利用酪氨酸激酶定向治療卵巢癌的報道，盡管它的療效還無定論[20]。這些結(jié)果暗示，酪氨酸代謝可作為卵巢癌診斷的指標(biāo)之一。

鞘脂類作為細(xì)胞膜的組成成分之一，其代謝產(chǎn)物神經(jīng)酰胺、神經(jīng)鞘氨醇和1-磷酸鞘氨醇在腫瘤的發(fā)生發(fā)展中發(fā)揮著重要作用，它們可以調(diào)節(jié)細(xì)胞的增值、存活及凋亡等。Babahosseini等[21]的研究發(fā)現(xiàn)，神經(jīng)酰胺和1-磷酸鞘氨醇能夠降低發(fā)生惡變的小鼠卵巢上皮細(xì)胞的彈性，而鞘氨醇的作用卻正好相反。Guillwrmet-Guibert等[22]發(fā)現(xiàn)，通過加入鞘氨醇激酶-1抑制劑、神經(jīng)酰胺類似物等提高細(xì)胞內(nèi)神經(jīng)酰胺的水平或者降低細(xì)胞內(nèi)鞘氨醇激酶-1的活性來提高神經(jīng)酰胺、1-磷酸鞘氨醇的比率，可有效提高胰腺癌細(xì)胞對胞苷的化學(xué)敏感性，從而殺死癌細(xì)胞。Anderson等[23]研究表明，通過添加外源性鞘氨醇能夠使小鼠上皮性卵巢癌細(xì)胞中檸檬酸合成酶活性降低、促進(jìn)三羧酸循環(huán)，減少膽固醇的合成和糖酵解，從而改變細(xì)胞代謝路徑。以上研究表明，鞘脂類代謝與腫瘤的發(fā)生、進(jìn)展密切相關(guān)，它在腫瘤的治療中具有巨大的潛力。Choi等[24]利用小鼠動物模型研究發(fā)現(xiàn)α-半乳糖神經(jīng)酰胺輔助腫瘤細(xì)胞疫苗治療卵巢癌，可打破腫瘤細(xì)胞的免疫耐性，明顯提高血清中γ-干擾素的水平。我們的通路富集分析結(jié)果表明，鞘脂類代謝在卵巢癌的發(fā)生發(fā)展中占有重要地位，與鞘脂類代謝相關(guān)的主要基因密切相互作用并都下調(diào)表達(dá)，這些基因的表達(dá)的異常可能是導(dǎo)致卵巢癌發(fā)生的重要原因。

綜上所述，本研究通過生物信息學(xué)方法對卵巢癌相關(guān)基因的表達(dá)和通路富集分析發(fā)現(xiàn)，在卵巢癌的發(fā)生發(fā)展中，有大量基因的表達(dá)發(fā)生了明顯變化（上調(diào)或下調(diào)），并引起了通路代謝的異常，p53信號通路、酪氨酸代謝通路和鞘脂類代謝通路可作為卵巢癌診斷的重要指標(biāo)。

1 Teneriello MG,Park RC.Early detection of ovarian cancer.CA Cancer J Clin,1995,45(2): 71-87.

2 Jacobs IJ,Menon U.Progress and challenges in screening for early detection of ovarian cancer.Mol Cell Proteomics,2004,3(4): 355-366.

3 Williams TI,Toups KL,Saggese DA,et al.Epithelial ovarian cancer: Disease etiology,treatment,detection,and investigational gene,metabolite,and protein biomarkers.J Proteome Res,2007,6(8): 2936-2962.

4 Fong MY,McDunn J,Kakar SS.Identification of metabolites in the normal ovary and their transformation in primary and metastatic ovarian cancer.PLoS One,2011,6(5): e19963.

5 Fan L,Zhang W,Yin M,et al.Identification of metabolic biomarkers to diagnose epithelial ovarian cancer using a UPLC/QTOF/MS platform.Acta Oncol,2012,51(4):473-479.

6 Tone AA,Begley H,Sharma M,et al.Gene expression profiles of luteal phase fallopian tube epithelium fromBRCAmutation carriers resemble high-grade serous carcinoma.Clin Cancer Res,2008,14(13): 4067-4078.

7 Tung CS,Mok SC,Tsang YT,et al.PAX2 expression in low malignant potential ovarian tumors and low-grade ovarian serous carcinomas.Mod Pathol,2009,22(9):1243-1250.

8 Gamwell LF,Gambaro K,Merziotis M,et al.Small cell ovarian carcinoma: genomic stability and responsiveness to therapeutics.Orphanet J Rare Dis,2013,8: 33.

9 Ma L,Robinson LN,Towle HC.ChREBP.ChREBP＊Mlx is the principal mediator of glucose-induced gene expression in the liver.J Biol Chem,2006,281(39):28721-28730.

10 Rifai N,Ridker PM.Proposed cardiovascular risk assessment algorithm using high-sensitivity C-reactive protein and lipid screening.Clin Chem,2001,47(1): 28-30.

11 Zhang L,Miles MF,Aldape KD.A model of molecular interactions on short oligonucleotide miroarrays.Nat Biotechnol,2003,21(7): 818-821.

12 Glaab E,Baudot A,Krasnogor N,et al.EnrichNet:network-based gene set enrichment analysis.Bioinformatics,2012,28(18): i451-i457.

13 Fu LJ,Wang B.investigation of the hub genes and related mechanism in ovaian cancer via bioinformatics analysis.J Ovarian Res,2013,6(1): 92.

14 Witkiewicz AK,Knudsen KE,Dicker AP,et al.The meaning of p16 (ink4a) expression in tumors:functional significance,clinical associations and future developments.Cell Cycle,2011,10(15): 2497-2503.

15 Schuyer M,van Staveren IL,Klijn JG,et al.SporadicCDKN2(MTS1/p16ink4) gene alterations in human ovarian tumours.Br J Cancer,1996,74(7): 1069-1073.

16 Bian Z,Yu Y,Yang T,et al.Effect of tumor suppressor gene cyclin-dependent kinase inhibitor 2A wild-type and A148T mutant on the cell cycle of human ovarian cancer cells.Oncol Lett,2014,7(4): 1229-1232.

17 Kim YM,Lee SW,Chun SM,et al.Analysis and comparison of somatic mutations in paired primary and recurrent epithelial ovarian cancer samples.PLoS One,2014,9(6): e99451.

18 Choudhuri S,Sharma C,Banerjee A,et al.A repertoire of biomarkers helps in detection and assessment of therapeutic response in epithelial ovarian cancer.Mol Cell Biochem,2014,386(1-2): 259-269.

19 Cho YR,Choi SW,Seo DW.Thein vitroantitumor activity of Siegesbeckia glabrescens against ovarian cancer through suppression of receptor tyrosine kinase expression and the signaling pathways.Oncol Rep,2013,30(1): 221-226.

20 Morotti M,Becker CM,Menada MV,et al.Targeting tyrosine-kinases in ovarian cancer.Expert Opin Investig Drugs,2013,22(10): 1265-1279.

21 Babahosseini H,Roberts PC,Schmelz EM,et al.Bioactive sphingolipid metabolites modulate ovarian cancer cell structural mechanics.Integr Biol (Camb),2013,5(11):1385-1392.

22 Guillermet-Guibert J,Davenne L,Pchejetski D,et al.Targeting the sphingolipid metabolism to defeat pancreatic cancer cell resistance to the chemotherapeutic gemcitabine drug.Mol Cancer Ther,2009,8(4): 809-820.

23 Anderson AS,Roberts PC,Frisard MI,et al.Metabolic changes during ovarian cancer progression as targets for sphingosine treatment.Exp Cell Res,2013,319(10):1431-1442.

24 Choi YS,Hoory T,Monie A,et al.Alpha-Galactosylceramide enhances the protective and therapeutic effects of tumor cell based vaccines for ovarian tumors.Vaccine,2008,26(46): 5855-5863.

Microarray data analysis and biomarker identification of ovarian cancer

Wenlong HE1,Xingling SUI2

1School of Public Health,Tianjin Medical University,Tianjin 300010,China;2Digital Library,Ji'nan Gene Cloud Museum,Ji’nan 250101,China

ObjectiveTo explore the gene expression pattern and identify the important metabolic pathways in ovarian cancer.MethodsThemicroarray data of ovarian cancer was downloaded from ArrayExpress database,and the differentially expressed (DE) genes in it were identified by the rank product method.The KEGG pathway enrichment analysis of these DE genes was performed by using the on line database STRING900,STRING,Bossi&Lehner and PSICQUIC.The R language was used to computing the intersection of these pathways and finding out the important pathways in ovarian cancer.Finally,the important genes and their relations in these pathways were analyzed by constructing the nodes network of the overlapped genes.ResultsA total of 1,732 DE genes and 176 metabolic pathways were involved in the occurrence.The p53 signaling pathway,tyrosine metabolism and sphingolipid metabolism were the most important pathways in ovarian cancer.The expression pattern of the genes that involved in these pathways was also greatly changed.ConclusionThe important pathways and change of the genes involved in them can be regard as the underlying biomarkers for the diagnosis of ovarian cancer.

Ovarian cancer; Differentially expressed genes; Pathway enrichment analysis; Network construction