【原】易学版m6A：肾癌预后和机制分析

Effect of m6A RNA Methylation Regulators on Malignant Progression and Prognosis in Renal Clear Cell Carcinomam6A RNA甲基化调控因子对透明细胞肾癌恶性进展及预后的影响

一.研究背景

m6A RNA甲基化调控因子在ccRCC的发生和发展中起着重要的作用。本研究旨在探讨13种m6A RNA甲基化调节因子在中的作用，并确定m6A RNA甲基化调节因子在ccRCC中的risk signature和预后价值。

二.分析流程

三.结果解读

1.m6A RNA甲基化调控因子的表达与ccRCC的临床病理特征相关

作者基于TCGA-KIRC样本RNA-Seq数据和临床病理特征，筛选了13个m6A调节因子。

• 图2A-C：所示为13个m6A 调控因子在不同组织样本中的热图(A图：肿瘤状态，B图：WHO分级，C图：病理阶段)。

• 图2D - F：所示为13个m6A 调节因子在不同组织样本中的表达(D图：肿瘤状态，E图：WHO分级，F图：病理阶段)。

• 图2D显示肿瘤状态中大多数m6A甲基化调节因子的表达较为显著。且与正常组织样本相比，ccRCC组织样本中ZC3H13、METTL14、YTHDF2表达下调，FTO、ALKBH5、WTAP、METTL3、YTHDC2、KIAA1429、RBM15表达上调。

• 图2E,F显示分级越高及晚期病理期的样本中FTO、METTLE14、YTHDC1、ZC3H13、KIAA1429显著低表达。

• 用GSE15641和GSE15641数据集验证这些差异表达的调控因子。验证结果表明，FTO、WTAP、RBM15、ZC3H13在ccRCC组织中明显表达异常(图2G、H)，与之前的结果近似一致。

图2.不同临床病理特征的ccRCC中m6A RNA甲基化调控因子的表达

2.ccRCC中13m6a RNA甲基化调控因子的相互作用及相关性研究

从STRING数据库中检索到13个m6A RNA甲基化基因之间的相互关系，构建PPI网络(图3A)，并使用R中的“corrplot”包对它们之间的相关性进行分析(图3B)。结果显示：

• 除METTL3和ZC3H13、METTL3和KIAA1429、ZC3H13和WTAP在ccRCC中的表达显著相关外，其余6个reader的表达均显著相关(图3B)。

• 在PPI网络中，5个阅读器之间的相互作用很少(图3A)。此外，FTO、KIAA1429、ZC3H13、METTL14、YTHDC1的表达呈高度相关(图3B)。

图3.m6A RNA甲基化调节因子之间的相互作用

3.m6A RNA甲基化调控因子的预后价值和基于两个m6A RNA甲基化调控一致建立的Risk Signature

接下来为了研究ccRCC中这13个m6A RNA甲基化调节因子的预后价值，作者根据TCGA-KIRC中调节因子的表达水平进行了单变量Cox回归分析(图4A)。

• 显示13个调节因子中有7个预后显著(p < 0.05)。

为了确定最优的预后m6A调节因子 ，作者对7个预后相关基因进行了LASSO Cox回归分析(图4B)，每个独立预后基因的系数如图4C所示。

• LASSO结果显示两个调节因子(METTL14和METTL3)是强有力的预后因素。

    然后作KM生存分析，对这两个调节因子的预后进行了验证。结果显示：

• METTL14高表达的ccRCC患者预后良好(图4D)，而METTL3高表达的ccRCC患者预后较差(图4E)，与LASOO结果一致。

    作者基于这两个强有力的预后因素(METTL14和METTL3)，构建risk signature。然后，根据LASSO分析得到的系数构建risk score。risk score(RS)的计算公式如下:

其中n为RNA模块的数量，Coef (i)为系数，X(i)表示经LASSO分析鉴定的每个基因的z-score-transformed相对表达量。每个基因的z-score-transformed相对表达量分别乘以各自的Coef系数后加和，即为RS。

    将TCGA-KIRC患者(n = 525)按中位risk score分为高、低风险组(表1)，来测试METTL14和METTL3基因risk signature的预后作用。生存分析的结果表明,低风险组的患者的OS较好(图4F)。

采用ROC曲线来评价预测的敏感性和特异性，结果显示AUC值为0.704(图4G)，预测效果良好。

图4.2个m6A RNA甲基化调节因子的risk signature

表1.本研究包括的患者的临床病理特征

4.risk score预后与ccRCC的临床特征密切相关

图5A：热图显示两种m6A RNA甲基化调控因子在高风险组和低风险组中的表达水平。结果显示，高风险组与低风险组中WHO分级(p < 0.001)、病理分期(p < 0.001)、T分期(p < 0.001)、M分期(p < 0.01)统计学上显著性差异。

还探讨了risk score与临床病理特征之间的关系：

• 图5B，C：显示按WHO分级和病理分期划分的患者的risk score存在显著差异。

• 图5D，E：ROC曲线显示，在预测5年生存率方面，risk score（AUC = 0.708） 预后优于WHO分级 (AUC = 0.655)；而与病理分期相近(AUC = 0.709)。

然后进行了单变量和多变量Cox回归分析：

• 在单变量分析中，年龄、WHO分级、病理分期、TNM分期、risk score等均与预后相关(图5G)，

• 而在多变量Cox回归分析中，显示只有risk score、诊断年龄与OS是独立的预后因素(p < 0.05)(图5H)。

作KM生存分析，分析不同分级及病理分期的risk signature对预后的影响。结果显示，低risk score组的具有更好的OS（图5I-J)。这些结果表明，risk score可以作为ccRCC的独立预后因素。

图5.risk score与临床病理特征的关系，以及按WHO分级和病理分期的患者risk signature的预后价值

5.两组独立预后m6A调节因子的一致聚类鉴定了两组具有不同临床结果的ccRCC

基于2个m6A 的METTL14和METTL3两个强有力的独立预后因子，作k = 2 - 9逐渐增加的聚类分析，发现k = 4为最优选择(图6A,B)。但当k = 2时，各子组之间的干扰最小。

所以作者使用k = 2进行一致性聚类分析，并确定了两个cluster：cluster1和cluster2。KM生存分析发现cluster1组患者的OS明显短于cluster2组(图6C)。

此外，使用主成分分析(PCA)来比较cluster1和cluster2的转录谱。结果表明，两个cluster之间存在显著差异(图6D)。

接着利用KEGG (图6E)和GO(图6F)对cluster1中不同表达的基因进行通路分析和功能注释。

• KEEG分析的结果表明,差异表达基因主要富集在恶性肿瘤相关的途径(图6e)。

• GO分析结果与癌症相关的生物学过程相关，包括“白细胞迁移”、“白细胞激活调节”、“免疫效应调节”等(图6F)。

• GSEA显示cluster1中高表达基因的恶性肿瘤特征显著富集(图7)。

• 以上结果表明，根据METTL14和METTL3识别出的两个cluster与ccRCC的恶性程度密切相关。

图6.cluster1和cluster2的ccRCC患者的总生存期

图7.对Cluster1的GEEA分析

            最后我们总结一下本文献的分析套路：首先根据TCGA-KIRC，分析m6A的表达情况，然后筛选关键的独立预后因子。然后建立risk signature，将ccRCC患者分为高危组和低危组。接着根据关键的独立预后因子的表达，通过一致性聚类方法。验证聚类后的结果与临床病理特征的关系。也对risk score与临床病理特征关系相关性验证。最后作GO、KEGG分析，验证与ccRCC癌症相关通路、关键的生物学过程和特征相关。