m 6 A signature与肿瘤免疫微环境对胶质瘤预后的预测价值

m6A Signature and Tumour Immune Microenvironment for Predicting Prognostic Value in Gliomas

m 6 A signature与肿瘤免疫微环境对胶质瘤预后的预测价值

首先让我们通过该文章的流程简单了解下其工作内容（图1），作者对TCGA-GBM / LGG数据集中的m 6 A RNA甲基化调控子进行了共识聚类。通过单样本基因集富集分析（ssGSEA）分析了癌症基因组图谱（TCGA）胶质瘤队列的肿瘤浸润免疫细胞。通过加权基因共表达网络分析（WGCNA）和最小绝对收缩和选择算子（LASSO）鉴定了四个hub基因。通过风险评分分析，Kaplan-Meier和ROC验证了诊断和预后价值。

图1.工作流程图
一.材料和方法
1.1 从TCGA中获取数据集
作者从UCSC xena下载了基于TCGA的低级别胶质瘤（LGG）和胶质母细胞瘤（GBM）队列的RNA-seq数据的基因表达数据集（665个样本）。其中的基因表达数据是FPKM值，并进行了删除批次效应处理。作者还从UCSC xena中下载了其对应的临床信息数据，包括年龄，性别，亚型等信息。
1.2选择m6A RNA甲基化调控子并使用ConsensusClusterPlus进行无监督分析
作者使用了已发表文献中12 个m 6 A RNA甲基化调控子的列表，然后，使用Gliovis系统比较了具有不同临床结果的神经胶质瘤中这些m6A RNA甲基化调控子的表达。为了研究胶质瘤中m 6 A RNA甲基化调控子的功能，作者使用“ ConsensusClusterPlus”将胶质瘤分为不同的组（50次迭代，重采样率为80％）并使用了R包“ PCA”进行了主成分分析，以研究不同胶质瘤簇中的基因表达模式。
1.3基于单样本基因集富集分析（ssGSEA）的m6A簇亚组功能分析和免疫浸润分析
使用R包“ gsva”进行了基因集变异分析（GSVA），以评估不同簇的通路富集情况。为了研究胶质瘤的免疫浸润情况，用R包“ gsva”技术，根据胶质瘤细胞特异性标记基因的表达水平，使用ssGSEA评估样本中的免疫浸润水平（记录为ssGSEA评分）。
1.4 COX回归分析
作者通过基于R包“ survival”和“ forestplot”的Cox比例风险回归分析评估了免疫细胞类型对生存时间和临床生存数据的影响。高危险比的细胞类型被认为是整体存活的危险因素。
1.5基于加权相关网络分析（WGCNA）的Hub基因与m6A RNA甲基化簇和免疫浸润相关
作者使用TCGA 中GBM / LGG数据集表达数据应用limma分析中提取了所有DEG（729个），以使用R包进行加权相关网络分析（WGCNA）。应用R包“ WGCNA”在其中找到与临床特征相关的模块和hub基因。这里面作者将邻接矩阵转换为拓扑重叠矩阵（TOM）。根据基于TOM的差异性度量，将基因分为不同的基因模块。在这里，作者将软阈值功率设置为9（无标度R2 = 0.85），切割高度设置为0.2，最小模块尺寸设置为30，以识别关键模块。hub基因被定义为基因显著性（GS）> 0.5和模块成员（MM）> 0.9的基因。
1.6 验证hub基因的预后价值
作者利用“Gliovis”在肿瘤基因组图谱（TCGA-GBM/LGG）数据集中验证了神经胶质瘤之间hub基因的表达水平及其与临床结果的相关性。为了预测具有中hub基因的胶质瘤的临床结果，作者对TCGA数据集中的5个hub基因应用了最小绝对收缩和LASSO   Cox回归算法。根据最低标准选择了四个基因来构建风险特征，并使用从LASSO算法获得的系数来计算每位患者的风险评分，如下所示：

其中n是预后基因的数目，是基因i的表达值，β是基因i在LASSO算法中的回归系数。使用中位风险评分，将脑胶质瘤患者分为高风险评分组和低风险评分组。此外，根据Rembrandt，GSE16011和CGGA数据集的三个外部队列研究了预后signature与总生存率之间的关系。使用R包survival，Kaplan-Meier评估了两组之间的生存时间差异。最后使用R包“ pROC”比较ROC曲线下的面积（AUC），使用ROC曲线来测量预后性能。
二.结果展示
2.1 m 6 A-RNA甲基化调控子的一致性聚类分析发现两组胶质瘤具有明显的免疫浸润
在这一部分，作者基于TCGA数据集中12个关键m 6 A RNA调控子的基因表达进行了共识聚类。基于m6A RNA甲基化调控子的表达相似性，聚类分析会将样本分为不同的聚类。在评估了累积分布函数，（CDF）曲线和热图下面积的相对变化之后，作者选择了一个三聚类解决方案（K = 2）。作者比较了由k = 2聚类的这两个亚组的临床病理特征，即cluster1和cluster2的亚组。结果显示m6A RNA甲基化调控子的大部分表达在两个簇亚组中表现出明显的区别和显著差异（图2a和b）。此外，cluster1亚组的总生存期（OS）明显短于cluster2亚组（图2c）。然后，作者评估了每个簇中的肿瘤突变负荷（TMB）；结果表明，cluster1的TMB明显高于cluster2（图2d）。

图2. 通过m 6 A RNA甲基化调控子识别一致性聚簇以及生存分析
此外，作者分析了cluster1和cluster2之间的DEG，并注释了它们对于生物学过程的功能基因集变异分析（GSVA）。结果表明，DEGs 富集于免疫相关的生物过程，包括IL2 / STAT5，IL6 / JAK / STAT3和干扰素-γ反应信号传导（表1）。这些发现表明，通过共识聚类确定的两个类别与胶质瘤的免疫浸润密切相关。

表1. cluster1和cluster2之间用GSVA对每个样本的通路活性评分的差异
2.2免疫景观与m 6 A RNA甲基化调控子显著相关
为了探讨免疫细胞对胶质瘤恶性进展的影响，作者分析了来自TCGA-GBM / LGG数据集的665名胶质瘤患者的RNA-seq数据，以评估免疫状况。在这里面免疫浸润的高低由欧氏距离和免疫细胞的ssGSEA得分定义。结果表明，与低免疫浸润胶质瘤相比，B细胞，T cm细胞和T辅助细胞富含高免疫浸润胶质瘤。相对而言，低浸润性神经胶质瘤的特征在于巨噬细胞，嗜酸性粒细胞，嗜中性粒细胞和aDC细胞（图3a）。为了分析m 6 A簇与免疫浸润之间的关系，作者进行了卡方检验（图3a）。然后使用Kaplan-Meier生存曲线分析来探讨免疫细胞浸润对胶质瘤患者预后的影响。观察到免疫浸润低的患者与免疫浸润高的患者相比，其总生存期较差（图3b）。作者还对TCGA数据集的免疫细胞进行单因素Cox回归分析，结果表明23/24细胞类型与OS显著相关（图3c）。

图3.胶质瘤的免疫景观
为了进一步确定m 6 A RNA甲基化调控子与免疫细胞浸润之间的关系，作者评估了免疫细胞浸润亚组与m 6 A RNA甲基化调控子表达之间的关系。结果表明，高免疫浸润与FTO，MELLT14，METTTL3，YTHDC1，YTHDC2和ZC3H13的高表达密切相关。低免疫浸润具有较高表达的ALKBH5，HNRNPC，WTAP，YTHDF1和YTHDF2（图4a（a-l））。然后作者计算了每个m 6 A RNA甲基化调控子与免疫细胞之间的关系，发现FTO，ZC3H13和YTHDC1与T cm细胞呈显著正相关，而巨噬细胞与FTO和ZC3H13呈负相关（图4b）。这些数据表明，m6A簇与免疫浸润高度相关。

图4. m 6 A RNA甲基化调控子与免疫浸润的关系
2.3 WGCNA和关键模块的标识
为了找到与胶质瘤中的m6A和免疫浸润亚型最相关的关键基因，作者对TCGA-GBM/LGG数据集进行了WGCNA分析。从该数据集中检索了胶质瘤样本信息，例如年龄，m6A簇亚组，免疫浸润亚组，OS（总生存）和OS状态（图5a）。通过软阈值功率和切割高度的设置作者确定了6个模块（图5b），从模块特征相关性的热图中发现黑色模块与临床特征相关性最高（图5c和d），尤其是免疫浸润和预后（图5d）。通过设置模块成员（MM）> 0.9和基因显著性（GS）> 0.5，作者从黑色模块中选择了5个hub基因（TAGLN2，PDPN，TIMP1，EMP3，CHI3L1）。这些基因彼此密切相关（图5e）。

图5. 通过WGCNA识别TCGA-GBM / LGG数据集中与
临床特征相关的关键模块
2.4 hub基因的表达与m 6 A RNA甲基化调控子和免疫浸润的表达相关
在这里，作者探讨了这五个hub基因的表达水平与m 6 A RNA甲基化调控子之间的关联，以阐明胶质瘤异常上调的潜在机制。相关分析表明，许多hub基因的表达与m 6 A RNA甲基化调控子显著相关。然后，我作者利用Spearman方法探索胶质瘤hub基因表达与免疫细胞浸润之间的潜在联系。结果显示，hub基因均与巨噬细胞正相关。相反，这五个基因与B细胞，Tcm细胞和Tem细胞的浸润之间负相关。这些结果表明，选定的五个hub基因与m 6 A RNA调控子和免疫浸润高度相关。
2.5 在四个数据集中验证hub基因
作者选择了上面筛选的五个hub基因，以验证其诊断和预后价值以及它们与临床特征的相关性。在预后部分，Kaplan-Meier曲线显示，这些基因的较高表达与外部数据集（GSE16011）的总体存活率差显著相关（图6a）。为了进一步预测具有hub基因的胶质瘤的临床结果，作者将最小绝对收缩和选择算子（LASSO）Cox回归算法应用于TCGA数据集中的五个hub基因（图6b）。

图6. WGCNA模块的hub基因的识别和验证并且与临床特征高度相关
最后，有四个基因与临床特征高度相关，随后将这四个基因用于建立风险特征。通过使用LASSO生成的系数，可以使用以下函数计算危险因素：

根据风险评分，将患者分为低/高风险组。使用Kaplan-Meier方法评估两组胶质瘤患者的生存时间差异，并采用对数秩检验确定两组之间观察到的差异的统计学意义。通过比较ROC曲线下的面积（AUC），使用ROC曲线来测量预后性能。结果显示，与低风险组相比，高风险组的胶质瘤患者的总生存期较短（图7a-d）。

图7. 在四个数据集中对胶质瘤患者进行风险评分分析
随时间变化的ROC曲线显示，在TCGA，Rembrandt，GSE16011和CGGA数据集中，用于4个基因标记的AUC的总体生存率分别达到0.80、0.79、0.67和0.72（图8a-d）。这些结果表明，该4-gene signature在预测胶质瘤患者的总体存活率方面有相当有效的表现。

图8. 四个数据集中的4-gene signature的Kaplan-Meier和ROC曲线
此外，与临床组织学、分级、IDH状态和年龄相比，风险评分显示OS的预后准确性更高（图9a-b）。

图9. ROC曲线比较风险评分与内部验证和外部验证队列（a，TCGA；b，CGGA）中临床组织学、分级、IDH状态和年龄的预后准确性。

转自生信人

分享