4+昼夜节律+分型+机器学习+单细胞，非肿瘤这样操作也很不错！！

‍

导语

今天给同学们分享一篇生信文章“Exploring the Role of Circadian Rhythm-Related Genes in the Identification of Sepsis Subtypes and the Construction of Diagnostic Models Based on RNA-seq and scRNA-seq ”，这篇文章发表在Int J Mol Sci 期刊上，影响因子为4.9 。

结果：

2.1. 获取和分析差异表达基因的结果

Figure 1 展示了本研究整体流程图。在本研究中，首先对 GSE134347 和 GSE69063 数据集进行了批次校正。从 GSE134347 数据集中提取了昼夜节律相关基因（CRGs），并进行了差异表达分析。十六个差异表达基因（DEGs）在表达上表现出显著差异（ Figure 2 A,B）。在可视化这些基因在染色体上的位置（ Figure 2 C）后，对这些基因进行了相关性分析。大多数基因之间的 Pearson 相关性显著（ Figure 2 D–E）。

整体流程图。

差异表达基因的批次校正与数据识别。(A)上调和下调差异表达基因的表达热图。(B)上调和下调基因的表达箱线图。(C)染色体不同位置差异表达基因的圆环图。(D)差异表达基因之间的相关性圆环图。(E)差异表达基因之间的相关性热图。*** p 小于 0.001；** p 小于 0.01。

2.2. 疾毒亚型的识别

在本研究中，作者使用共识聚类方法对 15 个差异表达基因的基因表达矩阵进行聚类。作者将聚类数量设置为 2-9。 Figure 3 图 A 显示了对应不同 k 值的相对一致性累积分布函数（CDF）曲线下面积的变化。 Figure 3 图 B 显示了不同 k 值下 CDF 的变化。此外，作者还构建了对应其他聚类数量的追踪图（ Figure 3 图 C）。当 k 设置为 2 时，颜色变化频率较低。经过仔细考虑，作者最终将聚类数量设置为 2（置信区间：0.5625-0.8750）。CRGs 在两种亚型之间的表达存在显著差异（ Figure 4 图 A,B）。作者使用 R 软件（4.2.3）中的“ConsensusClusterPlus”包评估了在不同参数组合下聚类结果的稳定性（ Table S2 ）。PCA 分析发现，两种亚型的样本在二维平面上可以被明显区分（ Figure 4 图 C）。使用“limma”算法（ Supplementary Materials “diff_subtype.txt”）进一步分析了两组在差异基因表达上调和下调方面的差异。在 p.adj < 0.05 和|logFC| > 0.5 的条件下，共鉴定出一个显著下调的基因（RORA）。最后，作者对两组样本进行了 GSVA 分析，旨在识别两种亚型所涉及的生物学通路差异（ Figure 4 D）。

基于共识聚类的脓毒症亚型识别结果。(A) k 值从 2 到 9 的共识聚类 CDF。(B) 从 2 到 9 个聚类数的共识聚类结果的累积分布函数(CDF)，每条线代表一个不同的聚类数，显示了在算法多次运行中始终聚在一起样本的比例。(C) 在不同聚类数下构建的跟踪图。图片底部的黑色条纹代表样本，显示了当 k 为 2-9 时样本的分类。不同颜色的块代表不同类别。(D) k 为 2 时的共识矩阵。

两种败血症亚型的表达和富集图谱。（A）两种亚型中 CRGs 的表达热图。（B）两种亚型中 CRGs 的表达箱线图。（C）通过 PCA 获得的两种亚型样本的散点图。（D）两种亚型的 GSVA 结果。*** p < 0.001；** p < 0.01；* p < 0.05。

2.3. WGCNA 分析结果

首先，作者基于 GSE134347 数据集中波动最大的前 25%的基因进行了 WGCNA 分析。为确保共表达模块的整体连通性，关键参数设置如下：“minModuleSize”设置为 50，软阈值功率设置为 6 ( Figure 5 A)。获得了 21 个共表达模块。 Figure 5 B 图展示了彩色编码的基因聚类树。 Figure 5 C 图显示，蓝色模块与基因的相关性最为显著 (p = 9 × 10 −107 , r = 0.93)。 Figure 5 D 图展示了该模块中基因表达与脓毒症样本之间的相关性散点图。

对照组和败血症组的 WGCNA 分析结果。(A) 从 WGCNA 分析中获得的无标度拟合指数和平均连通性线图。(B) 基因聚类树。图像顶部的分支代表基因，图像底部的矩形由各种模块组成。(C) 模块与标签之间关系的关联热图。(D) 蓝色模块与基因之间相关性的散点图。

同样地，作者从该数据集中提取了脓毒症样本，并选择了波动最显著的 25%的基因进行 WGCNA 分析。在将“minModuleSize”设置为 100 并将软阈值功率的关键参数设置为 9 ( Figure 6 A) 后，获得了 13 个共表达模块。 Figure 6 B 显示了两种亚型的基因聚类树。从 Figure 6 C 可以看出，红色模块与基因的相关性最为显著 (p = 4 × 10 −26 , r = 0.74)。 Figure 6 D 显示了该模块中基因表达与 C2 亚型样本之间的相关性散点图。在上述两个相关性最高的模块中交叉基因后，作者获得了 325 个基因 ( Figure 6 E)。

两种败血症亚型的 WGCNA 分析结果。(A) 从 WGCNA 分析中获得的无标度拟合指数和平均连通性线图。(B) 基因聚类树。(C) 模块与标签之间关系的关联热图。(D) 红色模块与基因之间相关性的散点图。(E) 上述两个相关性最高的模块中基因的交集。

2.4. 构建败血症诊断模型及验证结果

在本研究中，作者使用四种机器学习方法（RF、SVM、xgboost 和 GLM）筛选通过 WGCNA 分析获得的 325 个交集基因，并构建脓毒症诊断模型。 Figure 7 图 A 展示了这些算法在构建模型时选出的前 10 个最重要特征。其中，xgboost 算法的残差误差最低，证实该算法在构建诊断模型方面效果更好。在 GSE134347 数据集上，用于构建诊断模型的四种算法的 ROC 分别为 1、1、1 和 0.492（ Figure 7 图 B、 Figure 7 图 C、 Figure 7 图 D 和 Figure 7 图 E）。残差误差最低的 xgboost 算法在 GSE69063 数据集上的 ROC 达到 0.987，从而证明了其高预测准确性。为进一步确认 xgboost 算法在构建诊断模型中的有效性，对构成模型的五个基因（ARHGEF18、CHD3、PHC1、SFI1 和 SPOCK2）分别进行了分析。这些基因在两组间的表达差异显著（ Figure 8 图 A、B 和 Supplementary Materials 图 Figure S2 ）。在 GSE134347 数据集上的 AUC 分别为 0.999、0.991、0.992、1 和 1（ Figure 8 图 C–G）。在 GSE69063 数据集上，AUC 分别为 0.995、0.964、0.947、0.947 和 0.982（ Figure 8 H–L）。此外，这五个基因在两个数据集的两组间表达存在显著差异（p < 0.05）。在外部验证队列（ GSE13904 数据集）中，这五个基因的 AUC 值分别为 0.861、0.759、0.818、0.812 和 0.837（ Figure S3A 和 Figure S4 ）。基于这五个基因构建的诊断模型显示出 0.871 的 AUC，表明其诊断准确度较高。

使用四种机器学习方法构建的诊断模型评估结果。(A) 使用四种方法获得的顶部 10 个重要特征的贡献直方图。(B) 四个模型的累积残差分布图。(C) 四个模型的残差均方根分布的箱线图。(D) 使用四种方法构建的诊断模型在GSE134347 数据集上获得的 ROC 分析结果。(E) 使用 xgboost 方法构建的诊断模型在 GSE69063 数据集上获得的 ROC 分析结果。

在利用 xgboost 方法构建的诊断模型上对GSE69063 数据集进行分析的结果。(A,B)分别展示 GSE134347 和 GSE69063 数据集中对照组和败血症组诊断相关基因的表达热图。(C–G)展示 GSE134347 数据集中 ARHGEF18、CHD3、PHC1、SFI1 和 SPOCK2 的 ROC 曲线。(H–L)展示 GSE69063 数据集中 ARHGEF18、CHD3、PHC1、SFI1 和 SPOCK2 的 ROC 曲线。

2.5. 诊断模型构建结果

作者通过构建与五个基因相关的列线图（ Figure 9 A），直观地理解了败血症风险与五个 FRGs 表达水平之间的关系。使用校准曲线来测量列线图模型的预测能力（ Figure 9 B）。此外，作者还引入了 DCA（ Figure 9 C）。其横轴表示阈值概率，纵轴表示净收益。列线图模型和列线图均被证实具有良好的可靠性。

Nomogram 模型构建结果。(A) 使用诊断相关基因构建的 Nomogram 模型。(B) Nomogram 模型的校准曲线。(C) DCA 结果。

2.6. scRNA-seq 数据分析结果

作者还探索了不同细胞类型中诊断相关基因的表达谱（ Figure 10 ）。具体而言，作者对预处理后的 GSE167363 数据集（ Figure S3 ）进行了聚类和细胞类型注释。获得了七种细胞类型（T 细胞、NK 细胞、B 细胞、单核细胞、血小板、中性粒细胞和红系祖细胞）。CHD3、SFI1 和 SPOCK2 在 T 细胞中高度表达。

scRNA-seq 数据分析结果。(A) 不同细胞类型的 TSNE 降维可视化。(B) 五个诊断相关基因在不同细胞类型中的表达气泡图。(C) 五个诊断相关基因在不同细胞类型中的表达热图。

2.7. 诊断相关基因的表达验证

为确定脓毒症队列与对照组中十六个诊断相关基因的表达谱，作者采用了定量逆转录聚合酶链反应（qRT-PCR）。本研究相关的人口统计学数据提供于 Supplementary Table S1 。作者的研究发现，除了 CSNK1D、AANAT、CLOCK 和 RORA 外，其他基因在两组间的表达水平存在显著差异（ Figure 11 ）。这些结果与作者的生物信息学分析结果一致。因此，进一步探究这些基因在脓毒症中的作用，可能为该疾病的诊断和治疗带来重大进展。

来自脓毒症患者的全血和健康对照的 DBP（A）、CRY2（B）、CSNK1E（C）、PER1（D）、CSNK1D（E）、NR1D1（F）、NPAS2（G）、AANAT（H）、FBXL3（I）、TIMELESS（J）、BHLHE40（K）、PER3（L）、CLOCK（M）、CRY1（N）、RORA（O）和 NR1D2（P）的转录表达，相对于 GAPDH。均值±标准差。*** p 小于 0.001；** p 小于 0.01；* p 小于 0.05。ns：无显著差异。

总结

总之，通过差异表达 CRDRGs 的共识聚类识别出两种败血症亚型，为败血症患者的分层治疗提供了基础。使用五种 CRDRGs（ARHGEF18、CHD3、PHC1、SFI1 和 SPOCK2）构建的基于四种机器学习算法的稳健诊断模型表现出优异的预测性能，在验证集上的 AUC 为 0.987。这些发现不仅为败血症的潜在病理生理机制提供了有价值的见解，还可能为诊断标志物和治疗靶点。未来的工作应集中于在更大队列中验证这些标志物，并探索它们在败血症靶向治疗中的作用。。对这篇文章感兴趣的老师，欢迎扫码咨询！

往期推荐

纯生信选刊

非肿瘤生信

预后模型

单基因生信

单细胞系列