基于生物信息基因表达谱的乳腺癌转移基因标志物的分析研究
吴奇桥1 张树民2 郑婷婷3 刘 娟2 胡 永2 林登强4 戴梦婷1 孙 菁2▲
1.复旦大学附属中山医院厦门医院放疗科,福建厦门 361006;2.复旦大学附属中山医院放疗科,上海 200030;3.复旦大学附属中山医院厦门医院放射科,福建厦门 361006;4.复旦大学附属中山医院厦门医院泌尿外科,福建厦门 361006
[摘要] 目的确定潜在和乳腺癌转移相关的基因标志物,并对其进行生存分析。方法 使用基因表达综合数据库(GEO),筛选乳腺癌转移部位和非转移部位之间的差异表达基因(DEGs)。对DEGs 进行基因本体(GO)和京都基因与基因组百科全书途径富集(KEGG)分析,构建蛋白质-蛋白质交互(PPI)网络,并使用MCODE 进行模块分析。然后,通过Kaplan-Meier 工具对中枢基因进行总生存(OS)和远处无转移生存(DMFS)分析。通过GEPIA 对基因和临床生存数据的相关性进行验证。结果 本研究共获得295 个DEGs,主要与细胞外分泌、粘着相关通路、细胞分裂相关。KEGG 通路分析表明,重要的通路包括丝裂原活化蛋白激酶(MAPK)信号通路、Rap1 信号通路、细胞粘附分子(CAMs)、抗原加工和呈递。建立了具有280 个节点和357 个连接的蛋白-蛋白交互网络。从蛋白-蛋白交互网络中发现了10 个模块。结论 共筛选出20 个基因作为枢纽基因,其中TYMS、SKA1、ADCY7 的低表达和MX1 的高表达与OS 较差有关。POLR3H 的低表达和CDCA8、ASE1、KIF14、MX1 的高表达与较差的DMFS 有关。POLR3H和PRC1 可作为诊断乳腺癌转移或潜在药物靶标的新生物标志物。
[关键词]乳腺癌;转移;生物信息学分析;差异表达的基因;治疗药物;生存分析
乳腺癌是威胁女性健康的最常见恶性肿瘤,是美国第二大最常见的癌症相关死亡[1]。乳腺原位癌通常具有良好的预后,然而,乳腺癌若出现远处转移往往会导致危及生命的结果[2]。总体而言,乳腺癌的5年平均生存率为90%,但如果存在远处转移,则降至26%[3]
乳腺癌进展的分子机制尚未完全了解。鉴于其高死亡率,迫切需要弄清乳腺癌转移的潜在分子机制。既往的研究已经调查乳腺癌转移的相关基因,如Chen 等[4]揭示了ECM-受体相互作用可能有助于乳癌骨转移;Cai 等[5]表明CDCA8、CCNA2 与乳癌远处转移有关;Zheng 等[6]鉴定了几种与乳腺癌转移相关的基因。然而上述研究仅对单个数据集进行分析,目前仍然没有研究结合相似数据集进行基因分析。本研究分析了三个数据集中与乳腺癌转移相关的差异表达基因(differentially expressed genes,DEGs),目的是更好地了解潜在乳腺癌转移的机制,并找到潜在生物标志物和治疗靶标。
1 资料与方法
1.1 一般资料
GSE32489、GSE14776 和GSE103357[7]使 用 基 因表达综合数据库(Gene Expression Omnibus,GEO,http://www.ncbi.nlm.nih.gov/geo)得到的三个基因数据集,均使用Illumina HumanRef 平台芯片,根据平台中的注释信息将探针转换为相应的基因符号。GSE32489包含非转移尸检组织19 个,淋巴结转移组织90 个。GSE14776 包含8 个非转移细胞样本和6 个骨转移样本。GSE103357 包含2 个非转移细胞样本和3 个骨转移样本。
1.2 DEGs 的鉴定
Network Analyst(版本号:10.0)[8-9](http://www.networkanalyst.ca)用于提取数据集乳腺肿瘤样本和转移样本之间的DEGs。具有调整P 值<0.05 和|log2FC|>1.0的基因被认为是DEGs。通过维恩(Venn)工具(http://bioinformatics.psb.ugent.be/webtools/Venn/)检测重叠以得到共差异表达的基因(co different-expressed genes,Co-DEGs)。调整P 值、Benjamini 和Hochberg 假阳性发现率,用于在统计显著基因的发现和假阳性之间提供平衡。没有相应基因符号的探针组或具有多于一个探针组的基因被去除。识错率(false discovery rate,FDR)的计算公式如下:q-value(i)=p(i)length(p)/rank(p),表明当样本量增大时,检出假阳性率的概率增高,FDR 值越高,表明大样本数据富集后出现的假阳性概率越高。
1.3 基因本体和京都基因与基因组百科全书途径富集通路分析
基因本体(gene ontology,GO)分析用于对基因集进行功能研究[10]。京都基因与基因组百科全书途径富集(Kyoto encyclopedia of genes and genomes,KEGG)[11]是处理基因组和生物通路的数据库集合。注释可视化和 集 成 发 现 数 据 库(DAVID,http://david.ncifcrf.gov)(6.8 版)[12]用于解开已识别共基因的GO 和KEGG 途径。
1.4 DEGs 和模块选择的蛋白质-蛋白质交互网络分析
STRING(版本号:11.5)用于阐明蛋白质-蛋白质交互(protein-protein interaction,PPI)相互作用[13]。使用Cytoscape(版本号:3.8.2)[14]进行可视化PPI 网络。选择>0.4 的组合分数作为阈值。PPI 网络可以协助从蛋白层面识别参与乳腺癌转移的重要基因模块。此外,应用分子模块检测(molecular complex detection,MCODE)[15]插件对重要模块进行验证。
1.5 核心基因的热图绘制、GO 和临床分析
使用cytoHubba(版本号:3.8.2)插件和最大集团中心性(maximum group centrality,MCC)方法来识别前20 个中枢基因。核心基因的热图是通过使用R 中的热图包绘制的,使用的聚类方法为离差平和法(Ward),标准化方法为正态标准化(Z-score scaling)。使用Kaplan Meier-plotter(KM plotter,http://kmplot.com/analysis/)实现生存分析,Kaplan-Meier plotter 是一种生存分析软件[16]。通过Cytoscape 的BiNGO 插件(版本号:3.8.2)评估核心基因的GO 功能。
1.6 DEG 的生存分析
乳腺癌患者根据特定基因的表达分为两组(高表达与低表达)。通过Kaplan-Meier 生存图比较两个组患者的队列,并计算具有95%置信区间和对数等级P 值的风险比(HR),其中,HR=1 意味着基因低表达与高表达的等效性,若低表达处理优于高表达,则HR<1;若低表达处理劣于高表达,则HR>1。GEPIA[17]是基于TCGA 数据库的在线基因表达谱分析工具,用于验证枢纽基因与临床分期之间的相关性,其中,F value 代表单因素分析的F 值,这个值越大,表示组间差异越大,且当P<0.05 时,表示该基因在肿瘤不同分期中差异是有统计学意义的。
2 结果
2.1 三个基因数据集中DEG 和co-DEG 的鉴定
三个基因数据集中,韦恩图见图1A 及图1B,基因表达热图见图1C。结果显示,295 个基因被鉴定为co-DEGs。其中上调151 个,下调144 个。调整P<0.05和|log2FC|>1 被设置为截止标准。
图1 维恩图及基因热图
2.2 DEGs 的GO 富集分析及KEGG 通路分析
DAVID 在线工具用于阐明富集的co-DEGs 中的GO 和KEGG 通路。结果表明,对于BP,上调的基因主要富集在Ⅰ型干扰素信号通路、凋亡过程中,而下调的基因主要集中在DNA 修复、先天免疫反应等方面。KEGG 通路结果显示,DEGs 显着富集在与癌症发展的信号通路,包括丝裂原活化蛋白激酶(mitogen-activated protein kinase,MAPK) 信号通路和Rap1 信号通路(表1~2)。
表1 Co-DEGs 中上调基因的GO 分析
2.3 DEGs 和模块选择的PPI 网络分析
基于STRING 数据库得到PPI 网络分析结果,见图2A,结果显示,DEG 的PPI 网络由280 个节点和357 条边构成。使用MCC 方法,总共选择了20 个基因作为枢纽基因,分别为:CHEK1、POLR3H、IFITM1、XAF1、MCM5、ADCY3、ADCY6、KIF14、ADCY7、ADCY2、IFITM3、TYMS、IFI6、CDCA8、TIMP1、STAT2、SKA1、PKM、MX1、PRC1(ASE1)。
图2 PPI 网络分析结果和枢纽基因的鉴定结果
表2 Co-DEGs 中下调基因的GO 分析
枢纽基因的鉴定结果见图2B,结果显示,使用MCODE 从DEG 的PPI 网络中获得了3 个分值大于或等于6 的重要模块。
2.4 中枢基因的生存分析
根据每个基因的表达,绘制乳癌患者的无远处转移生存(distance metastasis free survival,DMFS)和总生存(overall survival,OS)曲线,结果分别见图3A 和图3B。结果显示,TYMS 的低表达[HR=0.52(0.29~0.91),P=0.021]、SKA1 的低表达[HR=0.57(0.32~1.01),P=0.049]、ADCY7 的低表达[HR=0.45(0.30~0.68),P<0.001]、MX1 的高表达[HR=2.07(1.17~3.66),P=0.011]与较差的OS 相关(图3A)。POLR3H 的低表达[HR=0.71(0.51~0.98),P=0.039]、CDCA8 的高表达[HR=1.67(1.37~2.04),P<0.001]、ASE1 的高表达[HR=2.1(1.72~2.57),P<0.001]、KIF14 的高表达[HR=1.86(1.33~2.59),P<0.001]、MX1 的高表达[HR=1.31(1.08~1.6),P=0.006]与较差的DMFS 相关(图3B)。
图3 枢纽基因的临床分析
选择数据集TCGA-BRCA 来验证20 个核心基因与乳腺癌临床分期之间的相关性。笔者比较了不同临床分期乳腺癌样本中核心基因的表达,结果见图3C,结果显示,在乳腺癌较晚分期中,TYMS 呈现低表达(P=0.0416),CDCA8 呈现高表达(P=0.003 66),PRC1(ASE1)呈现高表达(P=0.002 07),SKA1(P=0.000 909)呈现低表达,KIF14(P=0.000 711)呈现高表达。
3 讨论
本研究中分析了三个包含乳腺癌转移患者的基因数据集,进行了功能富集分析,表明了枢纽基因通过某些途径在转移的进展中发挥了作用。
结果表明,上调基因主要参与Ⅰ型干扰素信号通路、凋亡过程、粘着斑、蛋白质同二聚化活性Rap1 信号通路[18]和MAPK 信号通路[19-21],抗原加工和呈递,细胞粘附分子,而下调基因主要富集在p53 类介质、雌激素信号通路和趋化因子信号通路。这些发现与已报道的研究[19-29]高度一致,表明细胞凋亡、细胞迁移和粘附在乳腺癌进展过程中的重要作用。
先前的研究[19-29]已经为本研究所筛选的一些枢纽基因在乳腺癌症进展中的功能提供了大量证据。例如,TYMS 因其作为胸苷酸合酶的功能而被认为是5-氟尿嘧啶的靶标[23-24]。它与晚期乳腺癌[25]患者对化疗的耐药性和敏感性有关。CDCA8 是是有丝分裂的调节因子,SKA1 与有丝分裂有关,均被鉴定为乳癌远处转移的枢纽基因[5,26]。KIF14 通过负调节Rap1a-Radil 信号通路促进乳腺癌进展。参与细胞抗病毒的MX1与乳腺癌对淋巴结的侵袭有关[28-29]。尽管有报道称,PRC1(ASE1)促进了肺腺癌的发生和肝癌的早期复发[30-31],也和鼻咽癌的转移[32]和乳癌患者较差的无远处转移生存期相关[33],POLR3H 与促肿瘤作用相关[34-35]。它们在乳腺癌肿瘤转移中的确切功能仍然知之甚少,值得进一步深入研究。
综上所述,本研究确定了乳腺癌转移过程中的优势基因及其PPI 网络。其中一些基因从未被报道过影响乳腺癌的进展,因此可能作为潜在的药物靶点或生物标志物。然而,本研究仍具有局限性,即需要进一步的实验验证。需要检测更多基因来丰富网络,以便更全面地了解乳腺癌的转移通路。
[参考文献]
[1]Siegel RL,Miller KD,Jemal A.Cancer statistics,2019[J].CA Cancer J Clin,2019,69(1):7-34.
[2]Weigelt B,Peterse JL,Van′t Veer LJ.Breast cancer metastasis:Markers and models[J].Nat Rev Cancer,2005,5(8):591-602.
[3]Scully OJ,Bay BH,Yip G,et al.Breast cancer metastasis[J].Cancer Genomics Proteomics,2012,9(5):311-320.
[4]Chen X,Pei Z,Peng H,et al.Exploring the molecular mechanism associated with breast cancer bone metastasis using bioinformatic analysis and microarray genetic interaction network[J].Med (United States),2018,97(37):1-7.
[5]Cai Y,Mei J,Xiao Z,et al.Identification of five hub genes as monitoring biomarkers for breast cancer metastasis in silico[J].Hereditas,2019,156:20.
[6]Zheng T,Wang A,Hu D,et al.Molecular mechanisms of breast cancer metastasis by gene expression profile analysis[J].Mol Med Rep,2017,16(4):4671-4677.
[7]Edgar R,Domrachev M,Lash AE.Gene Expression Omnibus:NCBI gene expression and hybridization array data repository[J].Nucleic Acids Res,2002,30(1):207-210.
[8]Xia J,Gill EE,Hancock REW.NetworkAnalyst for statistical,visual and network -based meta -analysis of gene expression data[J].Nat Protoc,2015,10(6):823-844.
[9]Zhou G,Soufan O,Ewald J,et al.NetworkAnalyst 3.0:A visual analytics platform for comprehensive gene expression profiling and meta-analysis[J].Nucleic Acids Res,2019,47(W1):W234-W241.
[10]Ashburner M,Ball CA,Blake JA,et al.Gene ontology:Tool for the unification of biology[J].Nat Genet,2000,25(1):25-29.
[11]Altermann E,Klaenhammer TR.PathwayVoyager:Pathway mapping using the Kyoto Encyclopedia of Genes and Genomes(KEGG) database[J].BMC Genomics,2005,6:60.
[12]Huang DW,Sherman BT,Lempicki RA.Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources[J].Nat Protoc,2009,4(1):44-57.
[13]Szklarczyk D,Franceschini A,Wyder S,et al.STRING v10:Protein-protein interaction networks,integrated over the tree of life[J].Nucleic Acids Res,2015,43(D1):D447-D452.
[14]Shannon P,Markiel A,Ozier O,et al.Cytoscape:A software Environment for integrated models of biomolecular interaction networks[J].Genome Res,2003,13(11):2498-2504.
[15]Bader GD,Hogue CWV.An automated method for finding molecular complexes in large protein interaction networks[J].BMC Bioinformatics,2003,4:2.
[16]Gyorffy B,Lanczky A,Eklund AC,et al.An online survival analysis tool to rapidly assess the effect of 22 277 genes on breast cancer prognosis using microarray data of 1809 patients[J].Breast Cancer Res Treat,2010,123(3):725-731.
[17]Tang Z,Li C,Kang B,et al.GEPIA:A web server for cancer and normal gene expression profiling and interactive analyses[J].Nucleic Acids Res,2017,45(W1):W98-W102.
[18]Che YL,Luo SJ,Li G,et al.The C3G/Rap1 pathway promotes secretion of MMP-2 and MMP-9 and is involved in serous ovarian cancer metastasis[J].Cancer Lett,2015,359(2):241-249.
[19]Peng B,He R,Xu Q,et al.Ginsenoside 20(S)-protopanaxadiol inhibits triple-negative breast cancer metastasis in vivo by targeting EGFR-mediated MAPK pathway[J].Pharmacol Res,2019,142(151):1-13.
[20]Shao GL,Wang MC,Fan XL,et al.Correlation between Raf/MEK/ERKs signaling pathway and clinicopathological features and prognosis for patients with breast cancer having axillary lymph node metastasis[J].Technol Cancer Res Treat,2018,17(253):1-10.
[21]Yao Y,Fang ZP,Chen H,et al.HGFK1 inhibits bone metastasis in breast cancer through the TAK1/p38 MAPK signaling pathway[J].Cancer Gene Ther,2012,19(9):601-608.
[22]Lai YCC,Chang SJ,Kostoro J,et al.Vascular adhesion protein-1 as indicator of breast cancer tumor aggressiveness and invasiveness[J].Apmis,2018,126(9):755-761.
[23]Shan F,Liu YL,Wang Q,et al.Thymidylate synthase predicts poor response to pemetrexed chemotherapy in patients with advanced breast cancer[J].Oncol Lett,2018,16(3):3274-3280.
[24]Le K,Guo H,Zhang Q,et al.Gene and lncRNA co-expression network analysis reveals novel ceRNA network for triple-negative breast cancer[J].Sci Rep,2019,9(1):1-10.
[25]Guenter J,Abadi S,Lim H,et al.Evaluating genomic biomarkers associated with resistance or sensitivity to chemotherapy in patients with advanced breast and colorectal cancer[J].J Oncol Pharm Pract,2020,27(6):1371-1381.
[26]Fu Y,Zhou QZ,Zhang XL,et al.Identification of hub genes using co-expression network analysis in breast cancer as a tool to predict different stages[J].Med Sci Monit,2019,25:8873-8890.
[27]Ahmed SM,Thériault BL,Uppalapati M,et al.KIF14 negatively regulates Rap 1 a-Radil signaling during breast cancer progression[J].J Cell Biol,2012,199(6):951-967.
[28]Greenwood C,Metodieva G,Al -Janabi K,et al.Stat1 and CD74 overexpression is co -dependent and linked to increased invasion and lymph node metastasis in triple-negative breast cancer[J].J Proteomics,2012,75 (10):3031-3040.
[29]Aljohani AI,Joseph C,Kurozumi S,et al.Myxovirus resistance 1 (MX1) is an independent predictor of poor outcome in invasive breast cancer[J].Breast Cancer Res Treat,2020,181(3):541-551.
[30]Zhan P,Zhang B,Xi G min,et al.PRC1 contributes to tumorigenesis of lung adenocarcinoma in association with the Wnt/β-catenin signaling pathway[J].Mol Cancer,2017,16(1):1-15.
[31]Chen J,Rajasekaran M,Xia H,et al.The microtubule-associated protein PRC1 promotes early recurrence of hepatocellular carcinoma in association with the Wnt/β-catenin signalling pathway[J].Gut,2016,65(9):1522-1534.
[32]Yi L,Ouyang L,Wang S,et al.Long noncoding RNA PTPRG-AS1 acts as a microRNA-194-3p sponge to regulate radiosensitivity and metastasis of nasopharyngeal carcinoma cells via PRC1[J].J Cell Physiol,2019,234(10):19 088-19 102.
[33]Brynychova V,Ehrlichova M,Hlavac V,et al.Genetic and functional analyses do not explain the association of high PRC1 expression with poor survival of breast carcinoma patients[J].Biomed Pharmacother,2016,83:857-864.
[34]Mehra R,Vats P,Kalyana-Sundaram S,et al.Primary urethral clear -cell adenocarcinoma:Comprehensive analysis by surgical pathology,cytopathology,and next -generation sequencing[J].Am J Pathol,2014,184(3):584-591.
[35]Yán ez Y,Grau E,Rodríguez-Cortez VC,et al.Two independent epigenetic biomarkers predict survival in neuroblastoma[J].Clin Epigenetics,2015,7(1):1-14.
Genomic biomarker investigation for metastasis of breast cancer by combined bioinformatic gene expression profile analysis
WU Qiqiao1 ZHANG Shumin2 ZHENG Tingting3 LIU Juan2 HU Yong2 LIN Dengqiang4 DAI Mengting1 SUN Jing2▲
1.Department of Radiation Oncology, Fudan University Zhongshan Hospital (Xiamen Branch), Fujian Province, Xiamen 361006, China; 2.Department of Radiation Oncology, Fudan University Zhongshan Hospital, Shanghai 200030, China;3.Department of Radiology, Fudan University Zhongshan Hospital (Xiamen Branch), Fujian Province, Xiamen 361006,China; 4.Department of Urological Surgery, Fudan University Zhongshan Hospital (Xiamen Branch), Fujian Province,Xiamen 361006, China
[Abstract] Objective To identify potential genetic markers associated with breast cancer metastasis and to perform survival analysis.Methods The comprehensive Gene Expression Omnibus (GEO) was used to screen differentially expressed genes (DEGs) between metastatic and non-metastatic sites of breast cancer.Gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) pathway enrichment analyses were performed, the protein-protein interaction(PPI) network was constructed, and module analysis was performed using MCODE.Then, overall survival (OS) and distant metastasis free survival (DMFS) analysis of hub genes was performed by the Kaplan-Meier plotter online tool.The correlations between hub genes and clinical stages were validated through GEPIA online tool.Results A total of 295 DEGs were obtained, which were significantly enriched in the terms related to extracellular exosome, focal adhesion,cell division.KEGG pathway analysis showed that the significant pathways included mitogen-activated protein kinase(MAPK) signaling pathway, Rap1 signaling pathway, cell adhesion molecules (CAMs), antigen processing and presentation.The PPI network was established with 280 nodes and 357 edges.Ten modules were found from the PPI network.Conclusion A total of 20 genes were selected as hub genes, among which, low expression of TYMS, SKA1, ADCY7 and high expression of MX1 were related with worse overall survival.Low expression of POLR3H and high expression of CDCA8, ASE1, KIF14, MX1 were associated with worse DMFS.POLR3H and PRC1 may serve as new biomarker in diagnosing metastasis of breast cancer or potential drug target.
[Key words] Breast cancer; Metastasis; Bioinformatics analysis; Differently expressed genes; Therapeutic agent; Survival analysis
[中图分类号] R735.2
[文献标识码] A
[文章编号] 1674-4721(2021)12(c)-0022-07
通讯作者
(收稿日期:2021-08-04)