m6A修饰赋能肿瘤免疫疗法

今天给同学们分享一篇m6A+分型+免疫的生信文章“Thyroid cancer risk prediction model using m6A RNA methylation regulators: integrated bioinformatics analysis and histological validation”，这篇文章于2023年2月15日发表在Aging (Albany NY)期刊上，影响因子为5.2。

据报道，表观遗传学重编程在甲状腺癌的进展中起着至关重要的作用。RNA甲基化占所有RNA修饰的60%以上，而N6-甲基腺苷（m6A）是高等生物中最常见的RNA修饰。本研究旨在探讨m6A调控因子构建的相关修饰模式及其对甲状腺癌临床预后和治疗效果的影响。

1. 甲状腺癌中 m6A 调节因子的遗传变异

图1A显示了23个m6A调控因子的蛋白-蛋白相互作用（PPI）网络，包括8个writer因子、13个read因子和2个eraser因子。利用 STRING 数据库，用不同的颜色来区分调控因子的功能。进一步分析（图 1B）发现了 m6A 调控因子的 CNV 突变。YTHDC2、YTHDF2、RBM15、LRPPRC、METTL14和FMR1显示CNV频率增高，相反，IGF2BP2、ZC3H13和METTL16显示CNV频率降低。图 1C 显示了染色体上 23 个 m6A 调节因子的 CNV 变化。如图 1D 突变频率分析所示，487 个样本中有 7 个样本（1.44%）的 6 个 m6A 相关基因发生了突变。图1E显示了23种m6A甲基化调节因子在甲状腺癌组织和正常组织中的表达情况。在热图中，红色表示高表达，蓝色表示低表达。据观察，23个m6A调控因子在肿瘤组织和正常组织中具有不同的表达特征。为了确定这些基因变异是否会影响甲状腺癌患者体内m6A调节因子的表达，我们研究了肿瘤和正常样本中调节因子的mRNA表达水平。我们利用 TCGA 数据库分析了甲状腺癌中 m6A 基因的差异表达。图 1F 框图中，红色代表甲状腺癌组织，蓝色代表正常组织。从图中可以看出，大多数 m6A 基因在癌组织和正常组织中有差异表达。皮尔逊相关分析显示了 23 个 m6A 基因之间的相关性。图 1G 中，红色和蓝色分别表示正相关和负相关。如图所示，VIRMA 与 YTHDF3、METTL14 和 ZC3H13 高度相关，相关系数分别为 0.92、0.89 和 0.87。YTHDC1 与 METTL14 也高度相关，相关系数为 0.89。

图1 甲状腺癌中 m6A 调节因子的遗传变异

2. 甲状腺癌 m6A 相关基因与预后的关系

如图2A所示，利用TCGA数据库进行m6A生存分析，筛选出9个与甲状腺癌患者预后相关的m6A基因。HR为风险值，HR>1表示风险增加，HR<1表示风险减少，HR.95L和HR.95h为风险波动范围。根据m6A相关基因与预后的关系，图2B显示了预后网络，反映了m6A调控因子相互作用、调控因子联系以及对甲状腺癌患者预后意义的综合情况。我们发现，不仅同一功能类别的m6A调控因子在表达上有显著相关性，而且写者、读者和擦者之间也有显著相关性。如图 2C 所示，我们利用人类蛋白质图谱分析了甲状腺癌组织中与预后相关的 m6A 分子的蛋白质表达。免疫组化图像显示了肿瘤组织中RBM15、VIRMA、YTHDC2、YTHDF2、METTL14和IGF2BP2蛋白的表达水平。

图2 &nbsp;甲状腺癌 m6A 相关基因与预后的关系

&nbsp;

3.&nbsp;分子分型

利用类别发现工具 "sense clusterplus "对甲状腺癌患者的数据进行了基于 m6A 甲基化调节因子的统一聚类。图 3A 显示，在 k = 2 到 9 之间，当 k = 2 时可以得到最稳定的聚类结果。PCA（图 3B）表明，这种分类方法能有效区分样本。在图 3C 中，热成像分析用于显示两个 m6A 分组样本中不同临床特征的分布情况。图 3D 利用 GO 富集分析绘制了直方图和气泡图，以观察主要涉及差异基因的功能类型。从图中可以看出，功能主要集中在细胞周期调控和细胞有丝分裂上。在图 3E 框图中，横轴代表 m6A 相关基因，纵轴代表基因表达水平。图 3F 显示，基于 m6A 分组的 Kaplan Meier 生存曲线在各组间无显著差异（P = 0.934）。

图3 m6A 修饰模式的确定

4.&nbsp;构建 m6A 基因亚群

尽管基于m6A调控因子表达的一致性聚类算法将甲状腺癌患者分为两种m6A修饰表型，但这些表型的潜在基因变化和表达障碍仍不清楚。因此，我们应用经验贝叶斯方法筛选了两个m6A簇之间交集的差异表达基因（DEGs），并进行了无监督一致性聚类分析，得到了两个不同的m6A基因特征子集，分别定义为基因簇A和B（图4A）。图 4B 显示了 PCA 结果。图 4C 是根据基因型绘制的热图谱，其中横轴代表样本，纵轴代表基因，蓝色和黄色代表基因组。该组具有临床特征的患者人数较多。从图中可以看出，两种成分类型有不同的表现形式。图 4D 显示了甲状腺癌患者的生存曲线。结果显示，m6A 基因簇 A 和 B 两个亚型之间存在显著差异（P = 0.023）。在图 4E 框图中，横轴代表 m6A 相关基因，纵轴代表基因表达水平。从图中可以看出，一些 m6A 分子在不同基因型之间存在表达差异。

图4 m6A 基因亚群的构建

5.&nbsp;构建 m6A 得分系统

根据预后基因采用 PCA 方法得出每个样本的 m6A 评分，并选择最佳临界值将患者分为高分和低分，绘制生存曲线。如图 5A 所示，m6A 高分和低分之间存在显著差异，表明低 m6A 患者的预后相对较好。图 5B 直方图的横轴为 m6A 高分组和低分组，纵轴为 5 年内生存状况的百分比。如图所示，低分组患者的生存率相对较高。根据 5 年生存率，甲状腺癌患者被分为生存组和死亡组。比较两组患者的 m6A 评分。方框图图 5C 显示，死亡患者的 m6A 评分高于存活患者。从图 5D 可以看出，两组 m6A 亚型之间的 m6A 评分有显著差异，其中 m6A 亚型 a 的评分最高。图 5E 是根据 m6A 聚类、基因型、m6A 评分和生存状况绘制的。该图显示了不同基因型在其他基因型中的分布情况。

图5 m6A 评分系统的构建

6.&nbsp;m6A 评分可预测免疫疗法的疗效

图 6A 显示，m6A 评分与大多数免疫细胞呈负相关。根据 m6A 评分，绘制了两个瀑布图（图 6B）。红色表示 m6A 分数高，蓝色表示 m6A 分数低。对比两幅图可以看出，高分组和低分组各基因的突变频率通常不同，其中低分组甲状腺癌肿瘤标志物 BRAF 基因的突变频率较高。图6C显示了甲状腺癌患者m6A评分与TMB之间的显著性检验，以说明TMB与m6A评分之间的关系。相关性分析表明，m6A评分与TMB之间存在明显的负相关（R = -0.15，P = 0.00095）。我们对肿瘤突变负荷进行了生存分析。图 6D 中的结果显示，低 TMB 患者的生存率明显高于高 TMB 患者。结合 TMB 和 m6A 评分，图 6E 中的结果显示，低 TMB 和低 m6A 评分患者的生存率明显高于高 TMB 和高 m6A 评分患者。为了检测 PD-L1 在 m6A 评分中的表达差异并支持相关的免疫治疗，我们检测了 PD-L1 和 CTLA-4 在不同 m6A 评分组中的表达。结果显示，低m6A评分组的PD-L1和CTLA-4的表达明显高于高评分组（图6F）。图 6H 中，横轴为 m6A 评分，纵轴为免疫治疗评分。在阳性组（POS）和阴性组（NEG）中检测到了不同的 CTLA4 和 PD-1 抗体表达水平。可以看出，低 m6A 评分组的免疫治疗得分普遍较高，表明低 m6A 评分组在免疫治疗获益方面优于高 m6A 评分组（P < 0.05）。

图6 m6A 评分可预测免疫疗法的益处

7.&nbsp;m6A 评分可预测免疫疗法的益处

图 7 显示了 m6A 评分与临床特征之间的关系。图 7A 显示，在 65 岁以上的患者中，m6A 低分患者的预后优于高分患者。图 7B 显示，在男性和女性患者中，m6A 低分患者的预后均优于高分患者。图 7C 显示，在不同肿瘤分期的患者中，m6A 低分患者的预后优于高分患者。图 7D 显示，无论 T 分期如何，m6A 得分低的患者的预后均优于得分高的患者。图 7E 显示，无论 N 分期如何，m6A 低分患者的预后均优于高分患者。图 7E 显示，在 M1 阶段，m6A 低分患者的预后优于高分患者。

图7 m6A评分的临床评估

总结

总之，本研究基于大规模队列研究证明，m6A修饰在甲状腺癌的肿瘤发生过程中起着重要作用。m6A 评分能准确预测甲状腺癌患者的预后和临床特征，为探索甲状腺癌的潜在发病机制和确定患者治疗的新靶点提供了新的见解和方向。