今天给同学们分享一篇坏死性凋亡+预后模型+实验的生信文章“Analysis of necroptosis-related prognostic genes and immune infiltration in idiopathic pulmonary fibrosis”，这篇文章于2023年3月27日发表在Front Immunol期刊上，影响因子为7.3。
IPF是一种未确定的、进行性的肺部疾病。坏死程序性凋亡是一种参与COPD和ARDS等肺部疾病发病机制的程序性凋亡类型。然而，IPF中的坏死程序性凋亡尚未得到充分研究。本研究旨在调查IPF中的坏死程序性凋亡以及坏死性凋亡与免疫浸润之间的关系，构建基于坏死性凋亡相关基因的IPF预后预测模型。
图1 生物信息学分析步骤流程图
1. IPF和对照组之间的差异表达基因（DEGs）和非差异表达基因（NRDEGs）的鉴定
在进行差异分析之前，对GSE110147中的IPF组和对照组进行了主成分分析（PCA），并注意到两组之间的可比性明显。图2A显示了PCA的结果。随后，通过差异分析确定了GSE110147中IPF组和对照组之间的3253个差异表达基因（DEGs）。图2B通过热图可视化展示了DEGs的表达情况。将与GSE110147中的DEGs交集的坏死程序化细胞死亡相关基因（NRDEGs）进行了进一步分析，共有29个重叠的NRDEGs。图2C以韦恩图的形式展示了这29个NRDEGs。图2D和图2E展示了NRDEGs的表达情况，分别通过热图和火山图进行可视化。
图2 IPF中NRDEGs的鉴定
2. PPI网络和功能通路分析
NRDEG的PPI网络是根据STRING数据库构建的（图3A）。使用基于CytoHubba MCC的算法来识别前10个枢纽基因（图3B），包括CFLAR、BIRC2、BIRC3、CYLD、TRAF5、IL1B、STAT1、INFAR1、TICAM1和NLRP3。进行GO和KEGG富集分析以确定NRDEG的生物学功能和相关信号通路。总共鉴定了1805个相关的生物学过程和122个KEGG信号通路。GO富集分析显示，对病毒的防御反应富集、DNA结合转录的阳性调节、I-κB激酶/NF-κB信号传导的调节、肿瘤坏死因子介导的信号传导、坏死过程和生物过程中的程序性坏死细胞死亡（图4A）。此外，本研究还证明了细胞成分中质膜的细胞质侧、膜微结构域、NLRP3炎症小体复合体和膜筏的富集（图4B）。此外，在分子功能中还发现了肿瘤坏死因子受体超家族结合、参与凋亡过程的半胱氨酸型内肽酶活性和细胞因子受体结合的富集（图4C）。KEGG分析显示，NRDEG富集于坏死、NOD-like受体信号通路、TNF信号通路、NF-kappa B信号通路、甲型流感、冠状病毒疾病COVID-19、单纯疱疹病毒1型感染、脂质和动脉粥样硬化、C型凝集素受体信号通路和丙型肝炎（图4D）。
图3 通过STRING数据库和Cytoscape软件分析了PPI网络和关键基因
图4 使用NRDEG的GO和KEGG富集分析结果
3. 免疫浸润和NRDEGs与免疫细胞之间的相关性
通过富集分析发现，NRDEGs主要与免疫相关的信号通路有关。因此，利用ImmuCellAI算法研究了IPF中坏死程序性细胞死亡和免疫浸润之间的关联。结果显示，在IPF和对照组之间存在12种免疫细胞类型的显著差异（图5A-D）。值得注意的是，巨噬细胞、中性粒细胞、CD4 T细胞和nTreg细胞主要富集在IPF中。然而，在对照组中，富集的是DC细胞和NK细胞。此外，spearman相关性检验的结果显示，NRDEGs主要与巨噬细胞和中性粒细胞浸润呈正相关，与DC细胞和NK细胞浸润呈负相关，这些结果通过相关性热图进行可视化（图5E）。
图5 IPF免疫浸润的分析结果
4. 调查肺纤维化中II型肺泡上皮细胞中NRDEG的表达情况
为了进一步探索IPF中坏死性凋亡与细胞类型之间的关系，作者从GSE122960中提取了单细胞测序数据进行分析。结果显示，在纤维化肺组织的AT II中，坏死性凋亡相关基因的表达增加。作者观察到AT II中IPF与对照组之间存在4831个差异表达基因(DEGs)。将这159个坏死性凋亡相关基因与AT II中的DEGs进行叠加，筛选出31个重叠的NRDEGs，包括MLKL、IFNAR1和TRAF5。图6A以Venn图形式展示了这31个NRDEGs。NRDEGs的表达通过火山图可视化，如图6B所示。
图6 IPF中AT II的NRDEGs分析结果
5. IPF中AT II的NRDEGs分析结果
GSE70866包含了176名IPF患者的BAL细胞测序数据，并且通过"sva"软件包消除了批次效应。消除效果以PCA图表呈现（图7A、B）。提取了关键基因的表达矩阵和预后数据，并通过单变量Cox回归分析得到的生存差异基因展示在KM曲线上（图7C-F）。多变量Cox回归分析显示IFNAR1、TRAF5、NLRP3和CYLD与预后有关。结果以森林图的形式展示（图8）。采用Cox逐步回归方法从与坏死程序化相关的关键基因中筛选出4个最佳变量，然后根据多变量Cox回归模型中得到的4个特征关键基因的表达水平和相应的系数，为每个患者估计了一个风险评分。根据风险评分的中位数，将患者分为高风险组和低风险组。与低风险组相比，高风险组的预后明显较差（图9A）。ROC曲线显示评估模型对IPF患者的预后具有良好的预测价值，AUC为0.708（图9B）。本研究还基于Cox模型建立了一个总生存预测的示意图（图9C）。
图7
图8 多变量Cox回归分析中P <0.05的关键基因的森林图
图9
6. NRDEGs表达验证及预后预测模型效果
通过对GSE213001的分析，共发现了27个相同的NRDEGs。在这27个NRDEGs中，有12个基因的表达趋势与GSE110147中的一致，包括BIRC2、BIRC3、CYLD、DNM1L、HSP90AA1、IL1B、NLRP3、PLA2G4A、PLA2G4F、PYCARD、RBCK1和TRAF5（图10A、B）。总体而言，GSE213001中的验证结果与GSE110147比较一致。作者进一步验证了预后预测模型在外部数据集GSE93606中的效果。根据风险评分，作者将IPF患者分为高风险组和低风险组。在GSE93606中，生存比较显示低风险组的预后显著优于高风险组（图10C）。
图10 外部数据集验证
7. IPF A549细胞模型中与死亡相关的预后基因表达的初步验证
对肺纤维化的先前单细胞测序数据集进行的分析显示，坏死程序化细胞死亡与AT II密切相关。作者发现，坏死程序化细胞死亡标记物MLKL以及与坏死程序化细胞死亡相关的预后基因，如IFNAR1和TRAF5，在ATII中高度表达。因此，作者使用TGF-β1构建了IPF A549细胞模型，以验证预后预测模型中四个关键基因的mRNA表达。qPCR的结果显示，在TGF-β处理组中，CYLD、NLRP3、IFNAR1和TRAF5的表达水平较高，而在对照组中相对较低。对照组和IPF组之间的CYLD和TRAF5的表达差异具有统计学意义（图11）。
图11 基因的表达通过qPCR进行验证
8. 在博莱霉素诱导的肺纤维化小鼠模型中验证预后基因的表达
为了进一步验证体外细胞实验的结果，并尽可能模拟肺纤维化患者的病理过程，作者成功构建了一种致肺纤维化小鼠模型，用于体内实验。作者能够研究致肺纤维化小鼠模型中的坏死程序性细胞死亡事件的发生。Western blotting显示，在致肺纤维化小鼠模型的肺组织中，坏死程序性细胞死亡标志物p-MLKL以及四个与坏死程序性细胞死亡相关的预后基因（CYLD、NLRP3、TRAF5和IFNAR1）的表达水平较高，并且差异具有统计学意义（图12A、B）。体外细胞实验的结果通过蛋白质表达水平进一步验证。HE染色和Masson染色的结果显示，与对照组相比，致肺纤维化小鼠模型组的肺组织炎症和胶原沉积显著增加（图12C）。
图12 在小鼠肺纤维化的博莱霉素诱导模型中，通过免疫印迹法验证了基因的表达
总结
根据生物信息学分析，肺纤维化（IPF）患者和对照组的坏死程序化相关基因表达水平存在显著差异。此外，作者发现IPF患者肺组织中多种免疫细胞类型的浸润与坏死程序化之间存在相关性，并且正常肺组织与IPF肺组织之间的免疫细胞浸润也存在差异。通过单细胞测序数据分析，作者发现了AT II与纤维化肺组织中坏死程序化之间的关联，因为纤维化肺组织中AT II的坏死程序化标记物和相关基因的表达显著增加。此外，作者还探索了与坏死程序化相关的预后分子生物标志物，并基于多个标志基因构建了预后预测模型，为IPF的进展和预后预测提供了新的方法。因此，本研究为进一步探索IPF中的坏死程序化铺平了道路。此外，作者应该努力开发新的治疗靶点。
