Introduction
Sepsis is a life-threatening organ dysfunction caused by a dysregulated host response to infection, and it remains a clinical syndrome with a high mortality rate.1 Approximately 20–30% of patients in intensive care units develop sepsis.2 According to data from the Global Burden of Disease study in 2017, there were 48.9 million cases of sepsis worldwide.3 Although the age-standardized incidence and mortality rates have declined, sepsis continues to be a leading cause of global health loss, with the highest burden observed in sub-Saharan Africa, Oceania, South Asia, East Asia, and Southeast Asia.2,3 Current clinical criteria for sepsis diagnosis integrate the SOFA (Sequential Organ Failure Assessment) score with laboratory biomarkers. While the SOFA score is employed to quantify organ dysfunction and predict mortality, it suffers from delayed diagnosis requiring serial monitoring. In contrast, the qSOFA (quick SOFA) score enables rapid bedside assessment of disease severity but demonstrates reduced sensitivity and unreliability in immunocompromised patients.4 Supplementary biomarkers such as procalcitonin (PCT), C-reactive protein (CRP), lactate, and blood cultures exhibit inherent limitations. These include insufficient timeliness (eg, delayed elevation of PCT/CRP leading to treatment delays), limited specificity, and an inability to differentiate immune endotypes. Collectively, these issues elevate the risk of misdiagnosis. Recent studies further highlight the urgent need to develop novel sepsis biomarkers to enhance capabilities for early diagnosis and precision therapeutics.5 Novel biomarkers aim to address these unmet clinical needs by potentially offering enhanced specificity, facilitating the differentiation of sepsis from non-infectious syndromes, enabling early detection via rapid pathophysiological changes, and allowing for immune phenotypic stratification. This capability ensures more timely interventions and improved clinical outcomes. Consequently, this study focuses on characterizing such biomarkers to overcome the limitations inherent in current diagnostic approaches. At present, the treatment of sepsis mainly involves three aspects: infection control, hemodynamic management, and modulation of the host response. Infection control applies to all sepsis cases.6 Due to the heterogeneity of the sepsis process, the diagnosis and definition of sepsis remain key issues in clinical practice,6,7 therefore, underscoring the urgent need for improvements in sepsis diagnosis and treatment.
Parthanatos, a form of programmed cell death distinct from apoptosis, necrosis, and autophagy, is triggered by DNA damage and mediated by poly ADP-ribose polymerase 1 (PARP1).8 The main mechanisms include DNA damage, overactivation of PARP1, accumulation of PAR, depletion of nicotinamide adenine dinucleotide (NAD) and adenosine triphosphate (ATP), and nuclear translocation of apoptosis inducing factor (AIF).9 Parthanatos has been considered a promising therapeutic target in cancer.10 It has been reported in various diseases, including ischemic stroke and neurodegenerative diseases.11 In cancer, parthanatos has been linked to diseases such as hepatocellular carcinoma and triple-negative breast cancer.10 In the context of sepsis, Lorente L. et al investigated the relationship between parthanatos and sepsis mortality, finding an association between parthanatos and sepsis-related deaths,12 although the underlying molecular mechanisms remain unclear. Furthermore, it has been shown that caspase-11 parthanatos leads to pathological changes and shortens survival time in sepsis models. However, parthanatos has been less studied in the development of sepsis and further research is needed.
Mendelian Randomization (MR) is a novel epidemiological approach that uses single nucleotide polymorphisms (SNPs) as instrumental variables to infer causal relationships between exposures and outcomes based on genome-wide sequencing data.13 MR is built upon three key assumptions: (1) the genetic variant is strongly associated with the exposure (relevance assumption), (2) the genetic variant is not associated with confounders of the exposure-outcome relationship (independence assumption), and (3) the genetic variant influences the outcome solely through the exposure (exclusivity assumption).14 Compared to randomized controlled trials, MR offers statistical advantages by mitigating bias and reverse causation, while being more cost-effective.15
In summary, this study utilized sepsis transcriptome data from the Gene Expression Omnibus (GEO) database and 11 parthanatos-related genes (PRGs) from the literature to obtain biomarkers related to parthanatos in sepsis by differential expression analysis, weighted gene co-expression network analysis (WGCNA), MR analysis and machine learning. The diagnostic potential and underlying molecular mechanisms of these biomarkers were further analyzed, providing new insights into the clinical diagnosis, prevention, and treatment of sepsis.
Materials and Methods
Data Provenance
GSE65682 (GPL13667 platform) served as the training cohort, and GSE95233 (GPL570 platform) was used as the validation cohort, both of which were extracted from the GEO database (http://www.ncbi.nlm.nih.gov/geo/). GSE65682 comprised whole blood samples from 760 individuals diagnosed with sepsis and 42 samples from healthy individuals. Additionally, GSE95233 consisted of whole blood samples from 102 patients with sepsis alongside 22 samples from control subjects. The GSE167363 dataset was extracted from the GEO database for single-cell analysis, including 6 human peripheral blood mononuclear cell (PBMC) samples, comprising 2 healthy controls, 2 sepsis survivors, and 2 sepsis non-survivors. The parthanatos-related genes (PRGs) were collected from relevant literature.16 The sepsis dataset (ieu-b-4980) was gotten from Integrative Epidemiology Unit (IEU) Open Genome-Wide Association Study (GWAS) database (https://gwas.mrcieu.ac.uk/), including 486,484 samples of Europeans (11,643 sepsis and 474,841 control samples) and 12,243,539 SNPs. Acquisition of expression quantitative trait loci (eQTL) data of candidate genes as exposure factors via GWAS database.
Differentially Expressed Gene Identification
Differential expression analysis in training set was performed utilizing limma package (v 3.50.1)17 to derive differentially expressed genes (DEGs) between sepsis and control samples (padj < 0.05 and |log2Fold Change (FC)| > 0.5). Subsequently, DEGs were visualized by ggplot2 (v 3.4.4)18 pheatmap (v 1.0.12).19 The top 10 genes were annotated in each of the upregulated and downregulated DEGs, and these twenty genes were selected for heatmap.
Weighted Gene Co-Expression Network Analysis (WGCNA)
WGCNA package (v 1.7.1)20 was used to analyze gene modules associated with PRGs. Initially, we calculated the scores for PRGs in both sepsis patients and healthy samples via single-sample gene set enrichment analysis (ssGSEA) algorithm. Then, the GoodSamplesGenes function was utilized to cluster and exclude outlier samples. PckSoftThreshold function was employed to determine the soft value required to construct the scale-free network (R2 = 0.8620). Following this, the Dynamic Tree Cutting method was implemented to construct the hierarchical clustering tree. Finally, we established a correlation analysis by calculating the matrix of correlation coefficients between the eigenvectors of the modules and the ssGSEA scores.
Pinpointed Candidate Genes and Functional Evaluation
Candidate genes were pinpointed by intersecting DEGs with the genes from the key modules by the ggVenn package (v 0.1.9).21 Subsequently, in order to understand the possible pathways and biological functions of these candidate genes, we conducted Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis on candidate genes through clusterProfiler (v 4.7.1.003)22 (p < 0.05). Then, we investigated the interactions among the candidate genes on proteomic level, the genes were inputted into the STRING database for the acquisition of protein-protein interaction (PPI) data. (confidence > 0.4), and the network of candidate genes was then visualised using Cytoscape software (v 3.10.1).23
Mendelian Randomization (MR) Analysis
To identify candidate genes that had a causal relationship with sepsis, the TwoSampleMR package (v 0.6.4)24 extract_instruments function was harnessed for reading exposure factors and identify SNPs that were independently associated with these factors (p < 5×10-8, R²< 0.001, distance > 100kb and removed SNPs with F-statistics < 10). These SNPs were later applied as instrumental variables (IV) to determine the causality between the exposure and the resulting outcome. Subsequently, harmonise_data function was applied to standardize the effect alleles and effect sizes, and MR analysis was carried out combining 5 algorithms [MR Egger,25 Weighted median,26 Inverse variance weighted (IVW),27 Weighted mode,28 Simple mode].29 Among these, the IVW algorithm requires attention (p<0.05, OR ≠ 1) when assessing the results of the five algorithms. Following this, scatter plots (the negative slope suggested that the gene acts as a protective gene against Sepsis, while the positive slope suggested it as a risk factor for the condition), forest plots (the genes with effect size to the left of the dashed line are considered protective genes, whereas those to the right are identified as risk genes), and funnel plots (If the samples were evenly distributed around the IVW line,it would imply that MR adheres to Mendel’s second law) were utilized to further illustrate the MR results.
Heterogeneity test was conducted via mr_heterogeneity function. (Q_pval > 0.05), mr_pleiotropy_test function and the mr_presso function from the MRPRESSO package (v 1.0)30 were used to perform horizontal pleiotropy tests (p >0.05), and the mr_leaveoneout function performs leave-one-out analysis, which involves sequentially removing each SNP to see if the results change. Finally, the Steiger directionality test was employed to verify that the positive analysis structure was not affected by reverse causal effects, which was done by directionality_test function (steiger pval < 0.05). Genes analyzed by MR were defined as candidate signature genes for further research.
Identification and Validation of Biomarkers
After identifying the candidate signature genes, the least absolute shrinkage and selection operator (LASSO) and Boruta algorithms were instrumental into filter out candidate biomarkers. The intersection of the results from both methods was used to determine the candidate biomarkers. The LASSO was completed through glmnet package (v 4.1–4),31 and the Boruta analysis was accomplished with the Boruta package (v 8.0.0).32 Following this, within GSE65682 and GSE95233 cohorts, performed expression profiling and receiver operating characteristic (ROC) analysis on the candidate biomarkers, ROC curves were generated through pROC package (v 1.18.0).33 Candidate biomarkers with consistent expression trends, significant differences in expression profiles, and the area under the curve (AUC) greater than 0.7 were selected as biomarkers for subsequent analysis.
Constructed of Nomogram
To assess the predictive value of the biomarkers in Sepsis, a nomogram based on the biomarkers was constructed using the rms package (v 6.5–0).34 Subsequently, the accuracy of the nomogram’s predictive ability was evaluated through calibration plots and decision curve analysis. The accuracy of the nomogram was further verified using the ROC. The ggDCA package (v 1.2)35 was used for plotting the decision curves.
Gene Set Enrichment Analysis (GSEA)
To identify the signaling pathways in which the biomarkers are involved, based on the sepsis samples and healthy samples in GSE65682, the Spearman correlation was calculated between the biomarkers and all other genes. The genes were then ranked according to the correlation coefficients. Subsequently, the clusterProfiler (v 4.7.1.003) was used to perform GSEA, and the reference gene set chosen was “c2.cp.kegg.v2023.1.Hs.symbols.gmt” (p < 0.05 and false discovery rate (FDR) < 0.25). Selected the top 5 enriched pathways from the GSEA results for presentation.
Analysis of Immune Infiltration
Subsequently, we investigated the differences in immune cell infiltration in GSE65682 by utilizing the CIBERSORT algorithm to quantify the levels of 22 specific immune cell types, then applied to evaluate disparities in the expression of different immune cell (p < 0.05). Thereafter, we explored the relationships between biomarkers, immune cells, and among immune cells themselves, conducting a Spearman correlation analysis on biomarkers and immune cells to assess their interrelatedness. (p < 0.05, |cor| > 0.3).
Molecular Regulatory Networks of Biomarkers
Molecular regulatory networks are key elements in the process of gene expression regulation, and they are of great research value as they can uncover the intricacies and variety of gene expression control mechanisms. Utilizing the biomarkers obtained, we leveraged the miRNet database (https://www.mirnet.ca/miRNet/docs/RTutorial.xhtml) to identify microRNAs (miRNAs) and long non-coding RNAs (lncRNAs) that were involved in the regulation of these biomarkers. We then filtered for miRNAs that simultaneously target the biomarkers and predicted the related lncRNAs using these miRNAs, using these as the core to build a lncRNA-miRNA-mRNA regulatory network. In order to examine the interplay between biomarkers and transcription factors (TFs), we utilized the miRNet database to forecast TFs that have interactive relationships with the biomarkers and established a miRNA-mRNA-TF regulatory network. Both of regulatory networks were completed using Cytoscape (v 3.10.1).
Subsequently, we conducted research on the biological characteristics of these biomarkers, utilizing GeneMANIA to construct and analyze the interaction network among the biomarkers as well as their co-expressed genes.
Chromosomal Localization of Biomarkers
Gene localization is essential for exploring the structure, function, and interactions of genes, and it helps in understanding the impact of genetic factors on gene expression more profoundly. In order to understand the chromosomal locations of biomarkers, we employed the OmicCircos (v 1.36.0)36 to visualize the distribution of these biomarkers across the chromosomes and illustrate the expression profiles of biomarkers on various chromosomes.
Biomarker-Based Drug Prediction
Based on drug signature database (DsigDB) (http://dsigdb.tanlab.org/DSigDBv 1.0/), we utilized the enrichR package (v 3.1)37 to identify candidate drugs that may target the biomarkers (p < 0.05). Furthermore, we employed Cytoscape software to present a visualized network of interconnections between the biomarkers and a range of predicted drugs.
Single-Cell Analysis
The single-cell transcriptomic sequencing data (GSE167363) were constructed into Seurat objects using the R package “Seurat” (v 5.1.0),38 with cells having fewer than 300 genes and genes covered by fewer than three cells excluded. The parameters were set to 6000 > nFeature_RNA > 200 and nCount_RNA < 25000. The data were then normalized. The FindVariableFeatures function was used to identify highly variable genes based on the relationship between the mean and variance of gene expression, with the top 2,000 highly variable genes selected for further analysis by default. The results were visualized using the LabelPoints function, with the top 10 most variable genes labeled. To further confirm and analyze the cell populations of different cell groups, the R package “Seurat” (v 5.1.0) was used to normalize the single-cell transcriptome sequencing data (GSE167363) using the ScaleData function. Principal component analysis (PCA) was then performed on the highly variable genes of each sample for dimensionality reduction. The Jackstraw function was used to generate a Jackstraw plot, and re-clustering was performed using a permutation test algorithm to select the appropriate principal components. The ElbowPlot function was used to create an elbow plot to identify the usable dimensions. Based on the plateau in the PCA feature number shown in the elbow plot and the analysis from the Jackstraw plot, suitable principal components were selected for further analysis. The FindNeighbors and FindClusters functions in Seurat (v 5.1.0) were used to identify different cell clusters, and uniform manifold approximation and projection (UMAP) was employed for visualization (resolution = 0.1). For cell type annotation within the clusters, the CellMarker2.0 database (http://bio-bigdata.hrbmu.edu.cn/CellMarker) and manual labeling were used to annotate each cell cluster. Cell types were annotated based on the marker genes provided in the literature.39 Finally, UMAP plots of the annotated cell clusters were generated, and a DotPlot was created to display the distribution of biomarkers across different groups and cell types. In addition, the expression of biomarkers in different cells was verified to identify key cells.
Cellular Communication and Pseudotime Analysis
For cellular communication analysis, the CellChat package (v 1.6.1)40 was used, leveraging the CellChatDB database (https://www.cellchat.org/db/). This tool utilized cell expression data as input to model cell-cell communication through ligand-receptor and cofactor interactions. Additionally, pseudotime analysis was performed using the Monocle2 package (v 2.26.0)41 to explore key cellular transitions during different developmental stages, and DifferentialGeneTest was employed to analyze gene dynamics during the cellular differentiation process.
RNA Extraction and Reverse Transcription-Quantitative PCR (RT-qPCR)
Analyzed the expression of biomarkers in blood samples from sepsis patients and control subjects through RT-qPCR. Single blinding was applied during the validation phase. Total RNA was extracted from blood samples of 5 sepsis patients and 5 control subjects using TRIzol reagent (Ambion, Austin, USA). After standing for 15 minutes, 1 μL was taken to quantify RNA concentration using a NanoPhotometer N50. Then, reverse transcription to synthesize cDNA was performed using the SureScript First strand cDNA synthesis kit (Servicebio, Wuhan, China). The obtained cDNA was diluted and then subjected to RT-qPCR experiments according to a certain system. The results were analyzed using the 2-ΔΔCt method, with GAPDH serving as the reference gene to ensure the accuracy of the experiment. The specific primer sequences, reaction systems, and amplification conditions are shown in Supplementary Tables 1–3.
Statistical Analysis
All statistical analyses were conducted using R software (v 4.2.3). The Wilcoxon test was employed to assess differences between two groups, considering p < 0.05 as indicative of statistical significance.
Results
Screening of DEGs and Key Module Genes
Within the GSE65682, we identified 4,336 DEGs between sepsis patients and control groups, comprising 1,370 upregulated genes and 2,966 downregulated genes. The top ten genes from both the upregulated and downregulated DEGs were selected for marking and visualization through a heatmap (Figure 1A and B). Following that, we established a co-expression network for genes based on PRGs using WGCNA. At the outset, we analyzed the expression of PRGs in sepsis and control groups, selecting 9 PRGs with significant scoring differences (Figure 1C). A hierarchical clustering analysis was then conducted on all samples, from which we retained 801 samples (Figure 1D). We chose β= 8 to build a scale-free network (R² = 0.8620) (Figure 1E). After that, a hierarchical clustering tree was constructed to classify the genes and total of 17 co-expression modules were obtained (Figure 1F). Ultimately, the MEgreen module (cor = 0.63, p = 9e-91) as the key module for further analysis (Figure 1G). Finally, we retained 157 genes within the MEgreen module (Gene Significance (GS) > 0.2, Module Membership (MM) > 0.2) (Figure 1H).
|
Figure 1 Differential expression analysis. (A) Volcano map of differentially expressed genes; (B) Heat map of expression of differentially expressed genes; (C) Scoring violin plot for ssGSEA analysis. (D) Sample level clustering tree; (E) Selection of soft thresholds. (F) Identification of co-expression modules; (G) Heatmap of correlation between modular genes and single-sample Gene Set Enrichment Analysis (ssGSEA) score constructs; (H) Scatterplot of correlation between modules and traits.
|
Identification and GO and KEGG Analysis of Candidate Genes and Construction of PPI Networks
In the Venn diagram, 84 candidate genes were acquired through intersection of DEGs and MEgreen module (Figure 2A). Following this, GO analysis showed a significant enrichment of these candidate genes across 119 distinct categories, including 9 Cellular Components (CC), 61 Molecular Functions (MF), and 49 Biological Processes (BP), CCs featured complexes like the serine/threonine protein kinase complex and the protein acetyltransferase complex, MFs involved binding activities such as vinculin binding and the activity of cysteine-type endopeptidases in the apoptotic process, BPs included the regulation of transcription elongation from a DNA template and the positive regulation of signaling pathways involving the transforming growth factor beta receptor (Figure 2B). Additionally, the candidate genes were found to be enriched in KEGG pathways such as the cell cycle, the FoxO signaling pathway, and the Apelin signaling pathway (Figure 2C). Next, to explore the interactions among the candidate genes, we established a PPI network that comprised 60 nodes and 72 edges., where we excluded 24 genes that did not have any interactions (Figure 2D).
 |
Figure 2 Acquisition and functional analysis of candidate genes. (A) Venn diagram; (B) Results of Gene ontology (GO) enrichment analysis; (C) Results of kyoto encyclopedia of genes and genomes (KEGG) enrichment analysis. (D) Protein-Protein Interaction (PPI) network architecture.
|
BRD1, BTN2A1, MSL2, GCC1, NCOA6 and FOXJ3 Were Identified as Candidate Signature Genes
To identify candidate genes that have a causal relationship with sepsis, we conducted MR analysis. In the IVW algorithm, 8 candidate genes were significantly causally associated with sepsis, where BRD1 (β= 0.127, p = 0.004, 95% confidence interval (CI) = 0.040–0.214), CASP2 (β= 0.217, p = 0.000, 95% CI = 0.110–0.323), HNRNPF (β= 0.067, p = 0.044, 95% CI = 0.002–0.132), and GCC1 (β= 0.111, p = 0.000, 95% CI = 0.057–0.166) were risk factors, while BTN2A1 (β= −0.037, p = 0.003, 95% CI = −0.061 to −0.013), MSL2 (β= −0.109, p = 0.008, 95% CI = −0.189 to −0.029), NCOA6 (β= −0.077, p = 0.001, 95% CI = −0.122 to −0.032), and FOXJ3 (β= −0.190, p = 0.002, 95% CI = −0.308 to −0.072) were protective factors (Table 1). In the funnel plot (Figure S1), the SNPs associated with the candidate genes were found to be largely symmetrically distributed, suggesting that the MR analysis adheres to the principles of the second law of mendelian randomization. The outcomes of the scatter plot (Figure S2) and forest plot (Figure S3) further confirmed that BRD1, CASP2, HNRNPF, and GCC1 act as risk genes, whereas BTN2A1, MSL2, NCOA6, and FOXJ3 serve as protective genes. The p-value for the IVW heterogeneity test of Cochran’s Q for the 8 candidate genes was > 0.05 (Table 2), and except for CASP2 and HNRNPF, the p-value for the horizontal pleiotropy test of MR Egger for the remaining 6 candidate genes was also > 0.05 (Table 3). The LOO analysis was conducted for these remaining 6 candidate genes and the outcomes confirmed the reliability and stability of the MR analysis. (Figure S4). To verify that the results of the forward analysis were not confounded by reverse causal effects, the Steiger test was performed. The Steiger test for the 6 candidate genes all returned true (Table 4), and thus, these candidate genes were defined as candidate signature genes.
 |
Table 1 Mendelian Randomization Analytic Result
|
 |
Table 2 Results of the Heterogeneity Test
|
 |
Table 3 Results of the Horizontal Pleiotropy Test
|
 |
Table 4 Results of Steiger Test
|
BRD1 and FOXJ3 Were Identified as Biomarkers
After obtaining the candidate signature genes, LASSO regression and Boruta were utilized to identify candidate biomarkers. In the LASSO analysis, we retained five candidate signature genes (Figure 3A and B), and in the Boruta analysis, we retained six candidate signature genes (Figure 3C). By taking the intersection of these two results, we obtained five candidate biomarkers, which were BRD1, BTN2A1, MSL2, GCC1, NCOA6 and FOXJ3 (Figure 3D). Then, we proceeded to verify the expression levels of five candidate biomarkers within GSE65682 and GSE95233. Only BRD1, MSL2, NCOA6, and FOXJ3 exhibited significant decrease in expression levels in septic samples, with consistent trends across both cohorts (Figure 3E). The AUC for BRD1 and FOXJ3 were 0.940 and 0.953 in the training cohort and in the validation cohort, they were 0.946 and 0.816, consequently, BRD1 and FOXJ3 were as biomarkers for subsequent analysis (Figure 3F).
 |
Figure 3 Screening of key genes for diagnostic analysis. (A) Spectrogram of coefficients; (B)The ten-fold cross-validation; (C) Boruta algorithm to screen for characterized genes; (D) Venn plot for identification of key genes; (E) Expression difference of key genes in validation set and training set; (F) Receiver operating characteristic curve of key genes.
|
Constructing Biomarkers-Related Nomogram
We evaluated the prognostic value of the biomarkers in sepsis by developing a nomogram based on BRD1 and FOXJ3 (Figure 4A). With a p-value of 0.206 and a mean absolute error (MAE) below 0.05 (Figure 4B), it was demonstrated that the nomogram possesses high-precision in its predictive ability. The net benefit is positive in the decision curve further confirms the nomogram’s predictive accuracy (Figure 4C). Ultimately, the AUC of 0.96 substantiates the value of the nomogram in the prediction of sepsis (Figure 4D). These results indicate that the nomogram constructed based on biomarkers has a high predictive accuracy.
 |
Figure 4 Construction of the nomogram diagram model. (A) Nomogram; (B) Calibration curves for nomogram; (C) Decision curve analysis curve; (D) Receiver operating characteristic curves for nomogram.
|
Immune Infiltration Analysis and GSEA
The composition of immune cells in sepsis and control groups was shown in Figure 5A. We then examined the expression patterns of 22 types of immune cells and found significant differences in naive B cells, resting dendritic cells, CD8 T cells, and 15 other immune cell types (Figure 5B). Correlation analysis among these differentially expressed immune cells revealed significant correlations of various extents within the 18 identified cell types, such as there was a significant negative correlation between B cell memory and naive B cells. Additionally, monocytes show significant negative correlations with gamma delta T cells, macrophages, neutrophils, and CD4 naive T cells (Figure 5C). Lastly, the investigation into the relationships between the biomarkers and the differential immune cells uncovered that BRD1 is significantly negatively correlated with Macrophages M0 (p < 0.001, cor = −0.311) and Macrophages M1 (p < 0.001, cor = −0.304). Similarly, FOXJ3 exhibited a significant negative correlation with Macrophages M0 (p < 0.001, cor = −0.308) (Figure 5D).
 |
Figure 5 Immune infiltration and functional enrichment analysis. (A) Cell abundance map of immune infiltration; (B) Differences in immune cell expression between sepsis and control groups; (C) Correlation heat map of immune cells; (D) Heatmap of correlation between differential immune cells and biomarkers; (E) Gene set enrichment analysis of BRD1; (F) Gene set enrichment analysis of FOXJ3.
|
Afterward, we performed GSEA to investigate the potential pathways related to the biomarkers. BRD1 was found to be enriched in 365 pathways, while FOXJ3 was enriched in 354 pathways. Notably, both BRD1 and FOXJ3 were associated with functions like extracellular matrix organization, keratinization, and mRNA splicing (Figure 5E and F).
Exploration of the Regulatory Relationship for BRD1 and FOXJ3
Based on the miRNet database, we collectively predicted 562 miRNAs and 982 lncRNAs. Specifically, for FOXJ3, we predicted 513 miRNAs and 956 lncRNAs, and for BRD1, we predicted 189 miRNAs and 789 lncRNAs. There were 26 miRNAs that commonly regulated both BRD1 and FOXJ3. Using these 26 miRNAs as the core, we matched them to 299 lncRNAs and constructed a lncRNA-miRNA-mRNA network with 327 nodes and 1,434 edges (Figure 6A). In the JASPAR database, we identified 14 transcription factors (TFs) that regulate the biomarkers. For FOXJ3, we predicted 11 TFs, and for BRD1, predicted 3 TFs, then constructed a TF-mRNA regulatory network that includes 16 nodes and 14 edges (Figure 6B).
 |
Figure 6 Molecular regulatory networks. (A) lncRNA-miRNA-biomarker network. Orange represents biomarkers, green represents key miRNAs, and yellow represents lncRNAs. (B) TF-biomarker network. Orange represents biomarkers, and blue represents transcription factors (TFs). (C) GGI network. Different colored lines represent different interactions, and different colored blocks represent the potential functions involved.
|
To further identify genes related to the biomarkers, GeneMANIA was used to discover 20 genes associated with the biomarkers. Within the interaction network, physical interactions accounted for the highest proportion (77.64%), followed by co-expression interactions (8.01%) and the biomarkers and the genes that interact with them were functionally associated with acetyltransferase complexes, protein acetyltransferase complexes, and the process of protein acetylation (Figure 6C).
Chromosomal Localization and Drug Prediction
Chromosomal localization analysis results show that FOXJ3 is located on chromosome 1, and BRD1 is located on chromosome 22 (Figure 7A). To identify potential drugs for targeted therapy of the biomarkers, we conducted a search in the Dsigdb database and discovered 18 significant drugs, among which digoxin, doxorubicin and daunorubicin simultaneously target both BRD1 and FOXJ3. Subsequently, an mRNA-drug network consisting of 20 nodes and 24 edges was constructed (Figure 7B).
 |
Figure 7 Chromosome localization and drug prediction results. (A) The localization of biomarkers on chromosomes. The outer circle numbers represent different chromosomes, and the inner circle shows the position of the biomarkers on the chromosomes. (B) Drug-biomarker network diagram. Orange represents biomarkers, and blue represents drugs.
|
Single-Cell Analysis for FOXJ3
After rigorous quality control, 16,921 high-quality cells and 20,696 genes were retained for further analysis (Figure S5a and b). The top 2,000 highly variable genes were selected, and 10 of the most variable genes, including CXCL8, CCL4, HBD, and LCN2, were highlighted (Figure S5c). Principal component analysis (PCA) was then performed on these genes for dimensionality reduction, and the top 30 principal components, based on statistical significance (P < 0.05), were selected (Figure S5d and e). A total of 9 cell populations were ultimately identified (Figure 8A). Cell annotation was performed through marker gene analysis, identifying 9 distinct cell types: T cell, T cell precursors, Spermatogenic cells, B cells, Epithelial cells, Monocytes, Neutrophils, Macrophages, and Erythrocytes (Table 5, Figure 8B). The marker genes CD79A, MS4A1, and CD19 exhibited high expression in B cells (Figure 8C). Furthermore, the key gene FOXJ3 showed significant and high expression in Spermatogenic cells, thus Spermatogenic cells were identified as the key cell type (Figure 8D). Cell communication analysis revealed that the spermatogenic cells had strong interactions with neutrophils, monocytes, and B cells (Figure 8E–H). The key gene FOXJ3 was expressed in spermatogenic cells, and the proposed time-series analysis indicated that the expression of FOXJ3 was low throughout the different differentiation stages of spermatogenic cells. (Figure S6a and b).
 |
Table 5 Cell Annotation Types and Their Corresponding Marker Genes
|
 |
Figure 8 Single-cell analysis results. (A) Results of UMAP cell clustering analysis. (B) Manually annotated cell types. (C) Bubble plot of marker gene expression in each cell cluster. (D) Expression of FOXJ3 in different cell clusters. (E) Number of interactions between different cell clusters. (F) Interaction strength between different cell clusters. (G) Heatmap of the number of interactions between different cell clusters. (H) Heatmap of interaction strength between different cell clusters.
|
Experimental Verification of Biomarkers in Sepsis
As shown in Figure 9, BRD1 was significantly decreased in the blood samples of sepsis patients, which was consistent with the bioinformatics analysis. However, although the expression trend of FOXJ3 was consistent with the bioinformatics analysis, there was no significant difference in the expression levels between sepsis patients and control blood samples.
 |
Figure 9 Experimental verification of biomarkers in sepsis. * P < 0.05 and Not Significant (ns): P > 0.05.
|
Discussion
Sepsis is a hyper-heterogeneous syndrome, accompanied by systemic inflammatory responses throughout the entire course of the disease. Moreover, the inflammatory and immune responses change dynamically at different pathogenic stages.42 Parthanatos, a form of cell death distinct from apoptosis and necrosis, has been implicated in various diseases, including cancer and neurodegenerative disorders.11 However, few studies exist between sepsis and Parthanatos. In this study, we identified six candidate signature genes-BRD1, BTN2A1, MSL2, GCC1, NCOA6, and FOXJ3-via MR analysis. Afterwards, through machine learning, expression validation and ROC analysis, BRD1 and FOXJ3 were retained as biomarkers for further research. We also explored their roles in immune infiltration, molecular regulatory networks, GSEA, chromosomal distribution, and drug prediction.
Bromodomain protein 1 (BRD1) is a gene associated with several psychiatric disorders.43 It encodes a scaffold protein that interacts with epigenetic modifiers involved in histone acetylation and histone H3 N-tail cleavage, regulating the brain development and in a number of biological processes including cell proliferation, differentiation and development are also crucial.44 BRD1’s association with mental disorders, including depression, epilepsy, and schizophrenia, has been supported by various genetic studies.45 Kerstin Klein et al discovered that the inhibition of BRD1 can reduce the expression of lipopolysaccharide-induced tumor necrosis factor-α (TNF-α), suppress the proliferation of synovial fibroblasts, and exert beneficial effects on rheumatoid arthritis (RA). Thus, BRD1 is regarded as a potential therapeutic target for RA.46 Recent evidence indicates that silencing the BRD1 gene in rheumatoid arthritis (RA) patient-derived macrophages (MDMs) exerts a subtle anti-inflammatory effect. Specifically, BRD1 knockdown reduced TNF-α mRNA expression upon LPS stimulation but did not alter TNF-α-induced TNF-α mRNA expression. This suggests that BRD1-mediated regulation of TNF-α may be stimulus-dependent.46 These findings suggest that BRD1 may modulate the initiation or suppression of pyroptosis in immune cells by regulating key inflammatory mediators such as TNF-α, with effects contingent upon specific stimuli (eg, LPS vs TNF-α). We further propose that in the context of sepsis, BRD1 likely operates through distinct inflammatory signaling cascades. Notably, during bacterial infections (exemplified by LPS-induced inflammation), BRD1 may attenuate maladaptive inflammatory responses by downregulating TNF-α expression.Both sepsis and rheumatoid arthritis (RA) are inflammatory diseases. It is also possible that BRD1 could serve as a therapeutic target for sepsis by reducing the expression of inflammatory factors. Functionally, BRD1 plays a role in mitochondrial bioenergetics, impacting mitochondrial function.47 It is worth noting that, according to reports, multiple factors during sepsis can cause mitochondrial damage. The damaged mitochondria actively participate in the formation of the inflammatory environment through key signaling pathways (such as Toll-like receptors), which exacerbates the occurrence of sepsis.48 Therefore, it is speculated that BRD1 can affect mitochondria and thus be involved in the occurrence and development of sepsis. Mitochondrial-targeted intervention may be a prospective target for the treatment of sepsis. To date, no association between BRD1 and sepsis has been reported in the literature, and it has been found for the first time that we will continue to investigate its role in septic diseases and aim to elucidate its underlying mechanisms in future studies.
FOXJ3, a transcription factor belonging to the Forkhead box (FOX) family, is involved in processes such as cell proliferation, migration, and spermatogenesis. The Role of FOX Family Member FOXJ1 in Regulating T Cell Activation and AutoreactivityFOXJ1, a member of the FOX family, plays a role in regulating T cell activation and autoreactivity. The deficiency of this gene in T cells can lead to multi-organ systemic inflammation and an increase in NF-κB activity. Therefore, FOXJ1 may regulate the inflammatory response and prevent autoimmunity by antagonizing pro-inflammatory transcriptional activities.49 In addition, a cross-sectional study has revealed that seven single nucleotide polymorphisms (SNPs), including rs2455084, rs1393009, and rs7539485, in the FOXJ3 gene are significantly associated with the susceptibility to rheumatoid arthritis (RA). The polymorphisms of FOXJ3 are related to the diagnostic indicators reflecting the disease activity of RA.50 Although there is currently no direct evidence linking FOXJ3 to sepsis, previous studies have demonstrated that both FOXJ1 and FOXJ3 are closely associated with inflammatory responses. This suggests that FOXJ3 may participate in the pathogenesis and progression of sepsis by modulating immune responses, particularly through its effects on T-cell function. Notably, recent studies have revealed that FOXJ3 knockout ameliorates arthritis symptoms, further supporting its role in inflammatory regulation.51 Based on these findings, we hypothesize that FOXJ3 may attenuate excessive immune activation in sepsis by antagonizing pro-inflammatory transcription factors, thereby modulating systemic inflammatory responses. Although a direct association between FOXJ3 and sepsis remains to be established, these insights identify FOXJ3 as a potential therapeutic target and provide a novel perspective for investigating the immune mechanisms underlying sepsis pathogenesis.
The effectiveness, sensitivity, and specificity of biomarkers are typically evaluated using the area under the curve (AUC). Both BRD1 and FOXJ3 demonstrated excellent diagnostic performance with AUC values exceeding 0.8 in both training and validation cohorts.We identified BRD1 (AUC 0.946) as a significant diagnostic biomarker for sepsis, outperforming procalcitonin (PCT, AUC 0.79).52 Furthermore, FOXJ3 achieved an AUC of 0.816, surpassing the current gold standard for immunoparalysis detection – monocyte HLA-DR (AUC 0.73)53 expression.Given the current lack of direct multicenter comparative data for these novel targets, future prospective validation across different cohorts is warranted. Combined analysis with established biomarkers (eg, PCT, IL-6) could further confirm their diagnostic value.
Both BRD1 and FOXJ3 were found to be associated with pathways involved in extracellular matrix organization and mRNA splicing. The extracellular matrix is a complex network of macromolecules that play a critical role in cellular adhesion, migration, proliferation, and differentiation.54 During sepsis, inflammatory factors stimulate the degradation of the extracellular matrix, particularly collagen and elastin, leading to remodelling and loss of function mRNA splicing, the process by which introns are removed and exons are joined to form mature mRNA, is fundamental to gene expression regulation.55 Dysregulation of this process is linked to various genetic disorders, including cancers, neurodegenerative diseases, and muscular disorders.56 In sepsis, abnormal mRNA splicing not only affects the balance of the immune response, but may also exacerbate organ damage and influence the course of the disease An in-depth study of these signalling pathways could help with the development of new therapeutic strategies for the treatment of sepsis.
BRD1 exhibited a significant negative correlation with Macrophages M0 and Macrophages M1, while FOXJ3 showed a significant negative correlation with Macrophages M0. M0 macrophages represent an unpolarized naive state derived from bone marrow monocytes, possessing low-level antigen-presenting capacity and basal secretory activity of inflammatory mediators. Under microenvironmental stimulation, they can polarize into either M1 or M2 phenotypes, where M1 macrophages primarily mediate pro-inflammatory responses while M2 macrophages exert anti-inflammatory effects.51 During the early phase of sepsis, M0 macrophages rapidly differentiate into the M1 phenotype, triggering systemic inflammatory response syndrome (SIRS). In later stages, they may differentiate into M2 macrophages or remain in the M0 state, leading to immunoparalysis and increased infection risk.BRD1 may epigenetically suppress the TLR4/NF-κB pathway, thereby reducing the propensity of M0-to-M1 conversion, whereas FOXJ3 may negatively regulate IRF5 – a critical determinant of M1 polarization – to maintain macrophages in the M0 state. Together, these factors form a “molecular brake” for M0 homeostasis, whose dysregulation represents a key contributor to immune imbalance in sepsis.Targeted regulation of the macrophage polarization process may normalize the host immune response, thereby providing a new therapeutic approach for the treatment of sepsis.57 Macrophage polarization likely serves as a central target for immune modulation in sepsis. As the differentiation origin, M0 macrophages hold upstream regulatory value for intervention strategies. Compared to other inflammatory cells such as neutrophils, macrophages demonstrate superior plasticity and greater accessibility for both research and therapeutic purposes.Furthermore, Jamie R Privratsky et al have discovered that the macrophage-endothelial immunomodulatory axis is also a crucial regulatory axis for ameliorating acute lung injury in sepsis.58 Therefore, it is speculated that biomarkers (such as BRD1 and FOXJ3) may participate in the occurrence and development of sepsis by influencing the infiltration status of macrophages.Our findings demonstrate significant differences in regulatory T cell (Treg) abundance between the sepsis and control groups. Tregs play a crucial role in mitigating excessive inflammation during early-stage sepsis, while their upregulated frequency in late stages may induce immunosuppression.59 BRD1 potentially upregulates gene expression associated with immunosuppressive cells (including Tregs), leading to abnormal increases in both Treg quantity and function. These hyperplastic Tregs may suppress the activation and proliferation of effector T cells and natural killer cells, thereby impairing pathogen clearance capacity.The suppressive function of Tregs largely depends on high-level and stable expression of the transcription factor FOXP3.60 However, excessive FOXP3 activation may also cause abnormal Treg expansion, consequently weakening the host’s immune response against infections. Therefore, further investigation of BRD1 and FOXJ3 in immune regulation is warranted, which would help elucidate sepsis-related immune mechanisms and identify novel potential therapeutic targets for sepsis treatment.
Common Targets for BRD1 and FOXJ3: the drug prediction analysis identified digoxin, doxorubicin, and daunorubicin as common drugs targeting both BRD1 and FOXJ3. Digoxin, a cardiac glycoside derived from the foxglove plant, is primarily used to treat heart conditions and has potential preventive effects in osteoarthritis and intervertebral disc degeneration.61,62 Doxorubicin, an anthracycline antibiotic derived from streptomyces peucetius, has been widely used as a chemotherapeutic agent since the 1960s.63 It is reported that sialic acid-conjugate modified doxorubicin can treat neutrophil-related inflammation.64,65 However, studies have demonstrated that both digoxin and doxorubicin exhibit certain toxicity, particularly potentially inducing cardiotoxicity at excessive doses.66 Daunorubicin, another anthracycline antibiotic, is used to treat a variety of cancers, particularly leukemia. Similar to doxorubicin, it intercalates into DNA, inhibiting DNA and RNA synthesis to exert its antitumor effects.67 Given the therapeutic potential and safety concerns of these drugs in sepsis treatment, future studies should establish animal models to thoroughly investigate the mechanisms of action and therapeutic effects of digoxin and doxorubicin in sepsis. Through animal experiments, we can evaluate the pharmacological responses of these drugs during different sepsis stages, particularly their impacts on immune and inflammatory responses.The pathogenesis of sepsis involves complex and diverse mechanisms, and the prediction and development of related drugs still require further exploration in future research. As potential biomarkers, BRD1 and FOXJ3 could not only facilitate early sepsis diagnosis but also provide evidence for personalized treatment strategies.Precision diagnosis and drug prediction based on these biomarkers may offer novel therapeutic targets and approaches for sepsis patients in the future. This could lead to improved treatment outcomes, reduced medication risks, and advancement of personalized medicine.The pathogenesis of sepsis is complex and varied, and the prediction and development of its drugs still need to be explored in future research.
This study also has limitations in that the study was analyzed by MR, which has limitations, such as the validity and applicability of genetic tools and the possibility of genetic heterogeneity, which should be carefully considered when performing MR studies. Furthermore, the RT-qPCR results demonstrated that BRD1 was consistent with bioinformatics findings, confirming its potential as a biomarker for sepsis. Although FOXJ3 showed a similar expression trend, it did not exhibit significant differences between disease and healthy blood samples. This discrepancy may be due to the small sample size, batch effects, or differences in sample sources, necessitating further studies to explore its role in disease.
Conclusion
In summary, in this study, based on the public data download transcriptomics data and parthanatos-related genes, we highlight BRD1 and FOXJ3 as promising biomarkers for sepsis through a number of columns of bioinformatics methods. In addition, we performed pathway enrichment analysis, immune cell infiltration status and drug prediction of these genes, which provided new insights into the management of sepsis. Future studies will expand sample sizes and employ BRD1 knockout mouse models or RNA interference technology to investigate the role of BRD1 in programmed cell death within immune cells. Optimized disease models will be established to evaluate clinically relevant outcomes and compare them with in vitro experimental results, thereby validating the robustness of conclusions and laying the foundation for translational applications.The investigation will focus on elucidating the functions of these genes and related pathways to provide a theoretical basis for developing targeted therapies. Additionally, animal models will be utilized to assess the safety, toxicity, and optimal therapeutic windows of digoxin and doxorubicin, generating evidence for clinical applications. These studies will further verify the drugs’ specific effects on BRD1/FOXJ3, offering theoretical support for advancing targeted therapeutic strategies.Follow-up studies will further focus on the role of these genes and related pathways to provide a theoretical basis for the development of targeted therapies.
Abbreviations
PRGs, Parthanatos-related genes; DEGs, Differentially expressed genes; MR, Mendelian randomization; ROC, Receiver operating characteristic; RT-qPCR, Real-time quantitative polymerase chain reaction; NAD, Nicotinamide adenine dinucleotide; ATP, Adenosine triphosphate; AIF, Apoptosis inducing factor; SNPs, Single nucleotide polymorphisms; GEO, Gene expression omnibus; PRGs, Parthanatos-related genes; WGCNA, Weighted gene co-expression network analysis; GO, Gene ontology; KEGG, Kyoto encyclopedia of genes and genomes; PPI, Protein-protein interaction; IV, Instrumental variables; IVW, Inverse variance weighted; LASSO, Least absolute shrinkage and selectionoperator; AUC, Area under the curve; TFs, Transcription factors; MF, Molecular functions; BP, Biological processes; CC, Cellular components; BRD1, Bromodomain protein 1; TNF-α, Tumor necrosis factor-α; RA, Rheumatoid arthritis.
Data Sharing Statement
The original data presented in the study are openly available in [Gene Expression Omnibus] at [GSE65682, GSE95233, http://www.ncbi.nlm.nih.gov/geo/], [Integrative Epidemiology Unit (IEU) Open Genome-Wide Association Study] at [ieu-b-4980, https://gwas.mrcieu.ac.uk/].
Ethics Approval and Informed Consent
I certify that the research study titled Transcriptome Combined with Mendelian Randomization to Identify and Validate Biomarkers Associated with Parthanatos in Sepsis has been approved by the Shanxi Bethune Hospital Ethics Committee.The approval number and date of approval are as follows: YXLL-2023-269 and 2023.12.8.
Acknowledgments
We would like to express our sincere gratitude to all individuals and organizations who supported and assisted us throughout this research. Special thanks to the following authors: Suojuan Zhang and Jifang Liang. In conclusion, we extend our thanks to everyone who has supported and assisted us along the way. Without your support, this research would not have been possible.
Author Contributions
All authors made a significant contribution to the work reported, whether that is in the conception, study design, execution, acquisition of data, analysis and interpretation, or in all these areas; took part in drafting, revising or critically reviewing the article; gave final approval of the version to be published; have agreed on the journal to which the article has been submitted; and agree to be accountable for all aspects of the work.
Funding
This work was supported by Clinical Key Specialty Program of Shanxi Bethune Hospital and the Shanxi Provincial Department of Science and Technology,Shanxi Provincial Basic Research Program, Free Exploration Youth Scientific Research Project [grant number 202403021212186].
Disclosure
The authors declare no conflicts of interest in this work. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.
References
1. Piva S, Bertoni M, Gitti N, Rasulo FA, Latronico N. Neurological complications of sepsis. Curr Opin Crit Care. 2023;29(2):75–84. doi:10.1097/MCC.0000000000001022
2. Fleischmann-Struzek C, Mellhammar L, Rose N, et al. Incidence and mortality of hospital- and ICU-treated sepsis: results from an updated and expanded systematic review and meta-analysis. Intensive Care Med. 2020;46(8):1552–1562. doi:10.1007/s00134-020-06151-x
3. Rudd KE, Johnson SC, Agesa KM, et al. Global, regional, and national sepsis incidence and mortality, 1990–2017: analysis for the Global Burden of Disease Study. Lancet. 2020;395(10219):200–211. doi:10.1016/S0140-6736(19)32989-7
4. Qiu X, Lei YP, Zhou RX. SIRS, SOFA, qSOFA, and NEWS in the diagnosis of sepsis and prediction of adverse outcomes: a systematic review and meta-analysis. Expert Rev Anti Infect Ther. 2023;21(8):891–900. doi:10.1080/14787210.2023.2237192
5. Llitjos JF, Carrol ED, Osuchowski MF, et al. Enhancing sepsis biomarker development: key considerations from public and private perspectives. Crit Care. 2024;28(1):238. doi:10.1186/s13054-024-05032-9
6. Zampieri FG, Bagshaw SM, Semler MW. Fluid therapy for critically ill adults with sepsis: a review. JAMA. 2023;329(22):1967–1980. doi:10.1001/jama.2023.7560
7. Mirijello A, Tosoni A, Group OBOTIMSS. New strategies for treatment of sepsis. Medicina. 2020;56(10):527. doi:10.3390/medicina56100527
8. Huang P, Chen G, Jin W, Mao K, Wan H, He Y. Molecular mechanisms of parthanatos and its role in diverse diseases. Int J Mol Sci. 2022;23(13).
9. Koehler RC, Dawson VL, Dawson TM. Targeting parthanatos in ischemic stroke. Front Neurol. 2021;12:662034. doi:10.3389/fneur.2021.662034
10. Morana O, Wood W, Gregory CD. The apoptosis paradox in cancer. Int J Mol Sci. 2022;23(3):1328. doi:10.3390/ijms23031328
11. Elmore S. Apoptosis: a review of programmed cell death. Toxicol Pathol. 2007;35(4):495–516. doi:10.1080/01926230701320337
12. Lorente L, Martín MM, Ortiz-López R, et al. Parthanatos type programmed cell death and septic patient mortality. Medicina intensiva. 2023;47(12):691–696. doi:10.1016/j.medin.2023.04.016
13. Duan L, Xiao R, Liu S, Shi Y, Feng Y. Causality between cognitive performance and cardiovascular disease: a bidirectional Mendelian randomization study. Gene. 2024;891:147822. doi:10.1016/j.gene.2023.147822
14. Sekula P, M FDG, Pattaro C, Köttgen A. Mendelian randomization as an approach to assess causality using observational data. J Am Soc Nephrol. 2016;27(11):3253–3265. doi:10.1681/ASN.2016010098
15. Burgess S, Thompson SG, Collaboration CCG. Avoiding bias from weak instruments in Mendelian randomization studies. Int J Epidemiol. 2011;40(3):755–764. doi:10.1093/ije/dyr036
16. Fan S, Li H, Liu K. Molecular prognostic of nine parthanatos death-related genes in glioma, particularly in COL8A1 identification. J Neurochem. 2024;168(3):205–223. doi:10.1111/jnc.16049
17. Ritchie ME, Phipson B, Wu D, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47. doi:10.1093/nar/gkv007
18. Gustavsson EK, Zhang D, Reynolds RH, Garcia-Ruiz S, Ryten M. ggtranscript: an R package for the visualization and interpretation of transcript isoforms using ggplot2. Bioinformatics. 2022;38(15):3844–3846. doi:10.1093/bioinformatics/btac409
19. Gu Z, Hübschmann D. Make Interactive Complex Heatmaps in R. Bioinformatics. 2022;38(5):1460–1462. doi:10.1093/bioinformatics/btab806
20. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinf. 2008;9:559. doi:10.1186/1471-2105-9-559
21. Zheng Y, Gao W, Zhang Q, et al. Ferroptosis and autophagy-related genes in the pathogenesis of ischemic cardiomyopathy. Front Cardiovasc Med. 2022;9:906753. doi:10.3389/fcvm.2022.906753
22. Yu G, Wang L-G, Han Y, He Q-Y. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics. 2012;16(5):284–287. doi:10.1089/omi.2011.0118
23. Liu P, Xu H, Shi Y, Deng L, Chen X. Potential molecular mechanisms of plantain in the treatment of gout and hyperuricemia based on network pharmacology. Evidence-Based Complementary Alternative Med. 2020;2020:3023127. doi:10.1155/2020/3023127
24. Hemani G, Zheng J, Elsworth B, et al. The MR-Base platform supports systematic causal inference across the human phenome. eLife. 2018;7.
25. Xian W, Wu D, Liu B, et al. Graves disease and inflammatory bowel disease: a bidirectional Mendelian randomization. J Clin Endocrinol Metab. 2023;108(5):1075–1083. doi:10.1210/clinem/dgac683
26. Liu Z, Wang H, Yang Z, Lu Y, Zou C. Causal associations between type 1 diabetes mellitus and cardiovascular diseases: a Mendelian randomization study. Cardiovasc Diabetol. 2023;22(1):236. doi:10.1186/s12933-023-01974-6
27. Mounier N, Kutalik Z. Bias correction for inverse variance weighting Mendelian randomization. Genetic Epidemiol. 2023;47(4):314–331. doi:10.1002/gepi.22522
28. Xu J, Zhang S, Tian Y, et al. Genetic Causal Association between iron status and osteoarthritis: a two-sample Mendelian Randomization. Nutrients. 2022;14(18):3683. doi:10.3390/nu14183683
29. Dai Z, Xu W, Ding R, et al. Two-sample Mendelian randomization analysis evaluates causal associations between inflammatory bowel disease and osteoporosis. Front Public Health. 2023;11:1151837. doi:10.3389/fpubh.2023.1151837
30. Verbanck M, Chen C-Y, Neale B, Do R. Detection of widespread horizontal pleiotropy in causal relationships inferred from Mendelian randomization between complex traits and diseases. Nat Genet. 2018;50(5):693–698. doi:10.1038/s41588-018-0099-7
31. Li Y, Lu F, Yin Y. Applying logistic LASSO regression for the diagnosis of atypical Crohn’s disease. Sci Rep. 2022;12(1):11340. doi:10.1038/s41598-022-15609-5
32. Yue S, Li S, Huang X, et al. Machine learning for the prediction of acute kidney injury in patients with sepsis. J Transl Med. 2022;20(1):215. doi:10.1186/s12967-022-03364-0
33. Robin X, Turck N, Hainard A, et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinf. 2011;12:77. doi:10.1186/1471-2105-12-77
34. Sachs MC. plotROC: a tool for plotting ROC curves. J Statistical Software. 2017;79. doi:10.18637/jss.v079.c02
35. Vickers AJ, Calster B, Steyerberg EW. A simple, step-by-step guide to interpreting decision curve analysis. Diagn Progn Res. 2019;3:18. doi:10.1186/s41512-019-0064-7
36. Hu Y, Yan C, Hsu C-H, et al. OmicCircos: a simple-to-use R package for the circular visualization of multidimensional omics data. Cancer Inf. 2014;13:13–20. doi:10.4137/CIN.S13495
37. Kuleshov MV, Jones MR, Rouillard AD, et al. Enrichr: a comprehensive gene set enrichment analysis web server 2016 update. Nucleic Acids Res. 2016;44(W1):W90–7. doi:10.1093/nar/gkw377
38. Hao Y, Hao S, Andersen-Nissen E, et al. Integrated analysis of multimodal single-cell data. Cell. 2021;184(13):3573–87.e29. doi:10.1016/j.cell.2021.04.048
39. Li Y, Yu J, Li R, Zhou H, Chang X. New insights into the role of mitochondrial metabolic dysregulation and immune infiltration in septic cardiomyopathy by integrated bioinformatics analysis and experimental validation. Cell Mol Biol Lett. 2024;29(1):21. doi:10.1186/s11658-024-00536-2
40. Jin S, Guerrero-Juarez CF, Zhang L, et al. Inference and analysis of cell-cell communication using CellChat. Nat Commun. 2021;12(1):1088. doi:10.1038/s41467-021-21246-9
41. Qiu X, Hill A, Packer J, Lin D, Ma YA, Trapnell C. Single-cell mRNA quantification and differential analysis with Census. Nat Methods. 2017;14(3):309–315. doi:10.1038/nmeth.4150
42. Rubin J, Berra L. Electrical impedance tomography in the adult intensive care unit: clinical applications and future directions. Curr Opin Crit Care. 2022;28(3):292–301. doi:10.1097/MCC.0000000000000936
43. Rajkumar AP, Qvist P, Donskov JG, et al. Reduced Brd1 expression leads to reversible depression-like behaviors and gene-expression changes in female mice. Transl Psychiatry. 2020;10(1):239. doi:10.1038/s41398-020-00914-2
44. Paternoster V, Edhager AV, Qvist P, et al. Inactivation of the Schizophrenia-associated BRD1 gene in brain causes failure-to-thrive, seizure susceptibility and abnormal histone H3 acetylation and N-tail clipping. Mol Neurobiol. 2021;58(9):4495–4505. doi:10.1007/s12035-021-02432-8
45. Paternoster V, Svanborg M, Edhager AV, et al. Brain proteome changes in female Brd1± mice unmask dendritic spine pathology and show enrichment for schizophrenia risk. Neurobiol Dis. 2019;124:479–488. doi:10.1016/j.nbd.2018.12.011
46. Klein K, Kato M, Frank-Bertoncelj M, et al. Evaluating the bromodomain protein BRD1 as a therapeutic target in rheumatoid arthritis. Sci Rep. 2018;8(1):11125. doi:10.1038/s41598-018-29127-w
47. Paternoster V, Cömert C, Kirk LS, et al. The psychiatric risk gene BRD1 modulates mitochondrial bioenergetics by transcriptional regulation. Transl Psychiatry. 2022;12(1):319. doi:10.1038/s41398-022-02053-2
48. Hu D, Sheeja Prabhakaran H, Zhang YY, Luo G, He W, Liou YC. Mitochondrial dysfunction in sepsis: mechanisms and therapeutic perspectives. Crit Care. 2024;28(1):292. doi:10.1186/s13054-024-05069-w
49. Lin L, Spoor MS, Gerth AJ, Brody SL, Peng SL. Modulation of Th1 activation and inflammation by the NF-kappaB repressor Foxj1. Science. 2004;303(5660):1017–20.
50. Ban JY, Park HJ, Kim SK, et al. Association of forkhead box J3 (FOXJ3) polymorphisms with rheumatoid arthritis. Mol Med Rep. 2013;8(4):1235–1241. doi:10.3892/mmr.2013.1623
51. Yuan L, Jiang N, Li Y, Wang X, Wang W. RGS1 enhancer RNA promotes gene transcription by recruiting transcription factor FOXJ3 and facilitates osteoclastogenesis through PLC-IP3R-dependent Ca(2+) response in rheumatoid arthritis. Inflammation. 2025;48(1):447–463. doi:10.1007/s10753-024-02067-6
52. Ito A, Ishida T, Tachibana H, Arita M, Yamazaki A, Washio Y. Utility of procalcitonin for differentiating cryptogenic organising pneumonia from community-acquired pneumonia. Clin Chem Lab Med. 2019;57(10):1632–1637.
53. Mosayebi G, Alizadeh SA, Alasti A, et al. Is CD19 an immunological diagnostic marker for acute appendicitis? Iran J Immunol. 2013;10(4):216–228.
54. Karamanos NK, Theocharis AD, Piperigkou Z, et al. A guide to the composition and functions of the extracellular matrix. FEBS J. 2021;288(24):6850–6912. doi:10.1111/febs.15776
55. Rogalska ME, Vivori C, Valcárcel J. Regulation of pre-mRNA splicing: roles in physiology and disease, and therapeutic prospects. Nat Rev Genet. 2023;24(4):251–269. doi:10.1038/s41576-022-00556-8
56. Nikom D, Zheng S. Alternative splicing in neurodegenerative disease and the promise of RNA therapies. Nat Rev Neurosci. 2023;24(8):457–473. doi:10.1038/s41583-023-00717-6
57. Chen X, Liu Y, Gao Y, Shou S, Chai Y. The roles of macrophage polarization in the host immune response to sepsis. Int Immunopharmacol. 2021;96:107791. doi:10.1016/j.intimp.2021.107791
58. Privratsky JR, Ide S, Chen Y, et al. A macrophage-endothelial immunoregulatory axis ameliorates septic acute kidney injury. Kidney Int. 2023;103(3):514–528. doi:10.1016/j.kint.2022.10.008
59. Huang S, Liu D, Han L, et al. Decoding the potential role of regulatory T cells in sepsis-induced immunosuppression. Eur J Immunol. 2024;54(5):e2350730. doi:10.1002/eji.202350730
60. Vent-Schmidt J, Han JM, MacDonald KG, Levings MK. The role of FOXP3 in regulating immune responses. Int Rev Immunol. 2014;33(2):110–128. doi:10.3109/08830185.2013.811657
61. Wang K-D, Ding X, Jiang N, et al. Digoxin targets low density lipoprotein receptor-related protein 4 and protects against osteoarthritis. Ann Rheumatic Dis. 2022;81(4):544–555. doi:10.1136/annrheumdis-2021-221380
62. Meng Q, Liu K, Liu Z, et al. Digoxin protects against intervertebral disc degeneration via TNF/NF-κB and LRP4 signaling. Front Immunol. 2023;14:1251517. doi:10.3389/fimmu.2023.1251517
63. Carvalho C, Santos RX, Cardoso S, et al. Doxorubicin: the good, the bad and the ugly effect. Curr Med Chem. 2009;16(25):3267–3285. doi:10.2174/092986709788803312
64. Wang S, Lai X, Li C, et al. Sialic acid-conjugate modified doxorubicin nanoplatform for treating neutrophil-related inflammation. J Control Release. 2021;337:612–627. doi:10.1016/j.jconrel.2021.07.044
65. Peter S, Alven S, Maseko RB, Aderibigbe BA. Doxorubicin-based hybrid compounds as potential anticancer agents: a review. Molecules. 2022;27(14):4478. doi:10.3390/molecules27144478
66. Patocka J, Nepovimova E, Wu W, Kuca K. Digoxin: pharmacology and toxicology-A review. Environ Toxicol Pharmacol. 2020;79:103400. doi:10.1016/j.etap.2020.103400
67. Chiou J-T, Huang C-H, Wu T-H, et al. CREB/Sp1-mediated MCL1 expression and NFκB-mediated ABCB1 expression modulate the cytotoxicity of daunorubicin in chronic myeloid leukemia cells. Toxicol Appl Pharmacol. 2022;435:115847. doi:10.1016/j.taap.2021.115847