7+非肿瘤+WGCNA+分型+实验，筛选关键基因进一步分型以及表达验证

今天给同学们分享一篇非肿瘤+WGCNA+分型+实验的生信文章“Identification of molecular subtypes and immune infiltration in endometriosis: a novel bioinformatics analysis and In vitro validation”，这篇文章于2023年8月18日发表在Front Immunol期刊上，影响因子为7.3。

子宫内膜异位症是一种全球性的妇科疾病，影响了6-10%的生育年龄妇女。本研究的目的是调查子宫内膜异位症的基因网络和潜在的免疫浸润标志。

图1 探索子宫内膜异位症生物学特性的方法学流程图

1. 探索子宫内膜异位症生物学特性的方法学流程图

为了评估子宫内膜异位症患者的免疫浸润程度，作者使用CIBERSORT算法分析了子宫内膜异位症中免疫亚群的比例。作者构建了免疫细胞图谱（图2A）和这22种免疫细胞的相关性图谱（图2B）。作者使用GSVA计算了每种免疫细胞类型中所有基因的富集分数。作者观察到c2.cp.v7.4.symbols.gmt数据集中与10个免疫相关通路相关的特定基因集的变化，包括Toll样受体信号通路、Nod样受体信号通路、Rigi样受体信号通路、Jak-stat信号通路、NK细胞介导的细胞毒性、T细胞受体信号通路、B细胞受体信号通路、自身免疫性甲状腺疾病、移植排斥反应和移植物抗宿主病基因集（图2C）。HLA和KIR家族基因是与免疫相关的基因家族，家族中每个基因的表达都表现出不同的特征。HLA-A、HLA-B和HLA-C的表达相对较高，而KIR家族的表达相对较低（图2D）。

图2 免疫细胞浸润评估和免疫相关性分析

2. 鉴定两种子宫内膜异位症亚型并基于免疫特征探索差异表达基因（DEGs）

为了确定子宫内膜异位症中不同免疫亚型之间的生物学差异，作者根据免疫细胞浸润数据的共识聚类将子宫内膜异位症样本分为亚型1和亚型2。这两个亚型明显有区别（图3A）。使用limma鉴定了差异基因表达（P值<0.05且|log2FC|>0.3）。共获得了140个差异表达基因，其中包括29个低表达基因和111个高表达基因。火山图显示，上调表达基因的数量高于下调表达基因（图3B）。在热图中，亚型1中的基因呈总体下调趋势，而亚型2中的基因呈总体上调趋势（图3C）。GSVA的结果显示，B细胞受体信号通路和CD8 TCR下游通路是参与子宫内膜异位症的重要通路。这些发现表明，免疫系统在子宫内膜异位症的病理过程中起着至关重要的作用。

图3 基于免疫细胞亚型的差异分析

&nbsp;

3.&nbsp;共表达网络构建和核心基因鉴定

为了探索与子宫内膜异位症免疫浸润相关的核心基因，作者使用WGCNA算法构建了一个基因共表达网络。通过聚类分析评估了子宫内膜异位症免疫亚型的基因表达谱（图4A）。为了确保网络是无标度的，作者选择了β = 2的软阈值（图4B）。然后，作者将表示矩阵转换为邻接矩阵和拓扑矩阵，并使用平均连接层次聚类方法对基因进行聚类。根据混合动态剪枝树的标准，作者设置了每个基因网络模块中的最小基因数为50个。作者使用动态剪切方法确定基因模块，并计算了每个模块的特征基因值。

图4 加权共表达网络的构建和模块分析

作者随后对模块进行了聚类分析，并使用以下参数合并了彼此相近的模块：高度=0.25，深度分割=4，最小模块=50，共得到11个模块（图4C）。与免疫特征最相关的模块是绿色模块（r = 0.6，P = 5e-12；图4D）。作者在绿色模块中鉴定出了285个基因。然后，作者将这些基因与基于免疫特征的140个差异表达基因进行了比较，并鉴定出了52个共同的基因（图4E）。

4.&nbsp;对子宫内膜异位症关键基因进行GO和KEGG富集分析

为了进一步研究这52个常见基因的功能，作者进行了GO和KEGG分析。作者对这些分析使用了0.5的p值作为截断值。GO分析显示，这些常见基因主要与白细胞介导的细胞毒性、生长的负调控、NK细胞介导的细胞毒性、NK细胞介导的免疫以及蛋白质转运的正调控相关（图5A、C）。KEGG分析显示，这条通路主要与移植物抗宿主病、NK细胞介导的细胞毒性、抗原处理和呈递、补体和凝血级联反应、异体移植排斥以及静息和激活状态下的NK细胞相关（图5B、D）。通路富集分析还显示，这些常见基因在补体和凝血级联反应（图5E）以及NK细胞介导的细胞毒性（图5F）通路中显著表达。

图5 GO、KEGG和PATHWAY富集分析中的核心基因

5.&nbsp;GO、KEGG和PATHWAY富集分析中的核心基因

为了解决富集分析中交集基因潜在偏差的问题，作者对与子宫内膜异位症表达谱中免疫亚型相关的所有140个差异表达基因进行了基因集富集分析（GSEA）。作者在该分析中使用了C2.cp.v7.4.Symbols.gmt参考基因集。结果显示，差异表达基因在趋化因子信号通路、细胞因子与细胞因子受体相互作用、NK细胞介导的细胞毒性、细胞周期检查点、细胞有丝分裂、同源重组修复和反应组有丝分裂前期等通路中显著富集（图6A-G）。作者还使用C2.cp.v7.4.symbols.gmt参考基因集进行了GSVA，以探索子宫内膜异位症发病机制的潜在机制。根据免疫亚型，作者观察到10个基因集存在显著差异，包括先天免疫系统、biocarta csk通路、移植物抗宿主病、WP氧化损伤、中性粒细胞颗粒溶解、CSF3 G CSF信号传导、WP小胶质细胞病原吞噬通路、WP趋化因子信号通路、WP IL3信号通路和干扰素-γ信号传导。其中大多数差异表达基因集与免疫相关的通路有关（图6H）。

图6 GSEA和GSVA对关键基因进行富集分析

6.&nbsp;PPI网络的构建

为了探索52个候选基因的蛋白质功能，作者利用STRING数据库和Cytoscape软件构建了一个PPI网络。在置信度得分为0.55的情况下，作者确定了45个基因之间存在密切的相互作用，而只有8个基因的得分超过了0.7。作者利用0.55的相互作用得分来确定关键基因。利用cytoHubba插件，作者计算了PPI网络的值，并选择了前10个关键基因及其扩展：GZMB，PRF1，KIR2DL1，KIR2DL3，KIR3DL1，KIR2DL4，FGB，IGFBP1，RBP4和PROK1。为了进一步探索上游调控关系，作者利用multiMiR和TarBase v.8数据库，选择“Functional MTI”和“positive”证据的结果，预测了miRNA与关键基因的相互作用。作者发现关键基因中的8个mRNA与11个miRNA相互作用。

7.&nbsp;基于关键基因的无监督聚类鉴定三种子宫内膜异位症亚型

为了研究基于子宫内膜异位症特征的不同亚型之间的生物学差异，作者使用ConsensusClusterPlus软件包根据10个关键基因的表达谱构建了亚型。当k = 3时，分类是可靠且稳定的（图7A-C）。样本被分为亚型1、亚型2和亚型3。PCA证实了亚型1、亚型2和亚型3之间存在显著差异（图7D）。作者还分析了疾病亚型与免疫细胞比例之间的相关性，并发现三个组分在不同免疫细胞中表达不同。激活的NK细胞、CD8阳性T细胞和静息的NK细胞的表达差异较大。对于关键基因，作者计算了关键基因表达与免疫细胞比例之间的相关性（皮尔逊系数）。大多数关键基因与免疫标记物显著相关。作者观察到GZMB与激活的NK细胞之间存在显著正相关（P = 1.23e-11，R = 0.59），KIR2DL1与静息的NK细胞之间存在显著正相关（P = 2.51e-8，R = 0.50），KIR2DL3与NK细胞之间存在显著正相关。细胞休息（P = 1.55e-11，R = 0.59），KIR2DL4和休息的NK细胞（P = 4.26e-12，R = 0.60），PRF1和活化的NK细胞（P = 3.54e-10，R = 0.56），以及PRF1和休息的NK细胞（P = 1.15e-11，R = 0.59）。对中枢基因和免疫检查点基因之间的相关性进一步分析表明，中枢基因与CD40、IDO1、LAG3、TNF和TNFRSF18免疫检查点基因呈强正相关。这些结果表明，关键基因的表达与免疫特征显著正相关。

图7 基于子宫内膜异位症相关分子分型的特征

8.&nbsp;关键基因的表达验证

为了验证关键基因的转录和蛋白表达，作者使用qRT-PCR和Western blot分析检测了来自同一患者的卵巢子宫内膜异位组织和正常内膜组织中它们的表达情况（图8）。作者发现，除了FGB、KIR2DL1和KIR2DL3之外，所有其他关键基因在子宫内膜异位和正常对照组中的表达水平均显示出显著差异（图8A）。具体而言，KIR2DL4的表达明显下调，而PROK1、IGFBP1、RBP4、G2MB、KIR3DL1和PRF1的表达在卵巢子宫内膜异位中显著上调。此外，RBP4的蛋白表达在卵巢子宫内膜异位中明显高于正常对照（图8B）。

图8 关键基因的表达验证

总结

总之，这项研究从免疫浸润的角度为子宫内膜异位症的发病机制提供了全面可靠的证据。结果显示，大部分差异表达基因参与了与免疫相关的途径。与免疫标志物显著相关的十个中心基因可能成为子宫内膜异位症的潜在诊断和治疗靶点。