甲状腺乳头状癌潜在关键基因的生物信息学分析

赵珍，邢琰，张风华，刘妍，景尚华

（1.河北医科大学第四医院耳鼻咽喉-头颈外科，河北石家庄050000；2.河北省石家庄人民医院普通外科三病区，河北石家庄050000；3.河北省人民医院腺体外科，河北石家庄050000）

甲状腺癌是近年来增长最快的内分泌恶性肿瘤，发病率大于5%[1]。其中，80%～85%为甲状腺乳头状癌（papillary thyroid carcinoma，PTC）。PTC的治疗方式有手术切除，放射性碘消融等。尽管大多数PTC 患者预后良好，PTC 的5年总生存率为95%， 但仍有20%～30% 的患者会出现复发和进展[2]。

PTC 的形成是一个十分复杂的生物化学过程，其特征是各种分子的异常[3-4]。大规模基因组研究表明遗传改变在PTC 的发生发展中起着关键作用[5-7]。目前对其持续高发的原因尚未明确[8]。目前，我国多数地区将甲状腺超声筛查甲状腺癌作为常规体检项目。由于受检查医师主观认知水平限制较大，仍有很多的甲状腺癌误诊或漏诊，导致过度医疗和延误治疗。临床应用于诊疗PTC 的分子生物标志物在特异度和敏感度方面很低，对于临床的价值不大[9-11]。PTC 发病人数的迅速增长迫切需要开发新的有效的用于PTC 诊断和预后的分子标志物和治疗靶标。

本研究试图找到检测PTC 患者的新的高效的分子生物指标，使之能够提高临床的诊断和预后，成为潜在的治疗靶标。通过检测从Gene Expression Omnibus （GEO） 数据库（GSE60542，GSE33630 和GSE3467）中PTC 和正常甲状腺组织之间的差异表达基因，应用生物信息学方法分析基因的表达数据谱，对筛选出的基因进行基因本体论（gene ontology，GO）功能注释分析。建立蛋白质-蛋白质相互作用（protein-protein interaction，PPI） 网络来鉴别与PTC 相关的关键基因。这些关键基因的生存分析使用在线数据库Cbioportal 工具进行，并进一步研究验证。

1 材料与方法

1.1 数据来源

本研究中分析的基因表达数据集是从GEO 数据库（https://www.ncbi.nlm.nih.gov/geo） 获得的。从数据库中检索到总共2 079 个关于人类甲状腺癌的系列文章。经过仔细审查，选择了3 个基因表达谱（GSE60542，GSE33630 和GSE3467）。

1.2 差异基因的数据处理

使用GEO2R 在线分析工具（https://www.ncbi.nlm.nih.gov/geo/geo2r） 检测PTC 和正常样品之间的差异基因。对每个数据集进行统计分析，调整P值和|logFC|以符合临界标准。筛选出P＜0.05 和|logFC|≥2.0 的基因被视为差异基因，其中logFC≥2 代表基因表达上调，logFC≤-2 代表基因表达下调。Venn 图网络工具（bioinformatics.psb.ugent.be/webtools/Venn）确定重合的基因。

1.3 差异基因的GO功能富集分析

使用注释， 可视化和集成发现数据库（DAVID） 工具（https://david.ncifcrf.gov） 进行差异基因的GO 注释分析,供研究者挖掘大量基因背后的生物学价值。P＜0.01 和基因计数≥10 具有统计学意义。

1.4 PPI网络建设与关键基因识别

对于筛选出的差异基因，使用搜索工具（STRING） 数据库（http://string-db.org） 用于PPI 信息分析。以＞0.15 的组合分数提取PPI。通过Cytoscape 软件（www.cytoscape.org）将PPI 网络可视化处理。在本研究中，筛选出前10 个具有高连通性的节点基因，即为PTC 的关键基因。

1.5 关键基因对于PTC的生存分析

cBioPortal （https://www.cbioportal.org） 在线工具可以通过可视化的形式展示癌症研究样本的基因组数据，同时也可以帮助研究人员探索样本、基因和通路之间的遗传变化，并与临床结果相结合。本研究使用cBioPortal 评估筛选出的关键基因在PTC的预后价值。P＜0.01 表示具有统计学意义的结果。

1.6 关键基因在临床标本中的验证

1.6.1 临床标本收集 52 例PTC 组织和与之匹配的相邻正常甲状腺组织样本（距肿瘤边缘≥2 cm）是我院进行手术的患者中收集的。手术标本立即保存在液氮中直至使用。本研究中所有患者均签署知情同意书。我院医学伦理委员会批准此项研究。1.6.2 qRT-PCR 检测 TRIzol（碧云天）提取组织中的总RNA。NanoDrop ND-1000 分光光度计测量RNA 浓度及纯度。通过PrimeScript™RT 试剂盒（碧云天） 合成 cDNA。 采用 SYBR Green qPCR SuperMix-UDG（碧云天）和ABI 7500 实时PCR 系统（Applied Biosystems） 进行实时PCR。使用GAPDH作为内参，进行3 次独立的实验。KIT 的比较采用ΔCt 方法。引物序列由Sangon Biotech （上海，中国）合成。

1.7 统计学处理

用SPSS2 0.0 软件和GrapPad Prism 8.0 软件进行统计分析。计量数据用均数±标准差（±s）表示。关键基因在PTC 癌组织和癌旁正常甲状腺组织表达的比较采用配对t检验。PTC 癌组织中关键基因表达与临床病理特征相关性比较采用独立样本t检验。P＜0.05 为差异有统计学意义。

2 结果

2.1 筛选差异基因

在3 个基因表达谱中，GSE60542 包含33 个PTC 样本和30 个正常甲状腺正常样本，GSE33630和GSE3467 分别包含49 个、9 个PTC 样本和45 个、9 个正常甲状腺样本。根据标准，从基因芯片GSE60542中总共鉴定出307个差异基因，包括168个上调基因和139 个下调基因。在基因芯片GSE33630中，鉴定出了323 个差异基因，其中197 个基因上调，而126 个基因下调。从GSE3467 中鉴定出了155 个差异基因，包括93 个上调基因和62 个下调基因。通过将PTC 样本与正常甲状腺样本进行比较来识别所有的差异基因。随后，进行Venn 分析以获得差异基因轮廓的交集（图1）。最后，组中有102 个差异基因表达，其中62 个明显上调的基因和40 个明显下调的基因。

图1 3个GEO数据集共有的差异基因的Venn图A：上调的差异基因；B：下调的差异基因Figure 1 Venn diagram of the differentially expressed genes shared by the three datasetsA:Up-regulated genes;B:Downregulated genes

2.2 差异基因的GO功能富集分析

GO 是一个国际标准化的基因功能分类系统，全面描述生物体内基因和基因产物的性质。GO 共有3 个本体论，描述基因的分子功能、细胞的位置和所涉及的生物过程。使用DAVID 对差异基因进行GO 功能富集分析（图2）。

图2 差异基因的GO功能富集分析Figure 2 Significantly enriched GO terms of the differentially expressed genes

2.3 PPI网络建设与关键基因识别

Cytoscape 软件用于计算每个蛋白质节点的连通程度，通常有较高连通性的节点对整个网络的稳定性具有更重要的作用。STRING 检索先前确定的差异基因提取PPI（图3），并通过对提取的PPI网络进行可视化处理（图4）。此研究中，前10 个具有高连通性的节点基因被鉴定为关键基因。

图3 用差异表达基因构建的PPIFigure 3 PPI constructed based on the differentially expressed genes

图4 具有较高连通性的10个关键基因Figure 4 The 10 hub genes with high connectivity

2.4 关键基因对于PTC的生存分析

cBioPortal 是一种可视化和分析多维度癌症基因组数据的在线工具。其整合了多个肿瘤基因组研究的数据，还包括临床预后等表型的信息。使研究人员不仅能够获取大量的癌症基因组学数据，还能够将这些数据转化为图形化的结果，使复杂的癌症基因组学资料更易理解和接受。将10 个关键基因于cBioPortal 进行生存分析示，其中KIT 具有统计学意义（P＜0.01）（图5）。

图5 KIT与患者生存的关系Figure 5 Relationship between KIT and survival of the patients

2.5 临床标本验证结果

KIT 表达水平在PTC 组织较癌旁正常甲状腺组织表达明显降低（P＜0.001）（图6）。在PTC 组织中，KIT 与患者临床分期和淋巴结转移的差异有统计学意义（均P＜0.05），而与年龄、性别、组织分化方面的差异无统计学意义（均P＞0.05）（表1）。

表1 PTC中KIT的表达和临床病理参数的关系（±s）Table 1 Relationship between the expression of KIT and clinicopathologic parameters in PTC(±s)

表1 PTC中KIT的表达和临床病理参数的关系（±s）Table 1 Relationship between the expression of KIT and clinicopathologic parameters in PTC(±s)

参数年龄（岁）≤55＞55例数（n）值t P 30 26 8.64±1.53 7.29±2.48 0.958 0.709性别男女25 31 7.93±1.32 6.74±1.75 0.892 0.653分化程度高/中低29 27 7.84±1.38 6.36±2.52 1.112 0.132临床分析I/II III 38 18 5.18±1.05 8.26±1.53 4.273 0.008淋巴结转移无有19 38 4.83±1.32 9.25±1.46 3.753 0.023

图6 KIT在PTC及癌旁正常甲状腺组织中的表达Figure 6 KIT expressions in PTC and normal adjacent thyroid tissues

3 讨论

PTC 是甲状腺癌最常见的肿瘤亚型，主要发生在中老年人，发生的中位年龄是50 岁[12]。PTC 的主要特征之一是它能够侵入邻近的淋巴管等结构[13-15]。约有27%在就诊时出现同侧的颈侧区淋巴结转移[16]。手术是PTC 主要的治疗手段，通常具有良好的预后[17-19]。但是，有些临床病理和背景特征可能导致不良预后[20]。在分子水平上进一步阐明PTC 在肿瘤行为方面的发病机制，找到PTC 诊疗和预后的分子生物标志物有着重要的临床意义。

本研究通过公开的数据库资源进行基因表达和蛋白质-蛋白质表达的分析，筛选鉴定与PTC 相关的潜在的关键基因。依据GEO 数据库的基因表达数据，筛选出PTC 与正常甲状腺组织之间的差异基因，并确定了上调和下调的差异基因。进一步进行了差异基因的GO 功能富集分析。通过构建PPI 网络研究差异基因的相互关系，10 个基因被鉴定为关键基因。使用cBioPortal 在线工具来预测关键基因表达与PTC 患者预后之间的关系，发现PTC中KIT 基因的下调与PTC 患者的预后不良有关。

收集我院52 例PTC 手术患者样本，通过实时PCR 检测KIT 在癌组织和癌旁正常组织中的表达。结果表明KIT 在PTC 组织中的表达明显降低，统计具有统计学意义（P＜0.01），与以上生信分析结果一致。KIT 的表达与PTC 患者的临床分期和淋巴结转移相关，临床分期III 期的PTC 比I/II 期，淋巴结转移比未转移的PTC 患者KIT 表达下调的更明显，差异具有统计学意义（P＜0.05）。

KIT 是一种III 型受体酪氨酸激酶，位于染色体4q12，存在于Cajal 间质细胞、生殖细胞、造血祖细胞等多种细胞中，具有调节细胞的增殖和迁移等功能，在人体发生发育和生理调节中起着重要的作用[21]。KIT 在癌症的发生和增殖中起着至关重要的作用[22-23]。具有活化的KIT 的肿瘤是伊马替尼和其他选择性酪氨酸激酶抑制剂的潜在靶标[24-26]。在胃间质瘤中发现了突变的KIT，KIT 突变是胃间质瘤预后不良的标志，突变的KIT 的存在可预测疾病对伊马替尼的反应性[27-30]。

基于GEO 数据库获得的基因表达数据集进行数据分析，筛选PTC 与正常甲状腺组织的差异基因。鉴定出10 个基因为PTC 潜在的关键基因，其中KIT 下调与PTC 患者的不良预后相关,需要进一步的研究证实研究的结果。KIT 可能是治疗PTC 的治疗靶标和分子生物标志物。