卵巢癌發(fā)生發(fā)展關(guān)鍵基因的生物信息學(xué)篩選

余展鵬， 彭 亮， 姜 巧

卵巢癌是婦科最為常見的惡性腫瘤之一，其發(fā)病率逐年上升，病死率位居女性惡性腫瘤的第8位[1]。卵巢癌起病隱匿，大部分患者確診時(shí)已為晚期，預(yù)后不佳，5年生存率僅為30%～40%[2,3]。為此，尋找與其發(fā)生、發(fā)展相關(guān)聯(lián)的基因，對其診斷、治療和預(yù)后評估具有重要意義。本研究擬使用生物信息學(xué)方法，整合分析高通量基因表達(dá)(Gene Expression Omnibus，GEO)數(shù)據(jù)庫和癌癥基因圖譜(The Cancer Genome Atlas，TCGA)數(shù)據(jù)庫及KMplot(Kaplan-Meier Plotter)等多個(gè)數(shù)據(jù)庫中與卵巢癌的相關(guān)數(shù)據(jù)集，以篩選卵巢癌發(fā)生發(fā)展的關(guān)鍵基因。現(xiàn)報(bào)告如下。

1 資料與方法

1.1數(shù)據(jù)來源與差異基因(differentially expressed ge-nes,DEGs)篩選 (1)在GEO數(shù)據(jù)庫(https://www.ncbi.nlm.nih.gov/geo/)中，按以下設(shè)置檢索：(“ovarian neoplasms”[MeSH terms] or ovarian cancer[all fields]) and “homo sapiens”[porgn] and (“gse”[filter] and “expression profiling by array”[filter])。根據(jù)檢索結(jié)果下載得到GSE54388、GSE105437、GSE40595數(shù)據(jù)集，3份數(shù)據(jù)均由GPL570平臺完成、上傳。它們分別包含6例、5例及8例正常卵巢上皮樣本(正常組)和16例、10例及31例卵巢癌樣本(卵巢癌組)。使用R 3.5.2軟件進(jìn)行數(shù)據(jù)預(yù)處理：以log FC>1.5或log FC<-1.5且校正P<0.05為標(biāo)準(zhǔn)，篩選正常組與腫瘤組的DEGs，并選出其中的共有部分，作為本研究的候選基因。(2)從Xena數(shù)據(jù)庫[4](https://xena.ucsc.edu/)下載TCGA與GTEx(Genotype-Tissue Expression)的合并數(shù)據(jù)集(卵巢癌組=426，正常組=88)，分析兩組候選基因的表達(dá)情況，進(jìn)一步確定正常卵巢組織與卵巢癌組織的DEGs。

1.2DEGs與患者預(yù)后關(guān)系分析 使用KMplot數(shù)據(jù)庫(https://kmplot.com/analysis/)卵巢癌數(shù)據(jù)集分析DEGs對患者預(yù)后的影響，以其自動計(jì)算各DEGs的cut-off值將患者分為高表達(dá)組和低表達(dá)組，并通過數(shù)據(jù)庫缺省設(shè)置分別繪制各組Kaplan-Meier生存曲線，分析DEGs與患者生存預(yù)后的關(guān)系，篩選卵巢癌關(guān)鍵基因。

1.3卵巢癌關(guān)鍵基因的表達(dá)數(shù)據(jù)分析 從TCGA數(shù)據(jù)庫(https://portal.gdc.cancer.gov/)中檢索“Ovarian cancer”，下載380例卵巢癌組織的RNA-Seq數(shù)據(jù)集及與之對應(yīng)的臨床信息。剔除石蠟包埋標(biāo)本，整理數(shù)據(jù)并進(jìn)行標(biāo)準(zhǔn)化處理，分析卵巢癌關(guān)鍵基因表達(dá)水平與患者臨床特征的關(guān)系。

1.4卵巢癌關(guān)鍵基因調(diào)控機(jī)制的初步分析

1.4.1 基因集富集分析 選取TCGA下載的卵巢癌組織RNA-Seq數(shù)據(jù)，運(yùn)行GSEA 3.0(Gene Set Enrichment Analysis 3.0)軟件，以兩個(gè)卵巢癌關(guān)鍵基因[ARX(Aristaless Related Homeobox)，F(xiàn)AM150B(Family with Sequence Similarity 150 Member B)]的中位表達(dá)值為界，分為高表達(dá)組(190例)和低表達(dá)組(189例)，設(shè)置Hallmark Gene Sets為參照基因集，進(jìn)行KEGG(Kyoto Encyclopedia of Genes and Genomes)信號通路聚類分析。

1.4.2 卵巢癌關(guān)鍵基因與免疫相關(guān)基因的關(guān)聯(lián)分析 使用在線工具TIMER(https://cistrome.shinyapps.io/timer/)，選取“Ovarian Serous Cystadenocarcinoma”數(shù)據(jù)集，設(shè)置“tumor purity”為校正變量，分別使用GSEA富集得到的免疫系統(tǒng)相關(guān)信號通路的關(guān)鍵分子與ARX和FAM150B基因的表達(dá)水平進(jìn)行Spearman相關(guān)分析。

1.4.3 卵巢癌關(guān)鍵基因的相互作用網(wǎng)絡(luò)分析 使用STRING數(shù)據(jù)庫(https://string-db.org/)分別檢索“ARX”和“FAM150B”及與之相關(guān)的蛋白(即1.4.2步驟獲得的信號通路的關(guān)鍵分子)，構(gòu)建ARX和FAM150B的蛋白相互作用網(wǎng)絡(luò)。

2 結(jié)果

2.1卵巢癌DEGs的篩選結(jié)果 從GEO數(shù)據(jù)庫下載的GSE54388、GSE105437、GSE40595數(shù)據(jù)集中分別篩選到450、925、1097個(gè)DEGs。根據(jù)其上調(diào)、下調(diào)情況進(jìn)行韋恩圖分析，得到10個(gè)下調(diào)的候選基因及4個(gè)上調(diào)的候選基因。見表1。使用TCGA和GTEx數(shù)據(jù)集中的卵巢癌組織及正常卵巢組織的RNA-Seq數(shù)據(jù)對DEGs進(jìn)行驗(yàn)證，結(jié)果顯示有13個(gè)基因?yàn)槁殉舶┑暮蜻x關(guān)鍵基因。見圖1。

2.2卵巢癌DEGs與患者生存預(yù)后關(guān)聯(lián)性分析結(jié)果 使用KMplot數(shù)據(jù)庫分析以上13個(gè)DEGs與卵巢癌患者總生存期的關(guān)系，按檢驗(yàn)水準(zhǔn)α=0.01，結(jié)果顯示，在下調(diào)的9個(gè)DEGs中，ADH1B低表達(dá)組的生存預(yù)后顯著優(yōu)于高表達(dá)組(P<0.01)；而ARX、FAM150B低表達(dá)組預(yù)后顯著差于高表達(dá)組(P<0.01)。在上調(diào)的4個(gè)DEGs中，COL4A1高表達(dá)組預(yù)后顯著優(yōu)于低表達(dá)組(P<0.01)。見圖2。

表1 卵巢癌DEGs基因

圖1 卵巢癌DEGs驗(yàn)證結(jié)果圖

圖2 卵巢癌DEGs與患者生存預(yù)后的生存曲線圖

2.3ARX和FAM150B表達(dá)水平與卵巢癌患者臨床病理特征的關(guān)聯(lián)性分析結(jié)果 結(jié)合2.2的結(jié)果，發(fā)現(xiàn)ARX和FAM150B表達(dá)量的下調(diào)，可能會導(dǎo)致患者的預(yù)后不良，鑒此，本研究選取ARX和FAM150B作進(jìn)一步研究。結(jié)果顯示，ARX的表達(dá)量與卵巢癌患者腫瘤直徑具有關(guān)聯(lián)性(P<0.05)，腫瘤直徑越大，ARX表達(dá)量越低。FAM150B的表達(dá)量與患者的年齡及腫瘤直徑具有關(guān)聯(lián)性(P<0.05)，年齡或腫瘤直徑越大，F(xiàn)AM150B表達(dá)量越低。見圖3。

圖3 ARX和FAM150B表達(dá)水平與卵巢癌患者臨床病理特征的關(guān)聯(lián)性分析結(jié)果圖

2.4ARX及FAM150B調(diào)控機(jī)制的初步分析結(jié)果 使用TCGA卵巢癌RNA-Seq數(shù)據(jù)進(jìn)行GSEA PATHWAY富集分析，結(jié)果顯示，ARX低表達(dá)組較高表達(dá)組富集的通路有24條；FAM150B低表達(dá)組較高表達(dá)組富集的通路有14條。進(jìn)一步分析發(fā)現(xiàn)，ARX與FAM150B低表達(dá)組富集到的通路中分別有12條和8條與免疫系統(tǒng)相關(guān)。整理以上免疫系統(tǒng)相關(guān)信號通路的關(guān)鍵分子，并使用TIMER數(shù)據(jù)庫分析其與ARX和FAM150B表達(dá)量的相關(guān)關(guān)系。Spearman相關(guān)分析結(jié)果顯示，ARX與PIK3R2(Phosphoinositide-3-Kinase Regulatory Subunit 2)的表達(dá)量呈正相關(guān)(P<0.05)，與IL-6及PIK3CD的表達(dá)量呈負(fù)相關(guān)(P<0.05)。見表2。FAM150B與PIK3R1的表達(dá)量呈正相關(guān)，與VAV1、PIK3CA、RAC2、MAPK1的表達(dá)量呈負(fù)相關(guān)(P<0.05)。見表3。使用STRING數(shù)據(jù)庫分別構(gòu)建ARX與IL-6、PIK3R2、PIK3CD，以及FAM150B與VAV1、PIK3CA、RAC2、PIK3R1、MAPK1的互作網(wǎng)絡(luò)，預(yù)測結(jié)果顯示，ARX可能可以通過EGFR對IL-6、PIK3R2及PIK3CD進(jìn)行調(diào)控，F(xiàn)AM150B可能可以通過PLXNA2對VAV1、PIK3CA、RAC2、PIK3R1、MAPK1進(jìn)行調(diào)控。見圖4。

表2 ARX與免疫相關(guān)基因的Spearman相關(guān)分析結(jié)果

表3 FAM150B與免疫相關(guān)基因的Spearman相關(guān)分析結(jié)果

圖4 ARX和FAM150B與免疫相關(guān)基因的互作網(wǎng)絡(luò)圖

3 討論

3.1本研究通過聯(lián)合GEO、TCGA及GTEx數(shù)據(jù)庫分析發(fā)現(xiàn)ARX與FAM150B基因在卵巢癌中顯著下調(diào)，且它們的表達(dá)量與腫瘤體積及患者預(yù)后密切相關(guān)，表達(dá)量越低，腫瘤體積越大，患者預(yù)后越差(P<0.01)。ARX基因定位于Xp22.1-Xp21.3之間，其包含5個(gè)外顯子，編碼含562個(gè)氨基酸的ARX蛋白，該蛋白為轉(zhuǎn)錄調(diào)控因子。既往研究[5]表明其在神經(jīng)系統(tǒng)的發(fā)育中起重要作用，近來也有研究顯示ARX基因與胰腺神經(jīng)內(nèi)分泌瘤、前列腺癌、卵巢癌等疾病的發(fā)生發(fā)展密切相關(guān)。FAM150B基因定位于2p25.3，其包含11個(gè)外顯子，編碼含152個(gè)氨基酸的FAM150B蛋白，該蛋白為ALK(anaplastic lymphoma receptor tyrosine kinase)和LTK(leukocyte receptor tyrosine kinase)的配體。既往研究[6]表明其在細(xì)胞的增殖、分化中起信號轉(zhuǎn)導(dǎo)作用，但目前國內(nèi)外與之相關(guān)的研究較少，僅發(fā)現(xiàn)其與成神經(jīng)細(xì)胞瘤和鼻咽癌的發(fā)生、發(fā)展相關(guān)。

3.2本研究使用TCGA卵巢癌數(shù)據(jù)集中的數(shù)據(jù)，將獲取數(shù)據(jù)分為高表達(dá)組與低表達(dá)組進(jìn)行GSEA PATHWAY富集，結(jié)果發(fā)現(xiàn)兩個(gè)基因的低表達(dá)組均大比例地富集到與免疫系統(tǒng)相關(guān)的信號通路。因此，我們推測這兩個(gè)基因的下調(diào)可能會影響機(jī)體的免疫功能。鑒此，本研究進(jìn)一步使用TIMER數(shù)據(jù)庫分析ARX基因及FAM150B基因和與之相關(guān)免疫信號通路關(guān)鍵分子的表達(dá)情況。結(jié)果顯示，在ARX基因下調(diào)的狀態(tài)下，IL-6和PIK3CD基因的表達(dá)量上調(diào)，PIK3R2基因的表達(dá)量下調(diào)；而在FAM150B基因下調(diào)的狀態(tài)下，VAV1、RAC2、MAPK1及PIK3CA基因的表達(dá)量上調(diào)，PIK3R1基因的表達(dá)量下調(diào)。其中PIK3R1和PIK3R2基因編碼的蛋白為PI3K的調(diào)節(jié)亞基，PIK3CA和PIK3CD基因編碼的蛋白為PI3K的催化亞基[7]。調(diào)節(jié)亞基的下調(diào)或催化亞基的上調(diào)都會導(dǎo)致PI3K的活性增強(qiáng)，從而減弱抗腫瘤免疫效應(yīng)[8]。Vav1基因編碼的蛋白是一種鳥嘌呤核苷酸轉(zhuǎn)換因子，能促進(jìn)Rho家族的Rac2表達(dá)量增加[9]，而Rac2又可以進(jìn)一步上調(diào)MAPK1的表達(dá)[10]。MAPK1基因編碼蛋白ERK2作為MAPK信號通路上的重要分子，其上調(diào)可使MAPK信號通路激活，繼而上調(diào)PD-L1的表達(dá)量，使腫瘤細(xì)胞產(chǎn)生免疫耐受[11,12]。IL-6基因編碼的IL-6作為炎癥因子，在腫瘤微環(huán)境中不但能抑制樹突細(xì)胞的遷移能力，減少其與T細(xì)胞的有效接觸[13]，還能下調(diào)Th1細(xì)胞MHC Ⅱ類分子的表達(dá),并減少IFN-γ和IL-2的釋放，從而降低細(xì)胞毒性T細(xì)胞的活性[14],使腫瘤細(xì)胞能夠逃避抗腫瘤免疫反應(yīng)。

3.3為了探討ARX和FAM150B調(diào)控上述基因的機(jī)制，本研究使用STRING數(shù)據(jù)庫繪制了它們的互作網(wǎng)絡(luò)，結(jié)果顯示ARX可能通過TUBA1A及EGFR對IL-6、PIK3R2、PIK3CD進(jìn)行調(diào)控；FAM150B可能通過PLXNA2對VAV1、PIK3CA、RAC2、PIK3R1、MAPK1進(jìn)行調(diào)控。EGFR作為腫瘤發(fā)生、發(fā)展過程中的至關(guān)重要的基因，目前已被廣泛研究[15]，而對于TUBA1A和PLXNA2的研究則相對較少，但近來亦有研究[16,17]表明TUBA1A和PLXNA2也是腫瘤發(fā)生發(fā)展中的關(guān)鍵調(diào)控基因。TUBA1A基因編碼的α-微管蛋白是微管蛋白家族中重要的一員，是參與細(xì)胞分裂和運(yùn)動所必需的蛋白。Zhu等[18]研究發(fā)現(xiàn)，在宮頸癌中TUBA1A可作為miR-15a-5p和miR-15b-5p的靶點(diǎn)，參與腫瘤細(xì)胞增殖及侵襲過程的調(diào)控。另外，PLXNA2既往一直只被認(rèn)為在機(jī)體神經(jīng)系統(tǒng)的發(fā)育過程中發(fā)揮關(guān)鍵調(diào)控作用，但目前也有越來越多的研究強(qiáng)調(diào)了其在腫瘤細(xì)胞侵襲和遷移過程中的調(diào)控功能。Tian等[19]分析基因芯片發(fā)現(xiàn)，與原發(fā)或局限性的前列腺癌相比，在轉(zhuǎn)移性的前列腺癌組織中，PLXNA2基因顯著上調(diào)。該研究還通過體外實(shí)驗(yàn)證實(shí)了PLXNA2參與了前列腺癌細(xì)胞遷移和侵襲的調(diào)控。因此，通過對上述結(jié)果的分析，我們推測在卵巢癌中下調(diào)的ARX和FAM150B可能可以通過調(diào)控TUBA1A、EGFR、PLXNA2，進(jìn)而對IL-6、PIK3R2、PIK3CD、VAV1、PIK3CA、RAC2、PIK3R1及MAPK1進(jìn)行調(diào)控，從而使機(jī)體抗腫瘤免疫反應(yīng)減弱，最終導(dǎo)致腫瘤細(xì)胞出現(xiàn)免疫逃逸或者產(chǎn)生免疫耐受，促進(jìn)卵巢癌的發(fā)生發(fā)展。

綜上所述，本研究通過生物信息學(xué)分析發(fā)現(xiàn)ARX和FAM150B在卵巢癌組織中呈低表達(dá)狀態(tài)，且它們的表達(dá)量和患者的預(yù)后密切相關(guān)，經(jīng)過進(jìn)一步分析，我們還發(fā)現(xiàn)ARX和FAM150B能對一系列基因進(jìn)行調(diào)控，從而減弱機(jī)體的抗腫瘤免疫反應(yīng)。在此研究中，所有用于分析的基因芯片原始數(shù)據(jù)均來自于同一平臺，研究方法一致，數(shù)據(jù)客觀，因此認(rèn)為本研究能在一定程度上為ARX和FAM150B在卵巢癌中的進(jìn)一步研究提供理論基礎(chǔ)。