基于circRNA-miRNA-mRNA网络揭示circRNA在甲状腺乳头状癌中的作用

苗 欣，刘培发，菅雁兵，郗洪庆，陈 凛

解放军总医院第一医学中心 普通外科，北京 100853

甲状腺癌(thyroid cancer，TC)是临床上常见的一种内分泌恶性肿瘤，随着诊断技术的提升，其发病率在近几十年来持续上升[1-3]。甲状腺乳头状癌(papillary thyroid carcinoma，PTC)是TC中最常见的一种病理学类型，占所有TC病例的85%～95%，大部分患者预后良好，但高达30%的PTC病例在5年后复发[2,4]。尽管对PTC已有许多研究，但PTC是一种多因素疾病，其发生是一个以多分子异常为特征的复杂生物学过程，尚不清楚其高发病率的原因[5-6]。临床上对于PTC的诊断、术后治疗和预后评估也存在较大争议，因此寻找新靶点和特异性标志物以辅助诊断、判断预后是亟待解决的问题[7-9]。近年来，表观遗传调控在肿瘤发生发展中的重要作用已被广泛认识，大量证据表明ncRNA与肿瘤的生物学行为密切相关，如miRNA、circRNA和lncRNA的水平变化导致基因异常表达，可能导致PTC发生发展和转移[10-13]。本研究应用生物信息学方法对此进行探讨，首先差异分析基因表达综合数据库(Gene Expression Omnibus，GEO)来 源 的PTC相 关circRNA、miRNA和mRNA测序数据集，基于ceRNA理论并根据3种RNA之间的调控与被调控关系进一步构建circRNA-miRNA-mRNA网络，筛选出特定circRNA及其分子机制[14]。随后，利用 癌 症 基 因 组 图 谱(The Cancer Genome Atlas，TCGA)数据库进一步探讨具有上述差异mRNA变化的PTC患者生存情况，从而建立PTC预后相关circRNA-miRNA-mRNA调控子网络，为PTC的临床治疗及预后评估提供新靶标。

资料与方法

1 资料来源 在GEO数据库(http://www.ncbi.nlm.nih.gov/geo)下载PTC患者的三组相关数据集：circRNA数据集(GSE93522)、miRNA数据集(GSE113629)和mRNA数据集(GSE35570)。应用Limma和edgeR软件包分析差异表达的circRNA、miRNA和mRNA，将差异倍数(fold change，FC)＞2，即|logFC| ＞2且校正后P值(false discovery rate，FDR)＜0.05设定为筛选标准。登陆TCGA数据库(https://portal.gdc.cancer.gov/projects/TCGA)，下 载PTC患者的mRNA测序数据和临床数据。

2 筛选标准 选取2000年10月- 2019年10月收录在TCGA数据库中的PTC患者。纳入标准：1)mRNA测序数据完整；2)临床数据完整；3) 病理证实为PTC；排除标准：1)随访天数小于60 d；2)合并有其他恶性肿瘤。根据纳入和排除标准，本 研究共纳入PTC患者425例。

3 circRNA-miRNA-mRNA网络的构建 根据组织样本信息，将circRNA数据集GSE93522、miRNA数据集GSE113629和mRNA数据集GSE35570分别分为肿瘤组和对照组，分别在三个数据集中提取差异表达的circRNA、miRNA和mRNA。将差异表达的circRNA通过癌症特异环状RNA数据 库(cancer-specific circRNA database，CSCD；https://gb.whu.edu.cn/cscd)数据库预测相应的靶标miRNA，结合GSE113629中筛选得到的差异表达miRNA进一步筛选目标miRNA。然后，利用miRDB(MicroRNA Target Prediction And Functional Study Database；http://www.mirdb.org)、miRTarBase(http://mirtarbase.mbc.nctu.edu.tw/php/index.php)和TargetScanHuman 7.2(http://www.targetscan.org/vert_72)数据库预测目标miRNA对应的mRNA，再与GSE35570中差异表达的mRNA进一步取交集以识别目标mRNA。基于ceRNA理论，通过Cytoscape 3.7.1软件构建circRNA-miRNA-mRNA调 控网络。

4 预后相关circRNA-miRNA-mRNA调控网络的构建 从TCGA下载PTC患者的临床数据及mRNA表达谱数据，根据目标mRNA表达的中位数将患者分为高风险组和低风险组，以PTC患者死亡或疾病进展为观察终点，Kaplan-Meier法比较高低风险组无进展生存期(progress free survival，PFS)的差异，筛选出与患者预后相关的核心mRNA，从而构建与PTC患者预后相关的circRNAm iRNA-mRNA调控子网络。

5 统计学方法 使用R3.6.0软件进行统计学分析及相应图形绘制。采用Limma软件包对基于微阵列数据的circRNA、miRNA和mRNA表达水平进行校正，采用edgeR软件包对肿瘤组与对照组样本之间的差异表达进行分析，采用pheatmap软件包进行热图绘制，采用venn软件包对各差异表达的miRNA、mRNA与靶标取交集，采用survival软件包对患者进行PFS曲线的绘制与预后分析。P ＜0.05为差异有统计学意义。

结 果

1 PTC患者一般情况 根据纳入排除标准最终有425例患者纳入分析。其中男性103例、女性322例，平均年龄为(61.21 ± 7.92)岁；其中白种人 342例，黑种人54例，黄种人29例。见表1。

2 差异表达circRNA、miRNA和mRNA 通过分析PTC肿瘤组与对照组的circRNA数据集GSE93522，发现115个上调circRNA和21个下调circRNA(图1)。差异分析miRNA数据集GSE113629，发现19个上调miRNA和31个下调miRNA(图2A)。差 异 分 析mRNA数 据集GSE35570，发现242个上调mRNA和176个下调 mRNA(图2B)。

图1 PTC中差异表达circRNA的热图分析F ig.1 Heatmap of differentially expressed circRNAs in the PTC patients

图2 PTC相关差异表达miRNA(A)及mRNA(B)的热图分析F ig.2 Heatmap of differentially expressed miRNAs (A) and mRNAs (B) in the PTC patients

表 1 PTC患者一般情况Tab. 1 Clinicopathological characteristics of patients with PTC

3 circRNA-miRNA-mRNA网 络 的 构 建 通过CSCD数据库分别对上述得到的115个上调和21个下调circRNA进行靶向预测，筛选出2159个靶向miRNA。将这2159个靶标miRNA与miRNA数据集GSE113629得到的50个差异miRNA取交集，筛选出37个目标miRNA(图3A)。再利用TargetScan、miRTarBase和miRDB数据库对这37个目标miRNA进行靶向预测，筛选出1327个靶标mRNA。将这1327个靶标mRNA与mRNA数据集GSE35570得到的418个差异mRNA取交集，筛选出34个目标mRNA(图3B)。根据ceRNA调控原理并运用Cytoscape 3.7.1软件，我们最终筛选出29个circRNA、9个miRNA和12个mRNA构 建circRNA-miRNA-mRNA可视化网络(图4)。

图3 差异表达miRNA(A) 及mRNA(B) 的筛选F ig.3 Venn diagrams for differentially expressed miRNAs (A) and mRNAs (B)

4 PTC预后相关mRNA的筛选 从TCGA数据库下载PTC患者的临床数据及mRNA表达谱数据，根据纳入排除标准，最终有425例患者纳入分析。我们对上述用于构建circRNA-miRNAmRNA网络的12个mRNA进行预后分析，分别以各mRNA的表达中位数为界值，将425例患者分为高风险组和低风险组，利用R语言“survival”包的Kaplan-Meier方法进行预后分析，最终显示KIT、SPRY4、SFN这3个mRNA与PTC患者的PFS相关，患者信息见表2。SFN、SPRY4的高表达及KIT的低表达可能会导致PTC患者预后不良 。见图5。

5 预 后 相 关circRNA-miRNA-mRNA子 网 络 构建 根据ceRNA调控原理并运用Cytoscape 3.7.1软件，由KIT、SPRY4、SFN这3个mRNA我们最终筛选出10个circRNA、3个miRNA和3个mRNA构建PTC患者预后相关circRNAmiRNA-mRNA可视化子网络(图6)。其中10个c ircRNA的基本详细模式图见图7。

图6 PTC患者预后相关circRNA-miRNA-mRNA网络的构建Fig.6 Analysis of PTC prognosis related circRNA-miRNA-mRNA network

图7 筛选出的circRNA基本模式图F ig.7 Structural patterns of the ten circRNAs

讨 论

RNA是基因组编码的遗传信息的直接输出，在细胞功能中起着至关重要的作用。近年来越来越多研究表明ncRNA的失调与癌症发病机制密切相关，大规模基因组测序结果也揭示了基因改变在PTC发生发展中的关键作用[15-17]。circRNA近年来已成为肿瘤领域的研究热点，它是一类新的内源性非编码RNA，在序列上高度保守，并具有组织特异性表达特点[18]，这些特性使circRNA成为多种疾病的潜在分子生物标志物。miRNA是内源性保守类非编码小RNA，通过多种机制负反馈调节基因表达和细胞活性[19]。多种miRNA与PTC的发生、发展密切相关，其异常表达在临床应用中已作为生物标志物协助PTC良恶性的鉴别、淋巴结转移的判断、复发监测、治疗及预后预测，如miRNA-146、miRNA-222等[20-23]。circRNA可以作为miRNA分子海绵与内源RNA竞争，从而调节miRNA靶基因表达，但对其中ceRNA调控网络在PTC发生发展中的作用目前不清楚[24]。因此，本研究通过分析PTC非编码RNA高通量测序数据，基于ceRNA理论构建circRNA-miRNA-mRNA三元互作网络并进行预后分析，筛选参与调控PTC预后的潜在mRNA，为PTC临床治疗及预后评估提供新靶标。

本研究中首先从GEO数据库获得PTC相关circRNA、miRNA和mRNA测序数据集进行差异分析，基于ceRNA理论并根据3种ncRNA之间的调控与被调控关系进一步构建circRNA-miRNAmRNA网络，筛选出特定ncRNA及调控网络。在circRNA-miRNA-mRNA网络中，circRNA通过CircBase网站和CSCD数据库筛选，与差异miRNA取交集，而miRNA-mRNA的靶向关系来源于实验证实的miRTarbase 数据库，从而保证调控网络的可信性。随后，利用TCGA数据库进一步探讨可能影响患者预后的差异mRNA，发现SFN、SPRY4的高表达及KIT 的低表达可能会导致PTC患者预后不良。SFN基因编码细胞周期检查点蛋白14-3-3σ，作为一种肿瘤抑制蛋白在肿瘤中多呈低表达状态或活性降低，主要与基因启动子高甲基化和泛素化的增强导致蛋白量减少、自身活性调控因子异常引起的活性改变相关[25-28]。然而，由于14-3-3σ的生物学功能广泛、调控因子众多且涉及的调控机制复杂，在某些肿瘤呈现出相反的结果：SFN的过度表达刺激了癌细胞的增殖和生长且与较差的预后相关，如肺腺癌、胆管癌[27,29-30]。本研究发现SFN高表达患者的预后较低表达患者差，表明SFN在PTC中可能发挥着癌基因的作用，是潜在治疗靶点。受体酪氨酸激酶介导的信号通路与PTC的发生、发展密切相关，该通路抑制剂可能适用于靶向治疗带有这些基因突变的PTC患者，如药物索拉非尼和乐伐替尼。KIT基因的编码产物为一种酪氨酸激酶受体，与细胞的生长和生存密切相关，在临床上发现KIT突变与黑色素瘤、急性白血病及胃肠瘤的发生有关，PTC患者样本和混合淋巴结样本中也检测到KIT基因的多种突变，包括胚系突变和体细胞突变[31-33]。本研究结果显示KIT基因很可能是PTC发生、发展中一个非常重要的调控基因，并在临床诊治中具有指导意义。SPRY4基因编码的蛋白属于富含半胱氨酸和脯氨酸的蛋白质家族，是受体转导MAPK信号通路的抑制剂，SPRY4蛋白的失调参与了多种人类癌症的发生发展：在前列腺癌中受到DNA甲基化而下调；在人乳腺癌中异位表达SPRY4可抑制乳腺癌细胞系的细胞增殖和迁移；通过调节细胞增殖参与结直肠癌的发生和发展，可能是结直肠癌患者的潜在诊断和预后指标；还可以抑制非小细胞肺癌的细胞增殖和迁移[34-37]。MAPK信号通路是驱动甲状腺癌发病和发展的主要通路，然而SPRY4在控制PTC发生发展中作用机制的研究仍较少，本研究结果表明SPRY4可促进PTC的发展，高表达SPRY4的患者预后较差。

由KIT、SPRY4、SFN这3个mRNA建立PTC预后相关circRNA-miRNA-mRNA调控子网络，拓展了circRNA通过ceRNA网络调控PTC发生发展的分子机制，为PTC临床治疗及预后评估提供新靶标，未来还可能作为免疫治疗的新靶点。尽管本文系统揭示了circRNA-miRNA-mRNA在PTC患者预后评估中的重要作用，但依旧存在一定的局限性：首先，预测分子靶点仍需大量多中心的循证医学证据支持；其次，本研究是基于数据库中RNA高通量测序的分析结果，缺少临床、细胞、动物功能学实验，数据分析结果仍需要进一步验证。