Comprehensive Multi‐Omics Integration Reveals B Cell‐Derived ELL2 as a Novel Diagnostic and Prognostic Biomarker in Sepsis

We analyzed bulk transcriptome sequencing data from 10 cohorts of sepsis patients. Using unbiased patient clustering, we identified three subtypes with significantly different prognoses, which were consistently reproduced across all 10 cohorts. Through comprehensive multi-angle analyses, we revealed distinct differences among the subtypes in terms of inflammation, immune responses, and functional pathways. By integrating multiple machine learning algorithms with single-cell transcriptomic analysis, we identified ELL2 as an effective diagnostic and prognostic biomarker for sepsis. This study offers new insights into the mechanisms underlying sepsis progression and highlights the importance of continuous monitoring of ELL2 expression during early diagnosis and treatment. To the Editor, Sepsis is a life-threatening condition characterized by organ dysfunction resulting from the body's dysregulated response to infection [1]. Annually, sepsis affects over 31 million individuals globally, with a mortality rate of approximately 30% [2]. This condition is marked by its rapid progression, poor prognosis, and high mortality, posing significant challenges in intensive care medicine [3]. Early diagnosis and intervention are crucial to improving patient outcomes. However, the greatest challenge in clinical diagnosis and treatment arises from the delayed availability of microbial culture results. Existing biomarkers, although helpful, suffer from limitations such as insufficient sensitivity or specificity, leading to suboptimal treatment timing and increased mortality among sepsis patients. Therefore, there is an urgent need for new biomarkers that can accurately predict prognosis and diagnose sepsis. Recent advances in transcriptomics, proteomics, and metabolomics have enabled the identification of reliable biomarkers. Biomarkers like procalcitonin [4] and C-reactive protein [5] have proven valuable for the diagnosis, prognosis, and treatment of sepsis, aiding clinical decision-making. With the advent of bioinformatics, several studies have employed advanced machine learning methods to develop diagnostic models for sepsis, enhancing patient identification [6-9]. In sepsis research, there are generally two predominant strategies for patient subtyping. One approach relies on feature similarity, where patients are grouped based on molecular or clinical characteristics that resemble one another. This strategy provides valuable insights into the biological and clinical heterogeneity of sepsis [10]. The second strategy, developed under a causal inference framework, focuses on understanding how different sepsis subtypes influence the efficacy of interventions and clinical outcomes. By leveraging this framework, researchers aim to identify subtypes that may respond differently to specific therapeutic treatments [11]. This approach adds a layer of understanding to the management of sepsis, particularly in the context of personalized medicine. We conducted a comprehensive cohort analysis and identified a novel, highly reproducible classification of sepsis patients. Using weighted gene co-expression network analysis mediation analysis, and an integrated multi-machine learning pipeline, we developed an optimal prognostic prediction model. This model effectively identifies high-risk patients. Furthermore, through the integration of machine learning and single-cell sequencing analyses, we discovered the potential of B-cell-derived ELL2 as a novel biomarker for prognosis prediction and diagnostic differentiation in sepsis patients. This finding holds promise for addressing the limitations of current clinical diagnostic tests, improving early diagnosis, and enhancing patient outcomes in sepsis management. We began by processing the GSE65682 dataset, which comprises peripheral blood transcriptome sequencing data from 42 healthy individuals and 760 sepsis patients, with 468 patients having complete overall survival data and status. Using these 468 samples, we performed univariate Cox regression analysis on all genes (Table S1), identifying 2250 genes significantly associated with prognosis (p < 0.05). Based on the expression of these genes, we conducted unsupervised clustering of the 468 samples. We evaluated clustering results for k values ranging from 2 to 9. We determined that k = 3 provided the optimal clustering solution. This result was validated by the proportion of ambiguous cluster scores (Supporting Information S1: Figure S1A). The heatmap clearly showed distinct separation among C1, C2, and C3 (Figure 1A). Notably, Kaplan–Meier curves revealed distinct prognostic differences among C1, C2, and C3, with C2 having the best prognosis, C3 the worst, and C1 intermediate, identifying C3 as the high-risk sepsis patient cluster (Figure 1B). Subsequently, differential expression analysis was performed for each cluster, selecting the top 300 upregulated genes as characteristic genes for each cluster (Supporting Information S1: Figure S1B). We processed nine other sepsis cohorts as validation cohorts. By employing the Nearest Template Prediction method, we extended the subtypes to each validation cohort, consistently identifying clear distributions of C1, C2, and C3 subtypes in each cohort (Figure 1C). In summary, our results demonstrate high stability and scalability of the classification. Identification and Validation of Sepsis Patient Clusters. (A) Graph depicting correlation among patients when K = 3. (B) Kaplan–Meier curves for the three clusters. (C) Validation of subtypes across nine external cohorts using the Nearest Template Prediction algorithm. (D) GO enrichment analysis of top 300 genes specifically upregulated in C3. (E) Pathways specifically upregulated in each subtype based on GSEA analysis. (F) Differential analysis of transcription factor activity across each subtype. We delineated phenotypic and functional differences among various clusters. First, we conducted enrichment analyses on genes specifically upregulated in each cluster (Figure 1D, Supporting Information S1: Figure S1 C-D). The genes upregulated in C2 were enriched in multiple immune-related pathways, such as lymphocyte activation, MHC class II protein complex assembly, and type II hypersensitivity, indicating that C2 exhibits a highly immune-activated phenotype. In contrast, genes upregulated in C1 were enriched in pathways like nucleoside phosphate catabolic process and lipid biosynthetic process, suggesting a predisposition toward metabolic dysregulation in C1 patients. For C3, the upregulated genes were enriched in erythrocyte differentiation and homeostasis, typically occurring under anemic or hypoxic conditions, along with pathways related to porphyrin-containing compound metabolism, gas transport, and pigment biosynthetic process. This implies an increased metabolic demand for erythrocytes and a compensatory hematopoietic state, correlating with poorer prognosis in C3 patients. Gene set enrichment analysis (GSEA) results (Figure 1E) revealed activation of pathways related to response to type I interferon, T cell differentiation, and adaptive immune response in C2 patients. We also analyzed transcription factor activity differences across clusters, finding divergent transcription factor driving patterns (Figure 1F). For example, FOXP2 and GABPA activity increased in C1 but significantly decreased in C3 patients. In C2, IRF2, PFX5, MAZ, IRF4, and LYL1 showed elevated activity, while IRF2, PFX5, and MAZ activity decreased in C1, and IRF4 and LYL1 activity markedly decreased in C3. The increased activity of STAT2 in C2 aligns with its immune-activated state. Additionally, the elevated activity of SREBF2 in C3 is distinctive. In summary, immune suppression emerged as the most prominent feature in C3 patients, correlating with their poor prognosis. To identify genes highly associated with the C3 subtype, we performed causal weighted gene co-expression network analysis (CWGCNA) [12]. CWGCNA extends the traditional WGCNA framework, incorporating mediator models to establish causal relationships between WGCNA modules, module traits, and phenotypes, thus determining whether module changes lead to phenotypic alterations or vice versa. Causality is categorized as forward (phenotype -> module gene -> module) and reverse (module -> module gene -> phenotype). Initially, we conducted WGCNA, selecting a soft threshold of 6 (Figure 2A) to provide optimal power for analysis, leading to the clustering of genes into modules. We identified 29 gene modules, with modules ME1 and ME11 showing strong correlations with the C3 phenotype (Figure 2B), particularly ME11, which had the highest correlation coefficient of 0.78 and comprised 278 genes. Using the limma method, we performed differential analysis of each module between C3 and other clusters, revealing that ME11 exhibited the highest upregulation in C3 (Figure 2C) and the most significant differences (Figure 2D). Within the ME11 module, causal inference using a mediator model identified 96 genes with causal relationships to the C3 phenotype, all in the reverse direction (Table S2), indicating these genes contribute to the formation of C3 rather than being a result of it. Lastly, univariate Cox regression analysis of these 96 genes (Supporting Information S1: Figure S2A) identified 62 genes associated with patient prognosis, all being risk genes (HR > 1, p < 0.05). Identification of B Cell-Derived ELL2 as an Effective Biomarker in Sepsis Patients. (A) Selection of soft threshold for WGCNA analysis. (B) Correlation analysis between different gene modules and phenotypes. (C) p-values from correlation analysis between each gene module and phenotypes. (D) Differential analysis of each module's expression in the C3 phenotype using limma. (E) ROC curves of RS in the training and validation cohorts. (F) Survival analysis of high RS versus low RS patients in the training and validation cohorts. (J) AUC values distinguishing sepsis samples from healthy samples in each of the eight independent cohorts for every gene. (H) Box plot showing differential expression of ELL2 between sepsis and healthy samples across 8 cohorts. (I) Relative RNA expression levels of ELL2 in peripheral blood from 8 healthy individuals and 8 sepsis patients. (J) UMAP dimension reduction landscape in single-cell sequencing analysis. (K) Expression levels of ELL2 in different cell types between sepsis and healthy samples in single-cell sequencing. (L) Bar plot showing differential expression of ELL2 in B cells among healthy samples, surviving sepsis samples, and deceased sepsis samples. To accurately identify high-risk patients for early intervention, we used a prognostic model using a survival framework integrated with machine learning [13] on 62 genes. From the GSE65682 dataset, 468 samples were split into training and validation cohorts in a 7:3 ratio. In the training cohort, we employed 10 machine learning algorithms: RSF, Enet, Lasso, Ridge, stepwise Cox, CoxBoost, plsRcox, SuperPC, GBM, and survival-SVM. Additionally, combinatorial analysis was applied, where the first algorithm selected genes, and the second built the prognostic model. To ensure robustness, models with fewer than two genes were deemed ineffective. This process yielded 99 prognostic models, and the C-index was calculated for each gene signature. Ranking the C-index values in the validation sets, the Lasso+GBM combination achieved the highest C-index of 0.675 (Supporting Information S1: Figure S2B). This model included the five most important genes (ELL2, GABARAPL1, NFYB, C12orf29, GCLC). ELL2 belongs to the RNA polymerase II elongation factor ELL family and plays a role in immunoglobulin secretion in plasma cells [14],GABARAPL1 is involved in mitochondrial autophagy and autophagosome assembly [15], NFYB participates in gene transcription and steroid metabolism [16], C12orf29 is involved in maintaining RNA integrity [17], and GCLC is the rate-limiting enzyme in glutathione synthesis, associated with cellular ferroptosis [18]. Under the Lasso+GBM model, we calculated risk scores (RS) for each patient. Using overall survival and survival status, the ROC curves in the training cohort showed area under the curve (AUC)values for RS at 5, 15, and 25 days as 0.804, 0.829, and 0.829, respectively (Figure 2E). In the validation cohort, the AUC values at 5, 15, and 25 days were 0.773, 0.709, and 0.681, respectively (Figure 2E). Patients were classified into high and low RS groups based on the median RS from the training set. KM analysis indicated that patients in the high RS group had worse prognoses in both cohorts (Figure 2F). In summary, our results demonstrate that the RS is highly stable and effective for early clinical identification of high-risk patients, facilitating prompt diagnosis and treatment. To assess the diagnostic value of the five most important genes in constructing the RS, we computed the AUC values for each gene across 8 cohorts containing both healthy samples and sepsis patient samples (Figure 2G). Following ranking based on average AUC values, ELL2 exhibited the strongest diagnostic performance with an average AUC of 0.865. Across all 8 cohorts, ELL2 showed significant elevation in sepsis patients compared to healthy samples (Figure 2H). Peripheral blood samples were collected from 8 healthy individuals and 8 confirmed sepsis patients, with clinical details provided in Table S3. RT-qPCR results demonstrated significant upregulation of ELL2 in sepsis patients (Figure 2I). These findings collectively indicate that ELL2 expression exhibits accurate predictive and diagnostic capabilities for prognosis. However, bulk transcriptomic sequencing masks cellular heterogeneity. To determine the cellular source of ELL2, we further analyzed single-cell transcriptomic data from the GSE167363 dataset. This dataset included peripheral blood single-cell sequencing data from 2 healthy donors and 10 sepsis patients, with outcomes of 6 patients surviving and 4 deceased. After rigorous quality control, 58,463 cells were included in subsequent analyses. Cells were annotated based on common cell markers (Supporting Information S1: Figure S2B), identifying 8 distinct cell types: B cells, DC cells, T cells, monocytes, NK cells, neutrophils, platelets, and erythroid progenitor cells (Figure 2J). Interestingly, compared to healthy samples, ELL2 expression was significantly elevated in B cells of sepsis patients (Figure 2K). Additionally, ELL2 expression levels remained consistently higher in non-survivors compared to survivors (Figure 2L). In conclusion, integrating prognostic models with single-cell sequencing analysis, we identify B cell-derived ELL2 as a novel diagnostic and prognostic marker for sepsis patients. However, our study has inherent limitations, and future research efforts are needed to elucidate further the functional roles and regulatory mechanisms of ELL2 through additional foundational studies. Additionally, future studies could further integrate transcriptomic data with clinical scores (such as SOFA) to enhance the potential application of this molecular subtype in clinical risk prediction. Through large-scale cohort analysis, we have established a novel and stable classification of sepsis patients. Multi-machine learning integrated prognostic models prove highly effective in identifying high-risk C3 subtype patients, offering significant clinical utility. ELL2 derived from B cells as a potential new biomarker for prognosis prediction and diagnostic differentiation in sepsis patients. Please note: The publisher is not responsible for the content or functionality of any supporting information supplied by the authors. Any queries (other than missing content) should be directed to the corresponding author for the article.