基于GEO数据库的膀胱癌差异表达基因的生物信息学分析

孔钰琳，荣胜忠，王慧单，高 梦，李晓霞

(牡丹江医学院流行病与统计学教研室，黑龙江 牡丹江 157011)

膀胱癌是泌尿系统中最常见的恶性肿瘤之一。据WHO全球肿瘤流行病统计数据(GLOBOCAN 2020)显示，2020年世界范围内新发膀胱癌约有573278例，死亡约212536例，在男性人群中，膀胱癌的发病率高居第6位且居癌症相关性死亡原因的第9位[1]。目前，膀胱癌的治疗以手术为主，但其复发率和转移率较高，严重影响了患者的预后[2]。因此，探索通过研究膀胱癌在分子水平的发生发展机制，筛选具备诊断和预后潜力的分子标志物、研发新型基因药物及改善患者预后至关重要。

本研究通过从GEO数据库获取膀胱癌基因芯片数据进行分析，筛选出差异表达基因(differentially expressed genes,DEGs)，并利用基因本体论(gene ontology,GO)富集分析注释基因功能，采用京都基因与基因组百科全书(Kyoto Encyclopedia of Genes and Genomes,KEGG)通路分析对DEGs的信号通路进行富集，然后构建蛋白互作网络(protein-protein interaction network,PPI network)，同时结合Kaplan-Meier plotter数据库对关键基因的生存预后情况进行分析，以期发现膀胱癌潜在的治疗靶点。

1 材料与方法

1.1 基因芯片数据的获取进入GEO数据库(https://www.ncbi.nlm.nih.gov/geo/)，在搜索框内输入“Bladder Cancer”进行检索，纳入标准为：(1)样品来自人膀胱癌组织标本；(2)研究对象包含膀胱癌患者和正常对照；(3)研究类型为“Expression profiling by array”。从GEO数据库中筛选下载数据集(GSE7476,GSE37815和GSE13507)，其中GSE7476基于GPL570平台，GSE37815与GSE13507基于GPL6102平台。GSE7476包括9例癌组织，3例正常组织；GSE37815包括18例癌组织，6例正常组织；GSE13507包括188例癌组织，68例正常组织。

1.2 方法

1.2.1 数据处理及差异基因的筛选 利用GEO数据库自带的在线分析工具GEO2R(https://www.ncbi.nlm.nih.gov/geo/geo2r/)处理原始数据，将数据分为膀胱癌组和正常组进行分析。以|log FC(foldchange)|≥1.5且校正后P值<0.05为筛选标准分别对三组数据集进行DEGs筛选。log FC≥1.5的基因作为上调差异基因(UDEGs)，log FC≤-1.5的基因作为下调差异基因(DDEGs)。筛选出各数据集中DEGs后，通过Venny在线工具(http://bioinfogp.cnb.csic.es/tools/venny/)确定GSE7476、GSE37815和GSE13507三个数据集的相交部分。

1.2.2 DEGs的GO富集分析与KEGG通路分析 通过DAVID在线分析工具(The Database for Annotation,Visualization and Integrated Discovery,https://david.ncifcrf.gov/,version 6.8)对筛选出的DEGs 进行GO和KEGG分析。按照P<0.05为截取标准筛选GO分析和 KEGG通路分析结果。

1.2.3 PPI网络的构建和关键基因的筛选 通过String(Search Tool for the Retrieval of Interacting Genes,https://string-db.org/,version 11.0)在线工具构建DEGs的PPI网络。将PPI网络导入Cytoscape(version 3.7.1)软件进行可视化，利用Cytohubba插件，基于四种拓扑算法：最大聚集中心(MCC)、最大邻域组件(MNC)、度(Degree)、边缘渗透组件(EPC)分别筛选出PPI网络中前十个重要节点，最后将上述四种算法所获得的结果取交集，将重叠的差异表达基因作为关键枢纽基因。

1.2.4 关键基因的预后分析 通过Kaplan-Meier plotter(http://kmplot.com/analysis/)数据库进行评估关键基因的预后价值。每个基因根据 mRNA表达值自动将癌症患者分为高表达和低表达两组进行比较，生存分析采用Log-rank检验，以P<0.05为差异有统计学意义。

2 结果

2.1 DEGs的筛选对GSE7476、GSE37815和GSE13507三个数据集中的DEGs取交集，最终筛选出101个DEGs，其中95个UDEGs(表1)，6个DDEGs(表2)。不同数据集的DEGs筛选结果的火山图见图1。

图1 膀胱癌中的差异表达基因

表1 膀胱癌上调表达基因

表2 膀胱癌下调表达基因

2.2 DEGs的GO分析和KEGG信号通路分析GO富集分析，膀胱癌DEGs在细胞学过程(biological processes,BP)方面主要是参与肌肉收缩、细胞黏附、肌丝滑动等生物过程(图2A)。在细胞成分(cell component,CC)方面，主要分布于细胞外区、蛋白质细胞外基质、细胞外基质、黏着斑等组织(图2B)。分子功能(molecular function,MF)方面，DEGs主要有肌肉结构成分、细胞骨架蛋白绑定、肌动蛋白结合、钙调蛋白的结合等分子功能(图2C)。DEGs的KEGG信号通路分析显示，DEGs在心肌收缩、肥厚型心肌病、扩张型心肌病、黏着斑等通路中富集表达(图2D)。

图2 DEGs的GO分析和KEGG通路分析

2.3 PPI网络构建及关键基因筛选为了进一步确定筛选出的101个DEGs之间的关系，通过STRING11.0进行分析。PPI网络共涉及节点101个，边数157个，平均节点度为3.11，蛋白质相互作用富集P值小于1.0e-16。将PPI网络导入Cytoscape软件进行可视化，标注PPI网络的上调基因(橙色)及下调基因(绿色)(图3)，并利用MCODE插件筛选PPI网络中的关键功能模块，包括10个节点和45条边(图4)。利用Cytohubba插件，根据MCC，MNC，Degree和EPC这四种拓扑算法筛选出排名前十的关键基因(表3和图5)，确定了9个重叠差异表达基因为关键基因，即TPM1、ACTC1、ACTA2、TPM2、TAGLN、CALD1、LMOD1、MYH11、CNN1。

图3 膀胱癌差异表达基因编码蛋白作用网络图

图4 PPI网络中的关键模块

表3 基于Cytohubba插件筛选差异表达关键基因

图5 基于Cytohubba插件筛选差异表达关键基因的韦恩图

2.4 关键基因的预后分析利用Kaplan-Meier plotter数据库进行预后分析，对于每个基因根据mRNA表达值自动将癌症患者分为高表达和低表达两组进行比较，P<0.05表示对膀胱癌患者有显著影响。在Kaplan-Meier plotter数据库中共有405例膀胱癌患者用于分析总生存率，分析结果显示9个关键基因(TPM1、ACTC1、ACTA2、TPM2、TAGLN、CALD1、LMOD1、MYH11、CNN1)的上调差异表达均显著影响膀胱癌患者的总生存率，见图6。

图6 膀胱癌关键基因的预后分析

3 讨论

膀胱癌是我国男性泌尿生殖系统最常见的恶性肿瘤，其中根据是否侵犯逼尿肌分为肌层浸润性膀胱癌和非肌层浸润性膀胱癌，尽管有经尿道膀胱肿瘤切除术、放疗、化疗和BCG免疫疗法等多种治疗手段，但膀胱癌的生存率在近三十年里都未有改善[2]。膀胱癌作为一种异质性疾病，已造成巨大的社会负担。因此，膀胱癌的早期筛查以及早期诊断显得尤为重要，探索有效的肿瘤标志和治疗靶点对膀胱癌的诊疗十分必要。

本研究从GSE7476、GSE37815和GSE13507三个数据集中共筛选出101个DEGs，基于生物信息学对膀胱癌DEGs的功能和信号通路进行分析，并通过构建PPI网络以及利用Cytoscape软件获得与膀胱癌密切相关的9个关键基因，Kaplan-Meier plotter数据库探究关键基因与膀胱癌预后之间的关系，发现9个关键基因皆与病人的预后显著相关。基于以上结果，推测这9个关键基因可能作为肿瘤基因参与了膀胱癌的发生发展。

本研究的KEGG信号通路分析表明，DEGs在心肌收缩、肥厚型心肌病、扩张型心肌病、黏着斑等通路中富集表达。研究发现ACTC1、TPM1的突变与肥厚型心肌病的发展密切相关[3-4]，主要因为ACTC1、ACTA2、TPM1、TPM2、TAGLN都属于肌动蛋白家族，肌动蛋白是一类广泛分布于生物界的高度保守的蛋白质，对细胞增殖、细胞迁移、细胞凋亡具有调控作用[5]，且肌动蛋白装配动力学介导的细胞突起是细胞侵袭的关键步骤[6]。

TAGLN是一种肌动蛋白结合蛋白和肿瘤抑制因子[7]，此前已被报道参与了膀胱癌的发生发展过程，TAGLN在体外能抑制膀胱癌细胞的增殖和侵袭能力，在体内能抑制肿瘤发生[8-9]。

而ACTA2、ACTC1、TPM1、TPM2、CALD1、LMOD1、MYH11和CNN1虽在膀胱癌中尚未见报道，但它们的致癌作用在部分肿瘤中已被证实。ACTC1在胶质母细胞瘤中表达显著上调，可抑制癌细胞的迁移[10]。并且ACTC1与胶质母细胞瘤的预后相关，可作为神经胶质瘤中的一种新型的预后标志物[11]。CALD1作为一种细胞骨架蛋白[12]，它的错配导致的1-CAD表达上调会引发胶质瘤新生血管功能异常[13]。另外研究发现LMOD1可作为一种新的胃癌生物标志物和治疗靶点，通过调控FAK-Akt/mTOR通路诱导EMT[14]。

CNN1在乳腺癌中起肿瘤抑制作用，作为miR-106b-5p的下游靶标，其表达水平在乳腺癌组织和细胞系中下调，与miR-106b-5p的表达呈负相关。miR-106b-5p可以通过靶向CNN1和Rho / ROCK1通路促进乳腺癌细胞癌变[15]。

综上所述，本研究分析了TPM1，ACTC1，ACTA2，TPM2，TAGLN，CALD1，LMOD1，MYH11，CNN1这9个基因的上调差异表达均显著影响膀胱癌患者的总生存率。因此，本研究将为今后膀胱癌发生发展以及预后的分子机制研究提供初步的理论依据。但本研究尚有局限性，仍需要结合细胞实验和临床样本来进一步验证这些关键基因在体内外的作用机制。