最新见刊的坏死性凋亡多组学生信top文章

今天分享一篇今年7月发表在中科院二区top期刊Frontiers in Immunology上的一篇临床hallmark纯生信文章“Necroptosis identifies novel molecular phenotypes and influences tumor immune microenvironment of lung adenocarcinoma”。

发文单位：1. 武汉大学人民医院肿瘤科；2. 华中师范大学数学与统计学学院；3. 武汉大学医学院；4. 山东大学齐鲁医学院公共卫生学院。

第一作者：赵晨、熊可为；通讯作者：赵晨、李祥攀

流程图

摘要

本研究旨在研究肺腺癌（LUAD）中坏死性凋亡的免疫和表观遗传突变情况，识别新的分子表型，并开发基于坏死调节分子的预后评分系统，以更好地了解LUAD的肿瘤免疫微环境（TIME）。基于TCGA和GEO数据库，以及29个坏死性凋亡相关基因，利用无监督的共识聚类将患者分为不同的坏死表型。我们将这些表型与临床特征、免疫细胞浸润水平和表观遗传突变特征系统地联系起来。然后构建了一个新的评分系统，称为NecroScore，通过主成分分析对LUAD的坏死性凋亡水平进行量化，能够有效地预测生存和药物治疗。总之，坏死相关分子与LUAD的基因组多样性相关，在形成LUAD的TIME中发挥了重要作用。坏死性凋亡表型可以区分不同的TIME和分子特征，NecroScore是预测LUAD预后、免疫和化疗获益的一个有前景的生物标志物。

结果

表观遗传和转录特征

在本研究中我们通过R软件包 "DESeq2 "在TCGA-LUAD队列中确定了30个差异表达的坏死相关基因（补充表1）。这30个分子是 ALK, AXL, BACH2, BNIP3, CDKN2A, CFLAR, DIABLO, DNMT1, FADD, GATA3, ID1, IDH1, IDH2, ITPK1, KLF9, LEF1, MPG, MYCN, OTULIN, PANX1, PLK1, SLC39A7, TERT, TLR3, TNFRSF1B, TNFRSF21, TNFSF10, TRAF2, TRIM11, and ZBP1。补充图1的左上部分展示了坏死性调控的机制。为了揭示LUAD中坏死调节器的表观基因组学特征，我们分析了坏死调节分子的甲基化水平。结果表明TERT和GATA3基因的启动子甲基化和转录之间有显著的负相关关系，而BACH2的表达与其甲基化有弱的正相关（图1A）。由于CNV可以在arm-level水平影响整个基因组，或在focal-level改变一个广泛的区域，对体细胞突变产生更大的影响，我们进一步描绘了这些基因的CNV景观。在LUAD中，结果显示了30个基因的拷贝数异常，其中TERT和OTULIN具有广泛的扩增，而CDKN2A则表现为缺失（图1B）。在泛癌中，根据cBioPortal的突变分析，CDKN2A显示出最高的伴deletion突变频率，其次是TERT，TRIM和ZBP，其频率为11%。在561个LUAD样本中，有111个出现了与坏死相关的基因突变，其频率为19.79%。ALK的突变频率最高，其次是CDKN2A，而剪接位点只在CDKN2A和IDH1中被观察到（图1C）。
研究发现，大多数调节因子在癌旁和肿瘤的表达上有显著的差异。有几个基因在肿瘤组织中的表达水平较低，如AXL, AXL, AXL, AXL, AXL等；而其他基因则过度表达，如BNIP3，DIABLO，FADD和IDH1（补充图2A）。此外，BNIP3, DIABLO, DNMT1, FADD, LEF1, PANX1，PLK1，TNFRSF21和ZBP在不同的病理阶段（从I到IV），表达量有显著变化（补充图2B）。对这30个基因进行了层次聚类，得到了四个亚型。Cluster A（TERT、MYCN和ALK），Cluster B（AXL, BNIP3, CFLAR, DIABLO, DNMT1, FADD, ID1, ITPK1, KLF9, MPG, PANX1, TLR3, TNFRSF1B, TRAF2, TRIM11），以及Cluster D （BACH2, CDKN2A, GATA3, LEF1, PLK1和ZBP1）。为了充分了解LUAD中坏死相关基因的调控机制，我们收集了1,064个样本，使用 "ComBat "方法修正了批次效应（图2A），可以看出成分1和2的总和从90.9%下降到14.9%，表明批次效应被大致消除。我们用网络图全面地描绘了坏死性分子与遗传的相互作用情况的相似性及其对LUAD患者的预后价值的综合情况（补充图2C）。我们发现同一类别的坏死调节器在转录表达水平上呈现出明显的正相关关系，并且在四种修饰模式中表现出显著的相关性。根据log-rank检验，过度表达的DIABLO，FADD，ID1，MPG，TRIM11，PLK1，TNFRSF21和TERT与不良的总生存结果有关，即风险因素，而ALK被认为是一个保护因素。单变量Cox模型显示，FADD，PANX1，PLK1，TRIM11，TRAF2，TERT和ID1可能是风险因素；而ALK对预后有预防作用（补充图2D）。上述结果强烈表明整个坏死调节基因的不平衡与LUAD的发生有关。

图1. 坏死性凋亡相关基因在LUAD中的表观基因组变异特征. (A). 棒状图显示了DNA启动子甲基化和mRNA表达水平之间的相关性。圆圈大小代表相关系数。较蓝的点表示较低的P值。 (B). TCGA队列中坏死调节器的CNV频率。条形图的高度代表改变的频率。粉红色代表拷贝数扩增，绿色点代表拷贝数缺失。 (C). 左边部分显示了来自cBioPortal的28种癌症的突变频率和类型。每一列表示个别病人。上面的带子显示了癌症类型和生命状态。右边的数字显示每个基因的突变频率。右边部分说明了TCGA中肺腺癌的体细胞基因突变情况。上面的柱状图显示每个样本的总体突变频率。右边的叠加柱状图显示了每种变异类型的比例。瀑布图中只显示了有突变的分子，即突变频率>0%的分子。

无监督聚类识别LUAD的亚型

基于K-means的无监督学习方法被应用于 对29名患者进行分类（OTULIN在批次效应批量效应校正后的联合表达矩阵中没有检测到，因此被排除在外）的坏死调节分子进行分类。根据不同聚类的分布和累积分布函数，最终确定了最佳聚类和三种不同的坏死表型。包括Cluster A的210例，Cluster B的578例以及Cluster C的276例（补充图3A-C；补充表2）。PCA显示三种聚类亚型可以通过29种坏死相关分子的表达水平来区分（图2B）。热图直观地显示了坏死吞噬作用调节因子的表达变化。坏死性凋亡相关基因在Cluster C中表现出最高的表达水平，其次是Cluster B，然后是Cluster A（图2C）。通过Fisher's exact检验，分析了三个坏死性肿瘤群的临床病理表型的分布，包括年龄、性别、吸烟史、复发情况和病理阶段。我们发现，Cluster B和C有更多的女性患者（P < 0.001）。Cluster A的病人更有可能吸烟（P < 0.001）。早期阶段（I期或II期）的患者大多来自Cluster B（P < 0.001）。然而，患者的年龄在各表型之间的分布没有统计学意义。三种主要坏死表型的KM生存分析显示，Cluster B的预后优势相对突出（图2D）。单变量和多变量的Cox 比例危险回归模型表明，Cluster B与C可以预测LUAD的总生存期（单变量，HR = 0.616，95%CI = 0.485-0.782, p < 0.001；多变量，HR = 0.674, 95%CI = 0.529-0.860, p = 0.001；图2E）。此外，年龄和分期也是预测生存的基本因素。

图2. 平均联动的K-means无监督聚类识别三种坏死性凋亡表型。(A). 主成分分析（PCA）显示了29个坏死相关基因在肺腺癌队列中转录表达的分布（左部分）和批量效应校正后（右部分）。 (B). PCA表明29个坏死分子的表达可以区分TCGA和GEO-meta队列中的坏死性凋亡表型。 (C). 热图显示了三种坏死表型和临床病理特征之间的相关性，以及29个坏死相关基因的表达变化。星号代表P值（***p＜0.001，*p＜0.05，ns是无显著性的缩写），颜色代表统计学检验方法：蓝色表示Fisher精确检验，黑色表示Kruskal-Wallis检验，绿色表示Welch单因素方差分析检验。每个细胞的内部颜色显示了表达的变化。较红的与过度表达的水平有关，较蓝的与较低的水平有关。 (D). Kaplan-Meier曲线和log-rank检验显示，在TCGA和GEO-meta队列中，三种坏死性凋亡亚型之间的生存率差异显著。 (E). 单变量和多变量的Cox模型显示了坏死性凋亡聚类亚型对总生存的预后影响的结果。

肿瘤免疫微环境在三种坏死性凋亡表型中的变化

由于癌症免疫在肿瘤进展和细胞增殖中起着至关重要的作用，我们推测不同集群亚型的肿瘤免疫微环境可能是与众不同的。我们首先探索了坏死表型与多种免疫调节剂和免疫检查点之间的相关性。结果显示，Cluster B和C的趋化因子的表达水平明显上调，如CCL11、CCR1、CXCL10和PPBP（图3A，上部）。白细胞介素、干扰素、受体和其他一些细胞因子的表达也在Cluster B和C中出现，而Cluster A的免疫分子水平强烈下调（图3A，中间和底部部分）。与Cluster A和B相比，Cluster C表现出7种免疫分子的最高表达，这些分子是免疫治疗的潜在靶标，包括CD274（PD-L1）、PDCD1（PD-1）、CD247、PDCD1LG2（PD-L2）、CTLA-4、TNFRSF9和TNFRSF4（图3B）。
我们进一步分析了20种免疫细胞类型的浸润水平和免疫循环中每个步骤的富集评分。对基因特征的探索表明，免疫细胞的浸润在Cluster B中明显较高，其次是Cluster C，而在Cluster A中发挥的部分明显较低，如D4+ naïve T cells, follicular helper T cells和M1 macrophages（图3B）。值得注意的是，活化的记忆型CD4+T细胞在Cluster B中的丰度最低，而在Cluster A和C中的丰度较高。我们用ssGSEA算法计算了七个抗肿瘤免疫循环步骤的富集得分，所有的步骤在Cluster B和C中显示出较高的分数，但在Cluster A中显示出特别低的水平（图3B）。此外，Cluster A的MeTIL也很弱，这表明肿瘤浸润性免疫细胞的比例很低（图3B）。
我们还采用了TIDE算法和ESTIMATE方法来比较TME免疫亚型中免疫和基质成分的异质性，可以量化坏死部分和潜在的免疫治疗效果之间的相关性。研究发现，Cluster B的TIDE和功能障碍得分最高，其次是Cluster C，然后是Cluster A，而在排除分析中发现了一个相反的结果（P < 0.001，图3B上部；补充图4A-C）。聚类B和C的免疫和基质得分明显高于聚类A，表明聚类A显示出最高的肿瘤纯度（p < 0.001，图3B上部；补充图4D-F）。此外，我们分析了几个免疫途径的生物标志物。结果显示，功能的归一化富集分数明显较高，主要位于C簇，包括co-inhibition of APC，checkpoint inflammation-promoting，MHC class I，para-inflammation，coinhibition stimulation of T cells和type I IFN response（所有P < 0.001，图3C）。聚类B显示costimulation of APC，CCR，cytolytic activity，HLA和type II IFN response的得分相对较高，而聚类A的丰度最低（所有P < 0.001，图3C）。利用相关性策略，包括ALK，AXL，KLF9，LEF1，TNFRSF1B和ZBP1在内的几个坏死调节基因，与大多数的免疫细胞呈正相关，而负相关则表现在BNIP3、DIABLO、FADD、IDH1、IDH2、PLK1和SLC39A7等分子（图3D）。

图3. 三种坏死表型中的免疫调节剂和几种定量的肿瘤免疫微环境的改变。 (A). 热图显示了三种坏死吞噬亚型中趋化因子、白细胞介素、改变因子及其相应受体和其他细胞因子的mRNA表达水平的变化。星号代表P值（***p＜0.001，**p＜0.01，*p＜0.05，ns是无显著性的缩写），颜色代表统计学检验的方法。黑色表示Kruskal-Wallis试验，红色表示单因素方差分析试验，绿色表示Welch单因素方差分析试验。 (B). 综合热图显示了TCGA-LUAD队列中的肿瘤免疫微环境的变化。上面部分显示了几个 免疫检查点的变化，中间部分显示24种免疫细胞类型的富集分数，底部显示抗肿瘤免疫的七个步骤的富集水平。抗肿瘤免疫循环中的七个步骤。MeTIL的DNA甲基化水平展示在最下方。TIDE评分、功能障碍分数、排斥分数、肿瘤纯度、免疫分数和基质分数被注释在热图的顶部。 (C). 通过Kruskal-Wallis检验展示13个免疫途径的归一化分数的变化。 (D). 热图显示了29个坏死调节分子和免疫细胞之间的联系。免疫细胞之间的关系。红色表示正相关关系，蓝色表示负相关关系。星号表示P值（**p < 0.01，*p < 0.05）。空白的方格表示相关性没有统计学意义。

基因组在三个坏死性凋亡表型间的改变

我们进一步探讨了已确定的LUAD分子表型的基因组差异，包括突变特征，SNV，以及转录因子调控网络的活性。结果发现，Cluster A和ClusterC的TMB明显高于簇B（p < 0.001，图4A，B）。我们推断，即SBS1（与年龄有关）、SBS2（APOBEC活动相关的），SBS4（吸烟相关的）和SBS5（与ERCC2突变有关）。然而，观察到的统计学意义仅在SBS1中显示，Cluster B在年龄相关的特征中具有更多的突变（p < 0.001, 图4A；补充图5A)；与C组相比，A组在吸烟相关的特征中显示出更高的突变频率（图4A；补充图5B）。因此。这些突变特征可以作为解释
每个分子坏死性表型的扰动。在前20个经常突变的基因中，我们发现C组在TP53，TTN，MUC16的突变中更为丰富。CSMD3和RYR2的突变，而Cluster B显示出相对较低的突变频率。与ClusterA和Cluster C相比，Cluster B的突变频率较低（图4A）。拷贝数扩增和缺失在Cluster B和Cluster C中下降，但在群组A中显著增加。Cluster A在病灶和手臂水平上都有显著增加，而Cluster B的CNV负担明显降低（图4C）。我们探讨了不同染色体的CNV和GISTIC评分分布，验证我们之前的结果。在不同的染色体上，CNV和GISTIC分数分布的是最大的频率差异位于8号和14号染色体上（补充图6A，B）。
为了进一步分析集群之间的转录组差异，我们对总共8个LUAD特有的转录因子和15个与癌症染色质重塑相关的调节因子进行了分析。我们确定了三个坏死性凋亡的生物学差异，因为调控子的活动与这些集群密切相关（图4D上半部分）。Cluster A和C共享类似的regulon模式，而Cluster A有较高的FOXA2和ATOH8的较高活性，而Cluster C则与激活的HOXA4明显相关。Cluster B包含高活性的HOXA4，DACH1，EPAS1，ETV5，FOXA2，ATOH8和SMAD6调节子。与癌症染色质重塑相关的Regulon活性揭示了在三种坏死性疾病中的不同调节模式，表明表观遗传学驱动的转录网络可能是这些分析表性的重要区分因素（图4D下半部分）。我们发现，Cluster A和C也有类似的调控模式，而Cluster A表现出更高的SIRT1、SIRT4和SIRT5活性。同时，SIRT1，CLOCK，EP300，HDAC5，SIRT4、SIRT2、HDAC10和HDAC8的regulon活性在Cluster B中显著地富集。

图4. 三种坏死凋亡表型之间的基因组特征的区别。 (A). 根据坏死性凋亡聚类的体细胞突变情况亚型。从上到下分别代表肿瘤突变负担（TMB），四个突变特征（SBS1、SBS2、SBS4和SBS5）的相对贡献，选定的顶级突变基因。和最广泛的拷贝数改变（q值<0.05），以及选定的分子。变异转换的比例显示在右边的堆积条形图中。 (B). Kruskal-Wallis检验揭示了三个集群之间的TMB差异。 (C). Kruskal-Wallis检验揭示了坏死凋亡表型中的focal和arm-level拷贝数变化。星号表示P值（***p<0.001，**p<0.01和 ns是无显著性的缩写）。 (D). 热图显示了8个转录因子（上半部分）和与染色质重塑有关的潜在调节因子的调节活性与染色质重塑相关的潜在调节因子（下半部分）。

信号通路的功能异常

我们探讨了坏死的三种表型之间的基因表达差异，并分析了几种与癌症有关的信号通路的特点，这有助于了解潜在的潜在调节机制和可药用的通路。激活的通路如apoptosis, cell cycle A, EMT-A, and hormone AR-A与坏死性凋亡调节分子相关，但在泛癌中存在着抑制状态的cell cycle I, hormone AR/ER-I, and RTK-I等功能（图5A）。在LUAD中，大多数坏死性凋亡相关基因激活了cell cycle, EMT, and apoptosis等通路，而它们抑制了激素ER和RTK信号（图5B）。我们使用语义相似性来聚合紧密相关的类别，在30个调节分子中筛选出在坏死性凋亡中起作用的枢纽基因。结果表明，TNFRSF21的权重最高，其次是DIABLO和TRAF2（图5C）。为了分析这些不同的坏死表型的生物学行为，我们进行了GSVA富集的探索。如图5D所示，Cluster A在RNA polymerase，aminoacyl tRNA biosynthesis，homologous recombination和DNA replication中明显富集。Cluster B包括与viral myocarditis，cell adhesion molecules (CAMs)，natural killer cell-mediated cytotoxicity和hematopoietic cell lineage相关的生物途径。此外，Cluster C在primary immunodeficiency，the intestinal immune network for IGA production，cell cycle，DNA replication，mismatch repair和base excision repair方面有突出的富集。我们还评估了其他致癌通路在不同聚类亚型的差异。在Cluster B中NOTCH、NRF2、RAS和TGFb的得分明显偏高，而在Cluster C中细胞周期和TP53信号通路被激活（图5E）。此外，我们还对这三种表型与cGAS-STING信号通路之间的联系进行了研究。结果表明，cGAS（p < 0.001）、IRF3（p < 0.001）和TBK1（p = 0.001）在Cluster C中表达量最高，但在Cluster A中下调最多（补充图7A-D）。这意味着坏死性凋亡在塑造肿瘤免疫微环境和免疫反应中的重要性，并进一步表明坏死细胞在化疗和免疫治疗中的可行性。

图5. 坏死性凋亡和癌症相关通路的分析。（A）热图显示了29个坏死细胞分子的mRNA表达与公认的癌症相关信号通路之间的相关性。只有在至少5种癌症类型中的功能（被抑制或被激活）才被显示出来。红色的 "通路激活 "代表了某条通路可能被坏死性凋亡分子激活的癌症的百分比，而抑制则以类似的方式说明。（B）该网络通过线条连接显示了坏死性凋亡调节基因之间的关系。实线表示激活，而虚线表示抑制。(C). 坏死性凋亡相关分子的功能相似性分布被总结为复杂的圆点图。方框代表相似性的中间50%，上下边界表示第75和25百分位数。框中的线表示功能相似性的中值。连续分布被显示为盒子上的线。圆点代表原始数据。(D). 基因组变异分析富集显示了不同坏死性凋亡表型中生物途径的激活和抑制状态。热图被用于可视化生物学通路的活性。黄色代表激活，蓝色代表抑制。(E). 图中山形图显示了三种坏死细胞聚集亚型中十个重要的知名癌症信号通路。星号表示P值（***p < 0.001，ns是无显著性的缩写），其颜色代表不同的统计测试方法。黑色表示Kruskal-Wallis检验；红色 表示单因素方差分析检验，绿色表示Welch单因素方差分析检验。

NecroScore的临床价值

接下来我们开发了一个名为NecroScore的预后评分系统。在TCGA和GEO-Meta队列的最佳分界点上将患者分为低NecroScore组和高NecroScore组，其中TCGA-LUAD的cutoff为-3.80，GEO的为-0.22（补充表4）。通过Welch单因素方差分析检验比较NecroScore水平，结果显示，NecroScore在Cluster C中最高，其次是Cluster A，然后是Cluster B（图6A）。根据KM生存分析，在TCGA（p = 0.002，图6B）和GEO队列（p < 0.001, 图6C）中，NecroScore较高的患者比NecroScore较低的患者预后更差（p < 0.001）。不仅如此，不同临床阶层的患者，包括年轻人和老年人，女性和男性，以及早期和晚期的患者也具有生存差异（补充图8）。晚期患者也显著集中在高NecroScore组中（图6D）。单变量Cox模型显示，年龄、病理分期和NecroScore与生存率显著相关（图6E，上部），多变量Cox回归验证了年龄（HR=1.028，95%CI=1.016-1.041，p < 0.001），分期（HR = 1.750，95%CI = 1.540-1.986，p < 0.001）和NecroScore（HR = 1.086, 95%CI = 1.040-1.133, p < 0.001）可以作为独立的预后因素（图6E，底部）。通过t检验发现，低NecroScore患者与高NecroScore患者相比，具有明显的治疗优势（p<0.05，图6F）。在三个独立的预后因素中，分期呈现出最高的净效益，其次是年龄，然后是NecroScore（图6G）。此外，根据决策分析，三个变量的组合具有最高水平的3年、5年和10年净获益（图6G）。我们还建立了nomogram模型来预测LUAD患者的1年、3年和5年的生存率（补充图9A）。校准曲线证实了nomogram模型对总生存率的预测是良好（补充图9B）。我们还发现，低TMB和高NecroScore组合的患者预后不良，而高TMB和低NecroScore的患者预后较好（图6H）。
我们进一步研究了NecroScore和免疫细胞浸润之间的相关性。免疫反应的热图是基于TIMER，CIBERSORT，CIBERSORT-ABS，ESTIMATE，QUANTISEQ，MCP-counter和XCELL使用ssGSEA方法绘制的（补充图10A）。它表明，几种细胞类型的浸润程度与NecroScore呈负相关，如myeloid dendritic cells，cancer-associated Fibroblasts，hematopoietic stem cells和monocytes。从另一项研究中收集的免疫细胞的不同基因标记被用来验证这一结果。大多数的标志呈现负相关，特别是在endothelial cells，resting dendritic cells和resting mast cells，而在activated memory CD4+ T cells和Fibroblasts中可以发现弱的关联（补充图10B）。在免疫表型中，较高比例的C1（20%h vs. 6%l），C2（34%h vs. 19%l）和C6（6%h vs. 2%l）在高NecroScore组，而更多的C3（36%h vs. 66%l）和C4（4%h vs. 6%l）样本和C4（4%h vs. 6%l）位于低NecroScore组（补充图10C）。总的分布情况通过Fisher's exact test（p = 0.001）检验，总体分布具有显著性。来自TCIA的证据显示，低NecroScore患者的抗CTLA-4的阳性或阴性状态的IPS明显较低（补充图10D），表明抗CTLA-4抑制剂也可能对PD-1阴性状态的患者产生影响。此外，我们发现几个免疫检查点，如CD274、PDCD1、PDCD1LG2、LAG-3和TMIGD2在低NecroScore组中下调（补充图11A）。有趣的是，TIDE评分、MSI和
功能障碍评分在低NecroScore组中显著增加，同时，低NecroScore患者的排除分数较低（补充图11B）。免疫检查点阴性或弱表达的个体也可能表现出有效的免疫检查点阻断（ICB）治疗的结果。

图6 - NecroScore的临床意义。(A). Welch one-way ANOVA检验显示NecroScore处于三个不同的坏死性凋亡水平。框中的线代表中位数，黑点代表离群值。星号表示P值（***p<0.001）。(B, C). Kaplan-Meier曲线和log-rank检验显示，在TCGA-LUAD和GEO-meta队列中，低和高NecroScore组的总生存概率分别有显著差异。(D). 通过Fisher's精确检验，饼状图展示了低和高NecroScore组之间肺腺癌的临床病理特征的变化。(E). 单因素 (上部）和多因素（下部）Cox模型显示年龄、病理分期和NecroScore的独立预后价值，以森林图显示。(F). 初次治疗反应组（CR，完全反应；PR，部分反应；SD，疾病稳定；PD，疾病进展）后，NecroScore在不同组别中的分布情况用t检验。(G). 决策分析显示了NecroScore的临床有效性。X轴代表阈值概率阈值，Y轴代表3、5、10年的净获益。年龄、病理阶段和NecroScore的组合表现出稳定的最高净效益。(H). 不同肿瘤突变负荷和NecroScore分层的患者的生存分析是通过log-rank检验进行的，并通过Kaplan-Meier曲线显示。"H "代表高，"L "代表低。

高低NecroScore组的突变状态

为了分析LUAD的NecroScore相关机制，我们对TCGA-LUAD队列的体细胞突变进行了分析。通过比较两个评分组的突变频率，在高NecroScore组发现了更多的突变（r = 0.3, p < 0.001; Supplementary Figure 12A），包括非同义突变（r = 0.29, p < 0.001）和同义突变（r = 0.3, p < 0.001; 补充图12B, C）。同时，TP53、TTN、PCDH15、LRP1B、ZFHX4、NALCN、CSMD3、MUC16、XIRP2、APOB 、PAPPA2、PRDM9、ANK2和LRRC7的突变频率明显较高。在高NecroScore组中观察到LRRC7有更多突变，而KEAP1在低NecroScore组有更多的突变（补充图12D）。此外，这些基因之间存在着显著的共发生特征（补充图12E）。我们还发现了高突变基因变异对NecroScore的影响，NecroScore与坏死性相关的基因CDKN2A的突变频率呈正相关（补充图12F）。

药物发现及药物响应分析

我们分析了坏死性凋亡相关分子分子与LUAD治疗的临床效果的关系。根据GDSC的药物反应数据，显示一些基因如BAHC2, LEF1, DNMT1,TERT, TNFRSF1B, TNFSF10, MPG, 和ALK与药物有协同作用；而PANX1、AXL、SLC39A7、TNFRSF21和ID1与药物有拮抗作用（图7A）。ALK表现出与CH5424802（alectinib、ALK和RET抑制剂，图7A）的协同作用。这也，ID1与vorinostat（HDAC抑制剂）、PHA-793887（CDK2抑制剂），PIK-93（PI4K抑制剂），和SNX-2112（HSP90抑制剂）等药物有最强的负相关关系(图7A）。
鉴于坏死调节分子在化疗中的重要性，我们进一步判断了NecroScore能否准确预测LUAD患者的化疗敏感性。我们利用Ridge回归来预测IC50在低和高NecroScore组中的水平。结果显示较高的NecroScore对几种LUAD特异性药物的IC50较低，如顺铂、多西紫杉醇和紫杉醇。与此相反，低NecroScore的病人可能对厄洛替尼（酪氨酸激酶抑制剂/TKI）、Gefitinib（EGFR-TK抑制剂）和拉帕替尼（TKI）更敏感（图7B），这表明化疗是高NecroScore组的一个有效方案。相比之下，低NecroScore组的靶向治疗组可能会导致有利的结果。我们对基于坏死性凋亡表型的分类所得到的细胞系群内药物反应的AUC进行了比较，只有274个在90%以上的LUAD细胞系上测试的化合物被用于分析（补充表5-7）。一种名为SB-590885的一种B-RAF抑制剂和一种G9a/GLP抑制剂UNC0642在第3组LUAD细胞系中的AUCs显著降低，681640和KN001-27的信息并不存在（图7C）。
此外，我们确定了可以靶向干性特征和NecroScore差异的潜在化合物。总共分析了24个化合物的分析，其中19个MoA是由这些化合物共享的（图7D）。四个化合物（uphenazine、uspirilene、mesoridazine和pimozide）共享dopamine receptor antagonist的MoA；metixene和orphenadrine共享了acetylcholine receptor antagonist的MoA；amantadine和memantine共享了glutamate receptor antagonist的MoA；共享谷氨酸受体拮抗剂的MoA。还可以看出，pyrvinium-pamoate与AKT抑制剂的MoA 紧密联系。这些结果有助于根据NecroScore和参考突变特征，为每个LUAD组采取更有针对性的化疗策略（补充图12）。

图7. 基于坏死性凋亡亚型的药物敏感性以及化疗潜力和NecroScore之间的相关性。(A). 气泡图说明了坏死性凋亡调节分子对常见临床化疗药物的耐药性的相关分析结果。红色表示正相关，意味着坏死性凋亡相关分子的mRNA水平过度表达，对药物有抗性，反之亦然。颜色越深，相关系数越大。气泡大小表示错误发现率（FDR）的对数10的负值。(B). Wilcoxon测试显示了低NecroScore组和高NecroScore组之间的半最大抑制浓度（IC50）的差异。星号表示P值（***p < 0.001，*p < 0.05）。(C). Kruskal-Wallis检验揭示了三种肺腺癌细胞系聚类亚型之间药物反应的曲线下面积（AUCs）的差异。(D). 热图显示了每个共享作用机制的潜在化合物（MoA，以行为单位），并按右侧条形图中显示的共享MoA的化合物数量降序排列。

总结：本文基于坏死性凋亡行相关基因识别了三种LUAD亚型，分析了三种亚型间的临床、免疫、基因组、表观组的差异，并利用PCA构建了一个预后评分系统，能够有效预测病人生存和药物治疗响应。该文章内容比较丰富，在routine模式下又有一些新的分析，值得国内临床工作者们借鉴。