基于GEO甲状腺癌芯片数据的生物信息学分析

袁小艳，梁 韡，张 舒，周 萍

(成都大学四川抗菌素工业研究所抗生素研究与再评价四川省重点实验室，成都 610052)

甲状腺癌是甲状腺上皮细胞的恶性肿瘤。随着人们饮食结构、营养状况及生态环境等因素的变化，近年来，甲状腺癌发病率呈逐年上升趋势[1]。据统计，甲状腺癌已成为我国城市人群癌症患病率第4位、女性易发肿瘤第5位的恶性实体肿瘤[2]。目前，甲状腺癌的明确病因未见报道，探索甲状腺癌的发病机制及筛查生物标志物对其早期显示具有重要意义。

生物信息学是结合分子生物学与信息技术的新兴交叉学科[3]，对揭示肿瘤发病的分子机制具有重要意义。基因芯片作为一种高效、大规模获取生物信息的新兴技术，能够广泛采集疾病基因芯片表达谱数据。本文采用生物信息学工具对公共基因芯片数据库(gene expression omnibus，GEO) 中甲状腺癌表达谱芯片进行数据分析，为相关基因的生物学功能及其参与甲状腺癌发生、发展的分子机制提供理论依据。

1 材料与方法

1.1 材料 表达谱芯片数据来源于美国国立生物技术信息中心(NCBI)GEO数据库中甲状腺癌数据集，编号分别为GSE33630、GSE65144和GSE85457，包括76例甲状腺癌组织和61例正常组织。数据集由GPL570 芯片平台采集甲状腺癌及癌旁组织的基因芯片表达谱数据。

1.2 方法

1.2.1 数据处理及差异基因分析 3组原始数据集采用R软件及Affy、Limma、Pheatmap、Ggplot2等软件包进行数据处理，RMA算法进行背景校正、标准化及表达值计算。处理后数据采用Fold-change(FC)和T-test进行差异基因筛选，定义| log(FC)|≥1，P<0.05为有效基因。筛选3组数据集中共有有效基因作为最终差异表达基因进行后续分析。

1.2.2 差异表达基因的功能富集分析 生物学信息注释数据库(database for annotation，visualization and integrated discovery，DAVID)是在线生物信息学分析工具(https://david.ncifcrf.gov/)，可对基因或蛋白进行系统综合的生物功能分析和注释。差异表达基因利用David6.8进行基因本体论(gene ontology，GO)和KEGG(Kyoto Encyclopedia of Genes and Genomes) 通路的富集分析。

1.2.3 差异表达基因的PPI分析 String10.0数据库(http://string-db.org/)由已知或预测蛋白质相互作用的数据组成，包括直接和间接的蛋白质间相互作用。差异表达基因采用String10.0进行蛋白-蛋白相互作用(PPI)网络分析，设置有效结合分数为大于0.4。分析结果通过Cytoscape 3.5.1建立网络模型，以drgree得分筛选中心基因，并采用Cytoscape中MCODE插件对网络模型进行评价，选择scores>3且node>4模块中的基因进行通路富集分析。

2 结 果

2.1 差异表达基因筛选 经R软件分析，分别从数据集GSE33630、GSE65144和GSE85457中提取到差异表达基因3 844、2 822、1 724个。取3个数据集差异表达基因进行综合分析，得到相同差异表达基因383个(图1)；其中上调基因217个，下调基因166个。

图1 3组基因芯片数据共表达的差异基因关系图

2.2 差异表达基因GO富集分析 GO分析结果显示，差异表达基因的GO注释主要分为生物过程(BP)、细胞组成(CC)和分子功能(MF)3个部分。其中上调基因主要富集在细外基质组织、胶原纤维组织、调节细胞增殖等生物学过程，下调基因主要富集在调节脂肪细胞分化、肾发育、内分泌系统发育等生物学过程，部分结果见表1。

表1 GO富集分析部分结果

2.3 差异表达基因KEGG信号通路分析 KEGG分析结果显示，上调基因富集的信号通路(23条)多于下调基因富集的信号通路(6条)。上调基因主要集中在ECM受体相互作用、小细胞肺癌、PI3K-Akt信号通路、血小板活化、细胞周期等信号通路(图2)；下调基因主要富集在癌症转录失调、甲状腺激素合成、TGF-β信号通路等信号通路，见图3。

图2 上调基因KEGG通路分析结果

2.4 差异表达基因PPI分析 基于String数据库，筛选出CDC6、AURKA、FEN1、MCM4、MYC 5个degree得分较高的hub基因(表2)。Cytoscape软件 MCODE插件共筛选出显著模块3个(表3、图4)。模块1中涉及基因主要富集在DNA复制、细胞周期等信号通路，模块2中涉及基因主要富集在ECM受体相互作用、PI3K-AKT等信号通路，模块3中涉及基因主要富集在甲状腺激素合成、自身免疫性甲状腺疾病等信号通路。

图3 下调基因KEGG通路分析结果

基因名称关系强度分数拓扑系数MYC63 054300.072 463 77CDC617 620240.166 666 67AURKA13 140220.190 637 72FEN15 976210.251 984 13MCM43 426190.288 888 89

表3 3个模块中涉及基因通路分析结果

图4 模块1、2、3的蛋白互作关系

3 讨 论

甲状腺癌是内分泌系统和头颈部肿瘤中最常见的恶性肿瘤[4]。通过对肿瘤组织与正常组织基因表达谱的生物信息学分析发现，差异表达基因对血管发育、细胞增殖、转化生长因子β(TGF-β)反应及细胞生物合成调控等生物学过程都有一定影响。ECM受体相互作用、PI3K-Akt信号通路、各种癌症的疾病通路、细胞周期等信号通路是上调基因的主要富集通路，下调基因主要参与癌症转录失调、甲状腺激素合成、TGF-β等信号通路的调控。PI3K-Akt通路作为细胞内信号转导的重要通路之一，在多种肿瘤中异常活化[5]；研究表明，PI3K-Akt通路的异常活化可抑制肿瘤细胞凋亡，使肿其增殖、转移和侵袭能力显著增强[6]。肿瘤细胞的生长由各种生长因子、激素及细胞外基质等物质构成的微环境维持，这些因素的改变会使肿瘤的发生、发展过程及对药物敏感性发生变化。CAMs是包含整合素、选择素和免疫球蛋白等在内的一类黏附因子[7]，可通过介导细胞黏附作用及与相应配体结合，减少肿瘤细胞“失巢凋亡”和诱导肿瘤血管生成，促进癌细胞的生长和迁移[8]。TGF-β作为一种细胞生长因子，能够诱导肿瘤上皮细胞发生上皮间质转化及关键血管生成因子的表达[9]，促进肿瘤细胞的侵袭。甲状腺激素合成与甲状腺功能状态息息相关，其水平变化可作为恶性肿瘤患者其病情，以及预后判断的辅助指标[10]。因此，监测以上信号通路可能有助于预测甲状腺癌的进展及其对药物的敏感性。

同时，本文通过差异表达基因的PPI分析，筛选出MYC、CDC6、AURKA、FEN1和MCM4这5个关键基因。MYC基因家族及其产物可促进细胞增殖、永生化、去分化和转化等过程，在多种肿瘤形成过程中处于重要地位[11]；目前发现的MYC基因家族成员包括c-myc、n-myc及l-myc，其中c-myc是一种涉及细胞分化、生长、增殖以及凋亡的转录因子，与肿瘤的形成、发展和转移密切相关[12]，下调MYC表达的治疗可能会抑制甲状腺肿瘤细胞的增殖。CDC6是启动DNA复制的重要因子之一，在多种癌症组织和细胞中表达异常上调[13]，其主要通过调控起始DNA复制、激活S-M监测点及促进有丝分裂滑脱，影响肿瘤细胞的增殖及转移[14]；然而，目前尚无关于CDC6异常表达与甲状腺癌关系的研究报道，因此其可能成为未来甲状腺癌诊断的分子标志物和治疗的新靶点。AURKA基因编码的丝氨酸苏氨酸激酶，通过参与中心体的分离和成熟以及纺锤体两极的建立，确保细胞有丝分裂中染色体的正确分离和胞质分裂的顺利完成[15]；AURKA的高表达能够直接或间接地激活多种致癌蛋白或使多种抑癌蛋白失活，CHEN等[16]的研究表明AURKA在衰减P35和P73的肿瘤抑制功能上起着重要作用。FEN1是一种多功能蛋白，参与了细胞DNA 复制、修复和凋亡等过程[17]；研究发现，FEN1在多种肿瘤癌症中高表达，抑制其表达可降低肿瘤的增殖和侵袭能力，逆转肿瘤耐药性[18]。MCM4是启动真核基因组DNA复制所必需的微小染色体维持蛋白，据文献报道，MCM4 是与P53、KI-67、PCNA 等相比具有同等地位或更具临床指导意义的标记物[19-20]，可为研究甲状腺癌发生、发展，指导临床治疗及判断预后等方面提供更广阔的思路。此外，PPI网络模块分析显示，甲状腺癌的发展主要与DNA复制、细胞周期、PI3K-Akt信号通路、ECM受体相互作用、焦点黏附及甲状腺激素合成等信号通路有关，这些通路与肿瘤细胞的增殖、侵袭和转移密切相关，在甲状腺癌的发生、发展、浸润和转移等方面起着重要作用，为研究甲状腺癌的发病机制及肿瘤标志物的筛选提供了理论依据。

综上所述，本文通过生物信息学方法分析甲状腺腺癌及癌旁组织的基因表达谱芯片数据，发现影响甲状腺癌发生、发展的重要通路及关键基因，为阐述甲状腺癌的发病机制及诊断提供了全新视角，并为甲状腺癌靶向抑制剂的开发提供新的方向。然而，本文根据大数据预测的甲状腺癌相关信号通路及关键基因及它们是否也在其他癌症中发挥作用，还需通过分子生物学研究进行进一步验证。