6+的细胞焦亡+预后+免疫微环境新思路

The Pyroptosis-Related Signature Predicts Prognosis and Indicates Immune Microenvironment Infiltration in Gastric Cancer

胃癌焦亡相关特征预测预后并提示免疫微环境浸润

胃癌(Gastric cancer, GC)是癌症相关死亡的主要原因之一，具有高度的异质性。如果我们要改进治疗策略，开发一个特定的预后模型是重要的。细胞焦亡可因幽门螺杆菌（一种主要致癌物）以及化疗药物而发生。然而，GC对细胞焦亡的预后评估尚不充分。采用细胞焦亡相关调节因子的一致性聚类方法对来自4个GEO集的618例GC患者进行分类。对差异表达基因进行Cox回归后，采用LASSO-Cox分析建立我们的预后模型(PS-score)。采用TCGA-STAD队列作为验证集。ESTIMATE、CIBERSORTx和EPIC用于研究肿瘤微环境(TME)。采用阻断PD1/PD-L1的免疫治疗队列来研究治疗反应。基于焦亡相关调节因子的GC亚型可以根据不同的临床特征和TME对患者进行分类。本研究中确定的两种亚型之间的差异被用来建立我们命名为“PS-score”的预后模型。PS评分可预测胃癌患者的预后及总生存时间。PS评分低意味着炎症细胞浸润更大，PD1/PD-L1阻滞剂的免疫治疗反应更好。我们的研究结果为今后针对焦亡及其免疫微环境改善预后和免疫治疗应答的研究提供了基础。

1.胃癌数据集和预处理的来源

工作流程图(补充图1)描述了在统计分析的每个阶段使用了哪些样本。来自Affymetrix ®的微阵列数据来自基因表达综合数据库(Gene Expression Omnibus，GEO)。通过使用“SVA”包的“ComBat”算法校正非生物技术偏差的批量影响。这项研究中使用的所有临床信息都可以在GEO数据库中公开获得。对于TCGA中的数据集，基因表达的RNA测序数据(FPKM值)从UCSC下载。没有生存信息的患者被排除在进一步分析之外。应用“maftools”包中的瀑布函数来呈现突变景观图。免疫治疗数据集中的患者分别来自GSE78220、NCT01358721和IMURITY 210，分别接受派姆单抗、纳武单抗和阿特珠单抗治疗。

2.一致性聚类方法

通过k均值方法，应用一致性聚类来识别与焦亡调节因子表达相关的不同焦亡相关模式，聚类的数量和它们的稳定性是由使用“ConsensuClusterPlus”包的一致性聚类算法确定的，为了保证分类的稳定性，我们进行了1000次重复。

3.基因集变异分析(GSVA)

我们用“GSVA”R软件包在热图上进行了GSVA富集分析。我们从MSigDB数据库下载了“c2.cp.kegg.v6.2.symbols”来进行GSVA分析。经“LIMMA”软件包调整后P<0.05，表明不同亚组间差异有统计学意义。

4.差异表达基因(DEG)

我们使用“LIMMA”软件包的经验贝叶斯方法获得DEG。选择 DEG 的显着性标准设置为调整后的P < 0·05 和 Log2 FC 的绝对值 ≥ 0·8。

5.TME细胞浸润

我们用CIBERSORTx算法和EPIC对免疫细胞比例进行定量。对于CIBERSORTx，我们使用LM22 signature并置换1000次将标准化的基因表达数据上传到门户网站。EPIC是一个基于网络的分析和发现平台，用于以标准化方式分析来自免疫细胞的大量细胞质谱数据。肿瘤纯度评分由“ESTIMATE”包估计。

6.多变量Cox回归分析

我们使用单变量Cox回归分析和多变量Cox回归分析对前面描述的四个GEO数据集的总生存期(OS)进行了分析。以误发率(FDR)<0.05作为统计界限。应用“Forestlot”软件包对焦亡相关亚组进行多因素预后分析。

7.PS评分模型的建立及预后分析

我们用LASSO−COX分析建立了一个有效的预测模型。然后使用OS从训练集中获得最有用的预测特征，其中k、βi、Pi分别表示特征基因数、系数指数和基因表达水平。截止点是使用“survminer”软件包确定的。我们使用Kaplan-Meier生存曲线来识别模型区分不同亚型患者的能力，使用时间依赖的受试者工作特征曲线(ROC曲线)来确定模型的效率。C指数由“survcomp”软件包计算，并使用“cindex.comp”软件包进行比较。

8.RNA提取和Real-Time PCR

用Trizol试剂提取总RNA。随后，提取的RNA用带有gDNA Eraser的PrimeScript RT试剂盒逆转录。对 cDNA 进行基于SYBR Green的实时PCR分析。

9.统计分析与分界值

采用Spearman‘s相关分析和距离相关分析计算相关系数。采用对数秩检验确定生存曲线差异的显著性。本文提到的截止值是1258，这是“Survminer”包的最佳截止值。ROC曲线、时间相关ROC曲线和曲线下面(AUC)分别使用“pROC”和“timeROC”包导出。使用“iauc.comp”包进行曲线下积分面积 (IAUC) 的比较。

“RCircos”包使我们能够绘制23对染色体中细胞焦亡调节因子的拷贝数变异图。使用Wilcox检验进行两组之间的分析。Kruskal Wallis 检验也用于比较三个或更多组。基因表达数据和所有统计分析均在 R 4.0.0、GraphPad Prism 8 和 SPSS 26 软件中进行。在适当的情况下，通过卡方检验或 Fisher 精确检验比较临床特征。所有统计P值都是双面的，P< 0.05 表示统计显着性。

我们使用STRING平台分析了与 11 个细胞焦亡相关调节剂相关的潜在生物功能网络（图1A）。这些调节剂主要集中在免疫反应和焦亡的调节上。在基因水平上，433个样本中有59个(约13.63%)显示出与焦亡的调节子突变。其中，CASP5的突变频率最高。我们没有在任何GC样本中发现任何GSDME突变(图1B：TCGA-STAD队列中433例胃癌患者的突变图谱)。我们还在GSE62717的 11 个焦亡相关调节因子中的7个中发现了CNV；这些是常见的变化大多数都集中在拷贝数扩增上（图1C）。我们确定了染色体上具有CNV特征的7个调节子的变化(图1D)。在表达水平上，这11个调节子能够帮助我们区分胃癌患者的正常样本和肿瘤样本（图1E）。与正常样品相比，除GSDME和GZMA外，其余调节剂在GC样品中均有增加（图1F）。为了更有效地验证我们的发现，我们测试了从我院采集的 22对肿瘤组织和正常邻近组织样本中11个细胞焦亡相关基因的 mRNA 水平。同样的发现是这些调节因子在肿瘤中被上调（图1G）。

图1 GC中焦亡相关调节因子的特征和差异

我们使用来自同一平台的四个GEO数据集以及OS数据和临床信息创建了一个队列。基于11个与焦亡相关的调控因子的表达水平，我们用无监督聚类的方法确定了两种不同的调控模式，包括与焦亡相关的集群 1中的267例和与焦亡相关的集群2中的351例（图2A）。集群1的生存优势高于集群2(图2B)。为了探索这两种模式在生物学行为上的差异，我们进行了GSVA富集分析(图2C)，集群1 在凋亡途径中也显着富集。集群2 显示出致癌激活通路的富集。随后，我们还证实，通过应用11个与焦亡相关的调节剂，可以区分这两种调节模式。（图 2D）

图2 与细胞焦亡相关调节因子相关的胃癌亚组

接下来，我们分析了细胞浸润数据，发现激活的先天免疫细胞浸润在集群1 中很丰富，包括树突状细胞、M1巨噬细胞、B细胞以及激活的CD4 和CD8T细胞的存在，从而赋予显着的存活率优势。集群2富含内皮细胞、肥大细胞、M2 巨噬细胞和静息T4记忆细胞。（图 3A、B）为了进一步探索这两个不同集群的临床表现，我们重点关注ACRG队列，该队列代表了与我们的 618 名患者相关的临床信息的最全面研究。（图 3C）。我们对 ACRG 队列中具有复发数据的 282 名患者进行了统计分析。集群2模式与更多的复发病例相关（图 3D）。这些发现表明，不同的细胞焦亡相关模式代表不同的 GC 特征并具有不同的 TME 状态。

图3 不同的焦亡相关集群表现出不同的临床特征和TME细胞浸润

为了更好地将这些亚型应用于GC的临床治疗并确定每位患者的特定评分，我们接下来探索了两种模式之间的差异并确定了特定的基因特征。我们还量化了基因特征，以便将其应用于每位患者的诊断和治疗。首先，我们确定了113个DEG，其绝对值 Log2 FC < 0.8 和p< 0.05 与两种不同的调控模式相关（图 4A）。接下来，单变量和多变量 Cox 回归分析确定了 22 个可用作独立预后特征的基因。（图 4B）。10个基因被确定为与免疫反应最相关的基因，出现在第1类中，而12个编码癌症发生蛋白的基因倾向于更普遍地出现在第2类模式中（图 4C）。为了建立一个能够量化每个患者的模型，通过应用 LASSO-Cox 回归模型以最小的λ保留了22个DEG中的六个。我们使用这些来构建一个与焦亡相关的特征分数，我们将其命名为“PS-score”（图 4D）。接下来，我们试图通过预测患者的预后来进一步确定PS-score的价值。我们将患者分为高和低PS分数组，最佳截止值为1.258。我们发现低组比高组有明显的生存优势，为了证明 PS-score 的普遍指示值，我们还在更多的队列中验证了这个分数并得到了相同的结果（图4E：来自GSE62254队列的约300名胃癌患者的不同PS评分亚组的OS曲线（截断值为 1.258，图F：来自四个GEO队列（GSE15459、GSE34942、GSE57303和GSE62254）的618个胃癌样本中截断值为1.258的不同PS评分亚组的OS曲线））。我们进一步证明PS评分是GC患者3年生存率和5年生存率的良好指标。（图 4G）。PS-score 在预测GC患者预后方面优于其他模型（图4H）。

图4 生成基因表达特征以预测基于细胞焦亡相关簇的患者存活率

我们还发现，TCGA-STAD 队列中的高分组以相同的临界值预测了较差的生存特征（图 5A）为了验证PS-score对所有消化道癌症的普遍适用性，我们收集了包括胃癌在内的7 个消化道癌症样本的数据，并将这些数据分为具有相同截断值的两组。我们惊奇地发现低分组也表现出生存优势（图5B）为了确定具有不同临床表现的患者的 PS 评分的特异性，我们分析了 ACRG 队列中临床特征与PS评分之间的关系。我们发现，生存优势明显的cluster 1模式，得分明显偏低；这与我们之前的研究一致（图5C）。得分高的患者也经历了更多的复发（图5D），我们惊奇地发现，随着患者生存期的增加，PS-score呈逐渐升高的趋势，由此可见晚期胃癌患者的评分高于早期胃癌患者（图5E）。对治疗反应不佳的EMT亚型也表现出高分（图5F）。PS-score 低表明对后续治疗反应良好，因此表明PS-score可能代表GC治疗的良好标志物（图5G）。除了不同的临床表型外，患者体内不同焦亡相关调节因子的表达也会影响PS-score的特异性。我们发现，无论我们使用ACRG队列还是 TCGA-STAD队列，除了GSDME 外，其他调节因子的表达水平在高分组中均显着降低（图 5H）。然而，高分组主要富含与癌症发展相关的信号通路（图 5I）。这些发现表明，PS-score 在与不同信号通路激活模式相关的不同患者中发挥特异性，可用于评估GC患者的某些临床特征和治疗效果。

图5 PS-score 评分模型的特征

鉴于PS评分在预测GC患者预后方面的重要性，我们接下来尝试探索其临床应用价值。我们构建了一个具有7个临床特征的列线图，这些特征很容易获得，并且普遍认为对GC的预后和PS评分预测GC患者3、5和8年生存率的能力有一定影响。（图 6A）C-index 为0·764表明列线图具有良好的预测值（图 6B）。当我们使用 PS-score 作为单独的指标来区分 GC 患者时，我们发现 ACRG 队列的 300 名患者被分为高组和低组（图 6C）。冲积图用于可视化患者特征的变化（图 6D）。接下来，我们研究了 PS-score 是否对临床治疗，尤其是免疫治疗具有指导价值。我们分析了不同 PS-score 的 TME 细胞浸润情况（图 6E-G）。幸运的是，我们探索了不同 PS 分数的 TMB，发现较低组的 TMB 较高，这表明免疫治疗反应更好（图6H）。

图6 PS-score评分模型与不同PS-score亚组TME的临床应用价值

我们比较了 ACRG 队列和 TCGA-STAD 队列中免疫检查点基因表达水平的差异（图7A）。为了更好地说明PS评分对免疫治疗的预测能力，我们在ACRG队列中应用了TIDE。我们惊喜地发现TIDE和 PS-score 之间呈正相关（图7B）。此外，预测的反应表明PS评分可能是GC免疫治疗的良好预测指标（图7C）。我们还证实，较高的PS分数代表了乳腺癌和肾癌等“热肿瘤”的预后较差（图 7D、E）。我们检测了不同PS-score呈现的免疫检查点表达的变化，也发现低分组有更高的表达（图7G)，因此表明PS-score可能是一个重要的指示性指标，在其他癌症中发挥同样重要的作用。正如我们假设的那样，我们发现对免疫疗法有反应的患者也表现出较低的PS评分（图7H）。

图7 PS-score 评分模型在 PD1/PD-L1 免疫治疗中的强大作用

我们的分析最终形成了几个结论：（1）大多数与细胞焦亡相关的调节因子在集群1中表现出高表达水平；(2)集群2是一个单独的亚型，预后较差；(3)第1集群被确定为免疫炎症表型，第2集群被确定为免疫排除表型；(4)综合考虑临床信息和RNA数据更可能反映细胞表型。我们的研究旨在将GC患者分为亚型，识别DEG并建立预后模型，并将细胞焦亡与患者预后联系起来。尽管我们进行了多角度、多数据库的验证，但本研究仍存在需要考虑的局限性。创建的模型在接受免疫治疗的患者的所有时间生存阶段方面没有很好的预测价值。这可能是由不同癌症的特异性引起的，需要更广泛的研究。肿瘤异质性确实是一个不容忽视的问题，随着细胞焦亡相关研究的开展，可能会针对不同类型的肿瘤提出更有针对性的改进方案。