m5C甲基化相关基因预测肾透明细胞癌的预后

陈晓超，卜恒涛，梁 超，李 杰

南京医科大学第一附属医院泌尿外科，江苏 南京 210029

肾细胞癌是世界上第三常见的泌尿系恶性肿瘤，主要包括肾透明细胞癌、乳头状肾细胞癌、嫌色细胞癌、集合管癌和未分类肾细胞癌。其中肾透明细胞癌约占80%～90%［1］。尽管目前开发了许多靶向药物和免疫抑制药物，但肾根治性切除术仍然是金标准治疗方法［2］。此外，肾透明细胞癌的特点是复发率和转移率高［3］，这对患者的健康和生活质量造成严重影响。因此，寻找肾透明细胞癌中的相关基因以构建预测模型至关重要，同时这些基因也可以作为未来分子靶向治疗的潜在靶点。

最新研究表明，RNA修饰是一种新兴的基因调控机制，其主要包括m6A、m5C 和m1A，m6A 和m5C修饰技术是160多种化学修饰中最主要和最具代表性的两种RNA 转录后修饰类型［4-5］。其中，m5CRNA甲基化修饰主要依赖于甲基转移酶、去甲基化酶和结合蛋白，m5C-RNA甲基化修饰发挥多种生物学功能。如调节mRNA 的转运，增加RNA 的稳定性，调节蛋白质翻译，维持RNA 正常结构等，近年来许多研究表明，m5C-RNA 甲基化修饰参与调节多种肿瘤的发生发展、侵袭和转移。如在肺癌中，m5C 甲基化相关基因NSUN3、NSUN4 在肿瘤组织中的表达高于正常肺组织，且与临床病理特征及生存率有关，此外，还与6 种主要免疫细胞的浸润有关，从而调节肺癌的肿瘤免疫微环境［6］。在膀胱尿路上皮癌中，m5C 甲基化转移酶NSUN2 表达较正常组织高，从而通过增强癌基因肝癌衍生生长因子（HDGF）mRNA 的稳定性来促进疾病的发生发展。而对头颈鳞状细胞癌的研究显示，m5C RNA甲基化修饰不仅有助于促进疾病的进展，NSUN5、DNMT1和DNMT3A 的表达模式还有助于预测疾病的预后。因此，寻找m5C甲基化相关基因对肿瘤患者的早期诊断和预后评估具有重要意义

本研究基于肿瘤基因图谱（the Cancer Genome Atlas，TCGA）数据库下载的肾透明细胞癌患者的数据，通过生物信息学和统计学分析，构建了由m5C甲基化相关基因组成的生存预后模型，分析了与肾透明细胞癌生存有关的独立预后影响因素，同时用此模型将肾透明细胞癌患者分为高低分险组，通过GO 和KEGG 分析对高低风险组中的差异基因进行了分析，进一步揭示这些差异基因的生物学功能和潜在的信号通路。

1 资料和方法

1.1 资料

从TCGA 数据库中下载了572 个肾透明细胞癌样本和72 个正常样本的转录组数据用于进一步分析，临床数据也从TCGA 数据库下载（https：//www.cbioportal.org/）［9］。在cBioportal 数据库中研究了选定的4 种基因与拷贝数变异之间的关系（https：//www.cbioportal.org/）［10］。在GEPIA2 数据库找到这4 种基因的不同器官表达特征［11］。肾透明细胞癌和癌旁正常组织样本来自2008年4月—2019年4月南京医科大学第一附属医院接受根治性肾切除术的患者，本研究经过了南京医科大学第一附属医院伦理委员会授权（伦理批准号：2021-SR-430）。所有患者都签署了协议，允许他们的组织样本和相应临床信息用于进一步的研究。

1.2 方法

1.2.1 m5C-RNA甲基化相关基因分析

从已发表的文献中找到13 个m5C 甲基化相关基因，其中包括11 个甲基转移酶（NOP2、NSUN2、NSUN3、NSUN4、NSUN5、NSUN6、NSUN7、DNMT1、DNMT3A、DNMT3B、TRDMT1）、1 个去甲基化酶（TET2）和1 个结合蛋白（ALYREF）。然后系统地提取并分析了13 种调节基因的表达谱，以及临床病理参数。

1.2.2 m5C甲基化相关基因差异表达的数据处理

提取572 例肿瘤组织和72 例正常肾组织中m5C 甲基化相关基因的表达矩阵和临床资料。使用limma 软件包的R 版本（4.1.0）确定肿瘤和对照组之间差异表达的m5C甲基化相关基因。以P＜0.05和|log2FC|＞1 作为所有检验的显著性阈值。此外，用热图和小提琴图展示肿瘤组和对照组表达的差异。

1.2.3 m5C甲基化相关基因的共识聚类分析

用R（4.1.0）的“ConsensusClusterPlus”软件包对筛选出来的13 个m5C 甲基化相关基因进行共识聚类分析。此外，还对不同的聚类进行生存分析，以确定肾透明细胞癌样本中的最佳聚类。然后进行了Kaplan-Meier 生存分析，以揭示不同集群之间总体生存率的差异。比较临床病理参数，包括性别、分级、年龄、TNM 分期和不同聚类之间的关系。

1.2.4 蛋白质相互作用及功能分析

为了进一步筛选出中心基因，使用STRING数据库检索相互作用基因以分析差异表达的m5C甲基化相关基因，构建蛋白质-蛋白质相互作用（protein-protein interaction，PPI）网络（http：//stringdb.org/）［12］。

1.2.5 构建预测模型

首先，计算危险比（hazard ratio，HR，危险组个体和非危险组个体的比值）和95%置信区间（confidence interval，CI），通过单变量Cox回归分析识别适当的候选基因。其次，使用Lasso-Cox 回归分析，挑选出合适的能够预测肾透明细胞癌的预后基因，风险值计算公式为：风险值其中，n代表mRNA 的数目，coef 是系数，xi是基因的表达量）。接下来，根据风险值中位数，将TCGA 的肾透明细胞癌队列分为高风险组和低风险组，再用Kaplan-Meier生存分析/风险预测模型估计候选风险基因的预后价值。

1.2.6 GO和KEGG富集分析

使用R 软件包“clusterProfiler”、“enrichplot”和“ggplot2”对这13个选定基因进行GO 和KEGG 富集分析，以确定分子功能（molecular function，MF）、细胞成分（cellular components，CC）和生物过程（biological process，BP）并寻找潜在的信号通路。

1.2.7 细胞培养

肾癌细胞系（786-O、A498）和人肾小管上皮细胞系（HK-2）均购自中国科学院（上海）细胞库。其中，786-O 和HK-2 细胞在RPMI 1640 中培养，A498细胞在DMEM/F12中培养。所有细胞系均在37 ℃、含5%CO2的培养箱中培养。

1.2.8 总RNA的提取和qRT-PCR

使用TRIzol 试剂（美国赛默飞公司）从培养的细胞和组织样本中提取总RNA，并将其反转录为cDNA。然后使用StepOne Plus 实时PCR 系统和SYBR Premix Ex Taq 试剂进行qRT-PCR 实验，目的基因引物见表1。

表1 实时定量PCR的引物序列Table 1 The sequences of primers for RT-qPCR

1.3 统计学方法

统计分析使用R 软件（4.1.0）进行。此外，Perl编程语言（版本5.30.2）用于数据处理。Kaplan-Meier生存曲线分析用于分析总体生存率（overall survival，OS）。单变量和多变量Cox 回归分析用于评估预后意义。P＜0.05为差异有统计学意义。

2 结果

2.1 m5C 甲基化相关基因在肾透明细胞癌中的表达情况

为了明确肾透明细胞癌中m5C 甲基化相关基因的基本生物学功能的表达，首先从TCGA 数据库中下载并提取了有效的基因表达数据，如热图所示（图1A）。在572例肿瘤患者中，13个m5C甲基化相关基因中有12 个在肿瘤和正常组织中有差异性表达（P＜0.05）。其中包括NOP2、NSUN2、NSUN3、NSUN4、NSUN5、NSUN6、NSUN7、DNMT1、DNMT3A、TET2、DNMT3B 以及ALYREF。然而，与正常组织相比，肾透明细胞癌组织中TRDMT1 的表达没有明显差异（图1B）。在这12 个基因中，除了NSUN4 和NSUN7外，其余基因在肿瘤组织中表达较正常组织明显上调（P＜0.05）。说明m5C 甲基化相关基因在肾透明细胞癌组织和相应的正常组织中的表达有明显差异。

2.2 肾透明细胞癌中m5C 甲基化相关基因之间的相互作用

为了研究主要m5C 甲基化相关基因之间的关联，我们构建了PPI网络（图1C），以显示13个m5C甲基化相关基因之间的相互作用。TRDMT1似乎是相互作用网络的中心基因，与其余大多数m5C甲基化相关基因相关。由于PPI网络没有提供相关性的详细信息，对肾透明细胞癌进行了进一步的相关性分析（图1D）。TET2 和NSUN3 之间存在密切的相关性。除NSUN4外，NOP2与其他11个m5C甲基化相关基因相关。除NSUN7、NSUN4、NSUN3 外，NUSN2 与其他m5C 甲基化相关基因呈正相关。DNMT3A、DNMT3B、NSUN4、NSUN6亦如此。最后，TRDMT1和TET2与NSUN5负相关性最高，分别为-0.33和-0.31。

图1 m5C甲基化相关基因在肾透明细胞癌中的表达和相互作用Figure 1 Expression and interaction of m5C methylation-related genes in renal clear cell carcinoma

2.3 鉴定具有不同临床结果和特征的两组肾透明细胞癌样本

使用TCGA 数据库中572 份肾透明细胞癌样本的转录组数据进行共识聚类分析。从上述差异表达的m5C 甲基化相关基因中，使用13 个基因（NOP2、NSUN2、NSUN3、NSUN4、NSUN5、NSUN6、NSUN7、TRDMT1、DNMT1、DNMT3A、DNMT3B、TET2 和ALYREF）进行进一步研究。基于13个基因表达的相似性分析，结合共识聚类累积分布函数（CDF）和CDF曲线下区域的相对变化（图2A、B），当聚类稳定性数据集中k在2～10之间变化时，k=2被认为是合适的聚类数（图2C）。然后，根据生存参数将肾透明细胞癌样本分为两个亚群，可见亚群1中肾透明细胞癌患者的总生存期（overall survival，OS）较亚群2 缩短（P＜0.001，图2D）。此外，为了更好地预测这两个亚群肾透明细胞癌的临床病理特征，应用热图来显示分期的显著差异（P＜0.05）。而其他特征，如TNM分类、性别、年龄和临床分级等无显著差异（图2E）。从而，我们可以得出结论，13个m5C甲基化相关基因的表达特征与肾透明细胞癌患者的分期有关。

图2 共识聚类分析显示肾透明细胞癌两个亚群的预后差异Figure 2 Consensus clustering analysis revealed prognostic differences between two subgroups of clear cell renal cell carcinoma

2.4 生存预后模型的构建

为了研究m5C 甲基化相关基因在肾透明细胞癌预后中的作用，对13 个m5C 甲基化相关基因表达谱进行了单变量Cox 回归分析。结果显示，有4个基因（NSUN5、NSUN6、DNMT3B、TET2）可以用来预测肾透明细胞癌的预后（图3A）。在这4 个基因中，TET2 是HR ＜1 的保护基因，而NSUN5、NSUN6、DNMT3B 是HR ＞1 的风险基因。然后将筛选出的基因应用Lasso-Cox 回归建立生存预后模型，根据最低标准生成了单个候选基因的系数（图3B）。随后，TCGA 数据库中每个肾透明细胞癌患者的风险值计算：0.258 493×NSUN6 表达量+0.058 203×NSUN5 表达量+（-0.502 707）×TET2 表达量+0.070 293×DNMT3B表达量（图3C）。之后，根据风险值中位数界限，肾透明细胞癌患者被分成高风险组和低风险组，高风险组患者生存率显著低于低风险组（图3D）。这些结果表明，这4个m5C甲基化相关基因可以作为肾透明细胞癌的预测基因。

2.5 生存预后模型与肾透明细胞癌临床特征的相关性

热图显示，TET2 在低风险组中高表达，而NSUN5、NSUN6 和DNMT3B 在高风险组中高表达。此外，T、M 分期和分级存在显著差异（图3E，P＜0.05）。进一步运用单变量Cox 回归分析确定预后因素。结果显示：年龄（P＜0.001，HR=1.03，95%CI：1.017～1.044），分级（P＜0.001，HR=2.302，95%CI：1.858～2.851），分期（P＜0.001，HR=1.926，95%CI：1.678～2.209），T（P＜0.001，HR=1.987，95% CI：1.674～2.357），N（P＜0.001，HR=1.856，95%CI：1.528～2.255），M（P＜0.001，HR=4.388，95%CI：3.178～6.06）为单变量分析的独立预后因素（图3F）。再对这些参数进行多变量Cox回归分析，有价值的因素包括年龄（P＜0.01，HR=1.034，95%CI：1.018～1.049），分级（P＜0.001，HR=1.549，95%CI：1.218～1.968），分 期（P=0.027，HR=1.689，95%CI：1.062～2.689，图3G）。

图3 4种m5C甲基化相关基因所构建的生存预后模型Figure 3 A survival prognostic model constructed by four m5C methylation-related genes

2.6 4个模型基因在泛癌中的表达及其在肾透明细胞癌中的遗传变异

为了明确上述4个基因在人类相应肿瘤中的表达差异，在泛癌样本中对这4 个基因的转录水平进行比较，结果表明，TET2、NSUN5、NSUN6、DNMT3B在人类肿瘤中几乎均是过表达（图4A～D）。此外，在cBioportal 数据库中研究了这4 种基因与拷贝数变异之间的关系，结果显示：NSUN5、NSUN6、TET2、DNMT3B 的表达量和拷贝数变异之间的Pearson 系数分别为0.36、0.13、0.30、0.34（图5），我们认为TET2、NSUN5、DNMT3B 可能是由于其拷贝数的变异而对肾透明细胞癌的发生和进展产生影响。

图4 4个模型基因在泛癌中的表达及其在肾透明细胞癌中的遗传变异Figure 4 Pan-cancer expression of four model genes and their genetic variation in renal clear cell carcinom

图5 NSUN5（A）、NSUN6（B）、TET2（C）和DNMT3B（D）与拷贝数变异之间的关系Figure 5 Association between NSUN5（A），NSUN6（B），TET2（C）and DNMT3B（D）and copy number variation

2.7 4个模型基因对肾透明细胞癌预后产生影响的可能途径

通过limma软件包筛选出336个差异表达基因，其中281 个上调（P＜0.05，log2FC ＞1），55 个下调（P＜0.05，log2FC ＜-1），使用GOplot 和ggplot2 软件包进行了GO和KEGG富集分析，以说明差异表达基因的功能。KEGG分析结果显示：在高风险组中，这些差异表达基因的功能主要与补体和凝血级联相关（P＜0.05），而在低风险组中，则主要与癌症中的蛋白聚糖相关（P＜0.05，图6A、B）。GO分析结果显示：在高风险组中，这些差异表达基因在生物过程方面主要与体液免疫反应和水解酶活性的负调控相关（P＜0.05），在细胞成分方面主要与血液微粒相关（P＜0.05），在分子功能方面主要与酶抑制剂活性相关（P＜0.05）；而在低风险组中，这些差异表达基因在生物过程、细胞成分和分子功能方面分别与对类固醇激素的反应、含胶原蛋白的细胞外基质和碳水化合物的结合相关（P＜0.05，图6C、D）。这些结果为我们寻找m5C 甲基化相关基因在肾透明细胞癌中的潜在功能提供了新的途径。

图6 GO和KEGG富集分析结果Figure 6 GO and KEGG enrichment analysis results

2.8 4种m5C甲基化相关基因在肾透明细胞癌细胞系及组织样本中的表达

为了进一步验证4个选定的模型基因在肾透明细胞癌和组织样本中的mRNA 表达情况，本研究在肾透明细胞癌细胞系和组织样本中进行了qRTPCR实验。结果显示：与人肾小管上皮细胞系（HK-2）和癌旁正常组织相比，肾癌细胞系（786-O、A498）和肾透明细胞癌组织样本中的NSUN5、NSUN6、TET2、DNMT3B表达水平显著上调（图7）。

图7 肾癌细胞系及组织样本中4个m5C甲基化相关基因的验证Figure 7 Validation of four m5C methylation-related genes in renal cancer cell lines and tissue samples

3 讨论

本研究首先分析了这13 个m5C 甲基化相关基因在肾透明细胞癌中的表达情况，结果显示，NOP2、NSUN2、NSUN3、NSUN4、NSUN5、NSUN6、NSUN7、DNMT3A、DNMT3B、TET2和ALYREF在肾透明细胞癌组及正常组中差异表达。其中，除了NSUN4 和NSUN7 外，其余基因在肿瘤组中均高表达。接下来，根据m5C甲基化相关基因将肾透明细胞癌患者样本分为两个亚群，分析显示，亚群1中肾透明细胞癌患者总体生存率较亚群2 低，这预示着m5C 甲基化相关基因可能与肾透明细胞癌的预后相关。通过Lasso-Cox 回归分析建立了由NSUN5、NSUN6、DNMT3B和TET2 4个基因组成的生存预后模型，根据每个肾透明细胞癌患者风险值将样本分为高风险组和低风险组，观察发现，高风险组患者生存率显著低于低风险组。这一结果提示，这4 个基因可以作为肾透明细胞癌的预测因子。接着先制作热图评估了4个m5C甲基化相关基因的临床参数，然后通过单变量和多变量Cox 回归分析确立了分级和分期是肾透明细胞癌的独立预后因素。此外，GEPIA2 数据库分析发现这4 种基因在肿瘤组中的表达均上调，而在cBioportal 数据库中发现这4种基因的表达量与拷贝数变异之间有关，这一结果显示，TET2、NSUN5、DNMT3B可能是由于其拷贝数的变异而对肾透明细胞癌的发生和进展产生影响。最后，GO和KEGG富集分析验证这些m5c甲基化相关基因的主要生物学功能。总体而言，本研究结果表明m5C 甲基化相关基因对肾透明细胞癌具有相对稳定的预测能力。

目前，NSUN5、NSUN6、DNMT3B和TET2在肿瘤中的作用已被广泛研究。根据报道，NSUN5可以通过修饰真核rRNA来维持小鼠的蛋白质合成和正常的生长［13］。此外，NSUN5在结直肠癌肿瘤组织和细胞中表达上调，并通过促进细胞增殖，在体外触发细胞周期停滞进而促进肿瘤的生长［14］；NSUN6可以通过调节CDK10 参与有丝分裂纺锤体的组装进而影响胰腺癌细胞的增殖［15］，此外，NSUN6 可以通过NELFB 和RPS6KB2 相互作用来介导胶质母细胞瘤对替莫唑胺的反应［16］；另外有研究显示，DNMT3B能够与前列腺癌细胞中的转录抑制因子结合控制RAD9 的表达，进而调节前列腺癌细胞的进展［17］，DNMT3B还可以通过表观遗传学抑制miR-34a促进膀胱癌的迁移和侵袭［18］。TET2 对调节正常造血至关重要，并且可以作为肿瘤抑制因子以维持造血细胞的稳态，TET2的缺失可能会导致骨髓恶性肿瘤的发生［19］，此外，有研究表明，TET2可以抑制细胞的增殖和转移，诱导乳腺癌细胞的凋亡［20］。

总而言之，本研究构建了一个新的m5C相关基因的风险预后模型用于预测肾透明细胞癌患者风险程度和生存预后；根据该模型的风险值中位值对患者进行分组，生存分析结果显示：高风险组患者的分级和分期更高，这证明了所构建模型的可靠性。通过生物信息学技术对大规模的肾透明细胞癌患者数据进行了分析，并在细胞系和组织中进行验证，相较于以往的研究适用更广泛、可靠性更高，本研究结果对肾透明细胞癌患者的风险程度以及预后的预测具有一定意义，同时，本研究构建的模型基因可以作为肾透明细胞癌个体化治疗的治疗靶点。但是，本研究仍然存在着一些局限性，本研究数据主要来自TCGA 公共数据库，缺乏其他的数据库来源的数据，并且在使用模型之前需要进一步的多中心前瞻性临床研究进行验证。