基于GEO 和TCGA 数据库对肺腺癌差异表达基因的生物信息学分析

叶 汇, 孙 哲, 周丽婷, 齐 雯, 叶 琳

（吉林大学公共卫生学院劳动卫生与环境卫生教研室，吉林 长春 130021）

肺癌是中国发病率最高的癌症，约占全部恶性肿瘤的18%，肺癌导致死亡也是恶性肿瘤患者死亡的主要原因［1］。肺腺癌（lung adenocarcinoma，LUAD）是肺癌的主要组织学类型，约占肺癌总类型的40%［2-3］。LUAD 早期为原位腺癌和微浸润腺癌，转移可能性小，手术后复发率低；晚期LUAD转移和复发的风险高，严重者可危及生命。因此，早诊断和早治疗是提高LUAD 患者预后的关键。由于LUAD 患者早期无明显的临床表现，患者大多在LUAD 晚期确诊，尽管可使用靶向治疗和免疫疗法，LUAD 患者5 年生存率仍低于15%［4］。因此，发现影响LUAD 发生发展的关键基因有助于LUAD 的早期诊断和治疗。生物信息学技术、基因表达谱技术和微阵列技术的发展使癌症在诊断及治疗方面取得了巨大的进步。本研究通过对高通量基因表达（Gene Expression Omnibus，GEO） 数据库和癌症基因组图谱（The Cancer Genome Atlas，TCGA）数据库的基因芯片进行联合分析，筛选差异 表 达 基 因 （differentially expressed genes，DEGs），识别影响LUAD 发生发展及预后的关键基因并排序，分析前10 位关键基因的生物学功能及其对预后的影响，为发现LUAD 的生物标志物和治疗靶点及改善患者预后提供新的思路。

1 资料与方法

1.1 LUAD 基因芯片数据获取于GEO 数据库（https：//www. ncbi. nlm. nih. gov/geo/） 中以“adenocarcinoma”为关键词进行检索。筛选条件：①LUAD；②组织来源为人的组织；③研究类型为数组表达式分析（expression profiling by array）。通过文献筛选，采用GSE118370 和GSE136043 数据集。GSE118370 数据集为基于GPL570 构建的芯片，包含6 例正常肺组织和6 例LUAD 组织；GSE136043 数据集为基于安捷伦-026652 构建的芯片，包含5 例正常肺组织和5 例LUAD 组织。于TCGA 数据库（https：//www. cancer. gov/ccg/research/genome-sequencing/tcga） 中下载LUAD的转录组测序技术（RNA sequencing，RNA-Seq）数据，共包含594 例样本，其中59 例正常样本、535 例肿瘤组织样本和515 例病案样本。

1.2 DEGs 的提取和分析于GEO 和TCGA 数据库下载原始表达矩阵，基因组浏览器（https：//genome.ucsc.edu/）下载人类基因文件。采用Perl语言进行原始表达矩阵转换；使用R 软件中svaR包对GSE118370 和GSE136043 数据集进行合并和批次矫正；使用limmaR包筛选差异基因，计算倍数变化（fold change，FC）和错误发生率（false discovery rate，FDR），以|log10FC|>1 和|log10FDR|>1，矫正后P<0.05 为筛选条件；采用ggplot2 包绘制火山图，pheatmap 包绘制热图；采用venn 图网络分析工具确定公共的DEGs 进行后续分析。

1.3 DEGs 的基因本体（Gene Ontology，GO）功能富集分析和京都基因与基因组百科全书（Kyoto Encyclopedia of Genes and Genomes，KEGG）信号通路富集分析使用R 软件中clusterProfiler R 包对LUAD 的DEGs 进行GO 功能富集分析，使用DAVID 数据库（https：//david. ncifcrf. gov/） 对LUAD 的DEGs 进行KEGG 信号通路富集分析，基因集的截断标准为P<0.05。

1.4 DEGs 网络构建和关键基因筛选将上述获取的DEGs 导入STRING 数据库（https：//stringdb. org/） 进行蛋白- 蛋白相互作用（proteinprotein interaction，PPI）网络构建，评估潜在的蛋白质相互作用关系。采用Cytoscape 软件将PPI 网络可视化并去除游离蛋白质节点，CytoHubba 插件计算每个蛋白质节点的连接度并按照降序排列，排序前10 位的DEGs 为与LUAD 发病相关的关键基因。

1.5 关键基因mRNA 表达水平和蛋白表达情况分析在GEPIA（https：//gepia.cancer-pku.cn/）数据库中查询关键基因在LUAD 患者组织（n=483）和正常组织（n=347） mRNA 表达水平。以每百万转录本（transcripts per million，TPM） 计量关键基因mRNA 表达水平。于人类蛋白质图谱（The Human Protein Atlas， HPA） （https：//www.proteinatlas.org/）数据库中获取LUAD 组织和正常肺组织中关键基因蛋白的免疫组织化学染色图像并分析关键基因蛋白表达情况。

1.6 不同分期LUAD 患者肺组织中关键基因的表达分析使用GEPIA 数据库查询关键基因在不同分期LUAD 患者肺组织中的表达差异，以获取基因表达与肿瘤恶性程度的相关关系。

1.7 关键基因免疫浸润分析和生存分析使用TIMER 数据库（https：//timer. comp-genomics.org/）分析关键基因的表达与肿瘤纯度、B 淋巴细胞、CD8 T 淋巴细胞、CD4 T 淋巴细胞、巨噬细胞、中性粒细胞和树突状细胞之间的相关性。将关键基因导入Kaplan-Meier plotter 数据库（https：//kmplot.com/analysis/），按照关键基因表达量的中位数将患者划分为低表达组和高表达组，以总生存期（overall survival，OS）和患者不良预后生存期（poor prognosis survival，PPS） 为患者生存终点，筛选最佳探针，限制分析组织为腺癌组织，分析关键基因的表达与患者生存期的相关关系。以P<0.01 为差异有统计学意义。

2 结 果

2.1 LUAD 中DEGs 筛选在由TCGA 数据库中获得的数据集中，筛选出LUAD的DEGs为9 998 个，其中7 569 个DEGs 表达上调和2 429 个DEGs 表达下调（图1）。 将由GEO 数据库中获得的GSE118370 和GSE136043 数据集合并，筛选出LUAD 的DEGs 为1 425 个，其中601 个DEGs 表达上调和824 个DEGs 表达下调（图2）。将由GEO数据库和TCGA 数据库中获得LUAD 的DEGs 取交集，获得LUAD 的DEGs 954 个，其中428 个DEGs 表达上调和526 个DEGs 表达下调（图3）。

图1 TCGA 数据库中LUAD 的DEGs Fig. 1 DEGs of LUAD in TCGA Database

图2 GEO 数据库中LUAD 的DEGs Fig. 2 DEGs of LUAD in GEO Database

图3 DEGs 的Venn 图Fig. 3 Venn diagram of DEGs

2.2 DEGs 的GO 功能富集分析和KEGG 信号通路富集分析GO 功能富集分析结果显示：LUAD的DEGs 在生物过程（biological process，BP） 方面主要富集于外包被结构组织、细胞外基质（extracellular matrix，ECM）、细胞外结构组织、循环系统中的血管形成和上皮细胞的增殖和调节等方面；在细胞组分（cellular components，CC） 方面主要富集于含胶原的ECM、细胞交界处、顶端质膜、细胞外侧血浆膜和细胞基部等方面；在分子功能（molecular function，MF） 方面主要富集于氨基多糖结合、硫化合物结合、信号受体激活剂活性、化学趋避活性、肝素结合、G 蛋白偶联受体活性、蛋白质结合、受体配体活性、生长因子活性和ECM 结构形成等方面。见图4。

图4 LUAD 中DEGs 的GO 功能富集分析Fig. 4 GO fuctional enrichment analysis on DEGs in LUAD

KEGG 分析结果显示：LUAD 的DEGs 主要富集于蛋白质消化和吸收通路、ECM-受体相互作用通路、病毒蛋白与细胞因子和细胞因子受体的相互作用通路、轴突引导通路、血管平滑肌收缩通路及细胞因子受体相互作用通路等方面。见图5。

图5 LUADDEGs 的KEGG 信号通路富集分析Fig. 5 KEGG signaling pathway enrichment analysis on DEGs of LUAD

2.3 关键基因筛选结果LUAD 的DEGs 编码蛋白质PPI 网络由940 个节点和454 个边组成。连接度排名前10的关键基因分别为DNA 拓扑异构酶Ⅱα（topoisomerase Ⅱ alpha，TOP2A）、果蝇纺锤体异常基因（abnormal spindle microtubule assembly，ASPM）、细胞周期蛋白B1（cyclin B1，CCNB1）、人类细胞分裂周期相关基因8 （cell division cycle associated 8，CDCA8）、含杆状病毒IAP 重复序列蛋 白 5 （baculoviral IAP repeat containing 5，BIRC5）、苏氨酸激酶（aurora A，AURKA）、驱动蛋白超家族成员20A （kinesin family member 20A，KIF20A）、中心体相关蛋白55（centrosomal protein 55，CEP55）、着丝粒蛋白F （centromere protein F，CENPF） 和微管组织因子（targeting protein for Xklp2，TPX2）。见表1 和图6。

表1 连接度排名前10 位的关键基因Tab.1 Top 10 key genes in connectivities

图6 关键基因可视化图Fig. 6 Visualization map of key genes

2.4 关键基因mRNA 表达水平与正常肺组织比较，LUAD 组织中TOP2A、CCNB1、CDCA8、BIRC5、AURKA、KIF20A、CEP55、CENPF 和TPX2 mRNA 表达水平均增加（P<0.01）。见图7 和8。

图7 GEPIA 数据库分析2 组关键基因mRNA 表达情况Fig. 7 Expressions of key gene mRNA in two groups analyzed by GEPIA Database

图8 关键基因mRNA 表达情况Fig. 8 Expressions of key gene mRNA

2.5 2 组关键基因蛋白表达情况HPA 数据库中与LUAD 组织和正常肺组织有关的蛋白免疫组织化学染色结果显示：与正常肺组织比较，LUAD 组织 中 TOP2A、 CCNB1、 CDCA8、 BIRC5、AURKA、KIF20A、CEP55、CENPF 和TPX2 蛋白表达均增加。由于HPA 数据库中缺少ASPM 蛋白表达，因此该基因在肺组织中的表达情况无法验证。

2.6 不同分期LUAD患者肺组织中关键基因mRNA表达水平CCNB1、CDCA8、BIRC5、AURKA、KIF20A、 CEP55 和TPX2 mRNA 在不同分期LUAD 患者肺组织中的表达水平比较差异均有统计学意义（P<0.01）。与Ⅰ、Ⅱ和Ⅲ期LUAD 患者比较，Ⅳ期LUAD 患者肺组织中CCNB1、CDCA8、AURKA、KIF20A、CEP55 和TPX2 mRNA 表达水平增加（P<0.01）；与Ⅰ、Ⅱ和Ⅳ期LUAD 患者比较，Ⅲ期LUAD 患者肺组织中BIRC5 mRNA表达水平增加（P<0.01）。见图9。

图9 不同分期LUAD 患者肺组织中关键基因mRNA 表达Fig. 9 Expressions of key gene mRNA in lung tissue of LUAD patients with different stages

2.7 关键基因免疫浸润分析筛选的10 个关键基因表达与B 淋巴细胞浸润均呈负相关关系（－0.253≤r≤－0.104，P<0.01）；TOP2A、ASPM、CDCA8、BIRC5、CEP55、CENPF 和TPX2 表达与中性粒细胞浸润呈正相关关系（0.049≤r≤0.165，P<0.01）；CCNB1 和AURKA 表达与CD4 T 淋巴细胞、巨噬细胞和树突状细胞浸润呈负相关关系（－0.210≤r≤－0.100，P<0.01）。见图10和11。

图10 TOP2A、 ASPM、CCNB1、CDCA8 和BIRC5 基因免疫浸润分析Fig. 10 Immunoinfiltration analysis on TOP2A， ASPM，CCNB1，CDCA8，and BIRC5 genes

图11 AURKA、 KIF20A、CEP55、 CENPF 和TPX2 基因免疫浸润分析Fig. 11 Immunoinfiltration analysis on AURKA， KIF20A，CEP55， CENPF，and TPX2 genes

2.8 关键基因生存分析不良预后生存期曲线提示CDCA8 高表达会增加LUAD 恶化风险（P<0.01）。OS 曲线提示TOP2A、CCNB1、CDCA8、BIRC5、AURKA、KIF20A、CEP55、CENPF 和TPX2 高表达会增加患者死亡风险（P<0.01）。见图12 和13。

图12 关键基因不良进展生存曲线Fig. 12 Survival curves of adverse progression of key genes

图13 关键基因OS 曲线Fig. 13 OS curves of key genes

3 讨 论

可根据生长状况和病理特点等将LUAD 分为恶性程度逐级递增的原位腺癌、微浸润腺癌和浸润性腺癌。早期LUAD 多为原位腺癌和微浸润腺癌，转移的可能性小，手术后复发率极低，对患者生存影响较小。但浸润性腺癌转移和复发的风险高，严重者可能危及生命［5］。

本研究结果显示：LUAD 中DEGs 与ECM 组织、上皮细胞的增殖和调节及上皮-间质转化（epithelial-mesenchymal transition，EMT） 等BP、细胞-细胞连接和膜筏等CC 及糖胺聚糖结合和肽结合等MF 有关；DEGs 主要富集于ECM-受体相互作用通路。EMT 影响肿瘤细胞的迁移和侵袭，在肿瘤的发生发展中发挥重要作用。EMT 过程与ECM 受体相互作用通路有关。ECM-受体相互作用通路在肿瘤脱落、黏附、降解、运动和增生过程中起重要作用［6-9］。研究［10］显示：ECM-受体相互作用通路与LUAD 转移有关联。本研究筛选出影响LUAD 发生发展的前10 位关键基因TOP2A、CCNB1、 BIRC5、 TPX2、 AURKA、 ASPM、CENPF、KIF20A、CDCA8 和CEP55，均在EMT过程中发挥重要作用，可能通过影响ECM-受体相互作用通路促进LUAD 的进展。

不同分期LUAD 关键基因表达水平差异提示可以通过测定关键基因的表达水平确定LUAD 的分期，达到早诊断和早治疗的目的。肿瘤微环境中与肿瘤进展和预后不良相关的免疫细胞包括中性粒细胞、巨噬细胞和T 淋巴细胞等。本研究结果显示：关键基因可能参与LUAD 对肿瘤微环境的免疫应答，尤其是对B 淋巴细胞和中性粒细胞的免疫应答。B 淋巴细胞作为肿瘤微环境的重要组成部分，存在于肿瘤的各个阶段，B 淋巴细胞活化表明免疫治疗有效［11］。中性粒细胞是炎症通过破坏组织导致癌症的发生过程的重要成分，其浸润提示癌症转移［12］。本研究结果显示：关键基因表达抑制了B 淋巴细胞活化，进而促进中性粒细胞浸润，因此关键基因TOP2A 等高表达可能抑制细胞免疫并促进肿瘤转移。

TOP2A 是影响细胞分裂的关键基因，可作用于有丝分裂过程，促进染色体凝聚和分离并维持染色体结构。TOP2A 在宫颈癌、肝癌和胃癌等癌组织中高表达［13-15］。PEI 等［16］研究显示：TOP2A 是β-连环蛋白的共激活剂，可通过激活EMT 过程以促进肿瘤转移。WANG 等［17］发现：TOP2A 在早期肺癌患者中可改善其预后，TOP2A 表达与树突状细胞浸润水平有密切关联，与本研究结果一致。TOP2A 可能是LUAD 患者的预后生物标志物和潜在治疗靶点。

CCNB1 是细胞周期所必需蛋白，在细胞G2期向M 期的过渡阶段和完成M 期的过程中发挥作用。BAO 等［18］研究显示：CCNB1 过表达促进LUAD细胞增殖，微小RNA（microRNA，miR）-139-5p可以通过抑制β-连环蛋白的表达下调CCNB1 表达，进而抑制LUAD 细胞增殖、迁移和侵袭。研究［19］显示：巨噬细胞浸润与LUAD 患者的预后呈负相关关系，提示LUAD 患者预后不良。本研究结果显示：CCNB1 高表达与巨噬细胞浸润呈正相关关系，因此CCNB1 可用于提示LUAD 患者预后不良。

BIRC5 可调节胚胎细胞胚胎发生中的有丝分裂过程和细胞凋亡及恶性疾病的癌变。其可通过抑制EMT 进而抑制原位卵巢中的原发性肿瘤生长和原位卵巢癌小鼠模型中的腹膜转移［20］。WANG 等［21］研究显示：BIRC5 对肝癌有良好的预后价值，其高表达与免疫细胞浸润差和OS 短有显著关联。本研究结果显示：BIRC5 可能是LUAD 预后不良的生物标志物。

TPX2 是微管相关蛋白，是有丝分裂纺锤体组装和功能所必需的。TPX2 过表达可导致染色体不稳定，促进癌症发生［22］。TPX2 对微管的调控可能介导NSCLC 细胞的细胞骨架重塑，进而激活肿瘤细胞的EMT 过程和MMP 表达，使非小细胞肺癌细胞发生转移［23］。本研究结果显示：TPX2 高表达与LUAD 的临床分期、免疫细胞浸润和不良进展呈正相关关系。

ASPM 是果蝇纺锤体异常基因的人类同源基因，对胚胎神经母细胞有丝分裂过程中纺锤体功能至关重要。研究［24］显示：ASPM 通过磷脂酰肌醇3激酶（phosphatidylinositol-3 kinase，PI3K）/蛋白激酶B （protein kinase B，Akt） 信号通路促进EMT， 调节LUAD 细胞转移。 ASPM 表达与KIRC 和LIHC 中的肿瘤浸润、 B 淋巴细胞、CD8 T 淋巴细胞和M2 巨噬细胞有较强的相关性［25］。本研究结果显示：LUAD 预后不良和免疫细胞浸润水平增加提示ASPM 表达可作为LUAD临床结局和免疫细胞浸润的新预后生物标志物。

CDCA8 在人类胚胎干细胞和癌细胞中转录被激活，但在正常组织中低表达或不表达。研究［26］显示：在癌症发生发展过程中CDCA8 过表达。HU 等［27］发现：CDCA8 是LUAD 细胞中的致癌因子，其受上游miR-133b 调控。miR-133b 可以通过靶向CDCA8 抑制LUAD 细胞的增殖、侵袭和迁移。

AURKA 属于丝氨酸/苏氨酸激酶家族，其编码的蛋白是一种细胞周期调节激酶，在染色体分离过程中与纺锤体极的微管形成和稳定有关。AURKA 通过参与癌细胞的增殖、EMT 转移、细胞凋亡和癌症干细胞更新促进肿瘤发生发展［28］。基因扩增、转录激活和抑制蛋白质降解均可能导致癌症组织中AURKA 表达水平升高［29］。AURKA在LUAD 组织中高表达，与B 淋巴细胞和树突状细胞浸润水平呈负相关关系，同时B 淋巴细胞和树突状细胞浸润水平也与LUAD 患者的预后有关联。

KIF20A 是 驱 动 蛋 白 超 家 族 （kinesin superfamily protein，KIF） 成员，KIF 主要参与细胞器、蛋白质复合物、mRNA 的运输及有丝分裂和减数分裂过程中染色体及纺锤体的运动等。研究［30］显示：KIF20A 在纤维肉瘤组织中高表达，敲低KIF20A 调控增殖、迁移和侵袭并诱导G2/M停滞及细胞凋亡。在异种移植小鼠模型中，KIF20A 敲低抑制肿瘤生长和转移［30］。ZHAO 等［31］发现：LUAD 细胞中KIF20A 表达水平明显升高，敲低KIF20A 可大幅减少LUAD 细胞增殖，诱导肿瘤细胞凋亡，表明KIF20A 可能是LUAD 的治疗靶点。

CEP55 通过与参与转移、血管生成和增殖途径的蛋白质相互作用参与肿瘤发生。研究［32］显示：CEP55 通过激活P38、 丝裂原活化蛋白激酶（mitogen-activated protein kinase， MAPK） 和ERK 1/2 通路促进EMT 和肿瘤发生。FU 等［33］发现：Hsa-miR-195-5p/E2F7/CEP55轴可促进LUAD细胞凋亡，抑制增殖和自噬，降低LUAD 细胞对吉西他滨的耐药性，提示Has-miR-195-5p/E2F7/CEP55 轴可能是LUAD 的潜在治疗靶点。

CENPF 属于外层动粒蛋白，可构成纺锤体微管附着点，与染色体的运动和分离有关，并在染色体微管动力学作用和纺锤体检验电信号传导中发挥重要功能。研究［34］显示：CENPF 是多功能动粒蛋白，在细胞有丝分裂和分化中发挥重要作用。CENPF 可促进HCC 细胞系的迁移和EMT 的进展。敲低CENPF 在体外抑制HCC 细胞活性，在体内抑制异种移植肿瘤。SUN 等［35］发现：CENPF 通过激活PI3K-Akt- 雷帕霉素靶蛋白（mammalian target of rapamycin，mTOR） C1 信号传导促进乳腺癌骨转移。TANG 等［36］发现：CENPF 和雌激素受体β（estrogen receptor-β，ER-β） 2/5 的高表达与LUAD 的发展有关，敲低CENPF 可抑制由ER-β2/5 途径介导的LUAD 进展，提示CENPF 表达可作为LUAD 不良预后和免疫细胞浸润的新生物标志物。

综上所述， TOP2A、 ASPM、 CCNB1、CDCA8、 BIRC5、 AURKA、 KIF20A、 CEP55、CENPF 和TPX2 是参与LUAD 发生进展过程的关键基因，可能通过加速EMT 过程促进LUAD 发展，其高表达提示免疫细胞浸润、LUAD 进入晚期及LUAD 患者预后不良和死亡风险升高。