The microarray data set was screened from the GEO database (see text footnote 1) using the keywords “severe” and “viral infection “, resulting in the identification of GSE164805 (test set) and GSE199816 (validation set). The GSE164805 data set was generated using the GPL26963 platform, specifically the Agilent-085982 Arraystar Human lncRNA V5 microarray chip. Peripheral blood mononuclear cells (PBMCs) were collected from 5 patients with mild COVID-19, 5 patients with severe COVID-19, and 5 healthy controls. After standardization, annotation, and clearance of clinical information in dataset GSE164805, 5 samples from severe patients were selected as the coronavirus-associated viral sepsis group (test group) for bioinformatics analysis along with 5 samples from a healthy control group (control group). Dataset GSE199816 consisted of 9 patients with coronavirus-associated viral sepsis (test group) and 18 patients with bacterial sepsis (case control group). Ferroptosis-related genes including driver genes, suppressor genes, marker genes, and unclassified genes were obtained from FerrDb V2² within its database.

First, we preprocessed the raw file obtained from the GEO database by removing probes associated with multiple molecules. In cases where a probe was associated with the same molecule, we selected the average value of the probes as the gene expression value. Differential expression genes (DEGs) were analyzed between the coronavirus-associated viral sepsis group and control group using the limma package in R software. The DEGs condition was set as follows: p < 0.05, |logFC| > 1. We obtained 564 ferroptosis genes from FerrDb V2 database. We performed an intersection analysis between DEGs and ferroptosis genes and visualized the results using Venn diagrams created with VennDiagram R package, heat maps created with Heatmap package, and volcano maps created with ggplot2 package.

The gene ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) functional enrichment analyses were conducted on the aforementioned differentially expressed genes (DEGs) associated with ferroptosis resulting from coronavirus-associated viral sepsis using the clusterProfiler R package (8). Gene ontological (GO) enrichment analysis is a widely used approach for investigating functional enrichments, including biological process (BP), cell component (CC), and molecular function (MF), and plays a crucial role in exploring biological functions (9). The Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis is a network map that illustrates the relationships between individual genes or substances, aiming to explore potential signaling pathways (10).

The evaluation and analysis of protein–protein interaction (PPI) networks are utilized to elucidate crucial interactions among key cellular genes (11). In this study, a PPI network analysis was constructed for the differentially expressed genes (DEGs) associated with ferroptosis caused by viral sepsis. This analysis was performed using the online tool STRING,³ which integrates known and predicted associations between proteins, including physical interactions and functional associations. Subsequently, the results from STRING were imported into Cytoscape, and the cytoHubba plug-in was used to download interaction data and enhance the PPI network to extract critical subnetworks. Finally, based on the maximum correlation criterion (MCC algorithm), we identified the hub gene as the one with the highest score within these subnetworks.

This study utilized CIBERSORT (12) (Cell-type Identification By Estimating Relative Subsets Of RNA Transcripts), a deconvolution algorithm-based bioinformatics tool, to quantify the relative proportions of specific cell types from mixed tissue gene expression data. Through this methodology, we systematically evaluated immune cell infiltration in the samples, analyzing the distribution patterns of 22 immune cell subtypes (including CD8⁺ T cells, M1/M2 macrophages, neutrophils, and others). We further explored differences in immune cell marker expression between the coronavirus-associated viral sepsis group and healthy controls, thereby providing theoretical support for elucidating disease mechanisms and advancing precision therapeutics (13).

The experimental group consisted of 20 patients hospitalized from June 2022 to December 2023 with coronavirus-associated viral sepsis, while the control group included 20 healthy adults who underwent physical examination at the same period in the physical examination center. Inclusion criteria: (1) Experimental group: I. Individuals who met the diagnostic criteria for viral sepsis; II. Aged ≥18 years old; III. Provided informed consent. (2) Control group: I. Individuals whose physical examination indicated good health; II. Aged ≥18 years old; III. Provided informed consent. Exclusion criteria: (1) Age <18 years; (2) Patients who died on admission or were transferred to hospitals with incomplete clinical data; (3) End-stage patients with chronic and malignant diseases; (4) Individuals who refused to participate in the research or sign informed consent; (5) People with a history of mental illness or cognitive dysfunction. This study was approved by the Clinical Research Ethics Committee of the hospital. Peripheral venous blood samples were collected from each subject, centrifuged at 4°C, and then stored immediately in a refrigerator at −80°C for future use. Table 1 shows the basic characteristics of the patients recruited from the control group and the test group.

The total RNA was extracted from the serum of the subjects. Reverse transcription was performed under the following conditions: 11 μL of RNA (1 μg) was mixed with 1 μL of random primer (0.2 μg/mL) and incubated at 65°C for 5 min. Subsequently, a mixture containing 4 μL of 5 × Buffer, 3 μL of dNTP (10 mmol/L), 1 μL of RNase inhibitor (20 U/uL), and 1 μL of reverse transcriptase (20 U/uL) was added to the reaction. The reaction proceeded at 25°C for 10 min, followed by incubation at 42°C for 1 h and finally at 72°C for 15 min. The real-time PCR reaction system consisted of 10 μL FastStart Universal SYBR Green Master (ROX), 0.5 μL upstream primer (15 μM), 0.5 μL downstream primer (15 μM), 2 μL cDNA, and 7 μL DNAse and RNase-free water, making a total volume of 20 μL. The primer sequence for human β-actin was as follows: upstream primer -GTCAGGTCATCACTATCGGCAAT-3′ and downstream primer -AGAGGTCTTTACGGATGTCAACGT-3′. Additionally, the ferroptosis hub gene primer sequence was included. PCR reaction conditions: The PCR reaction was initiated by pre-denaturation at 94°C for 10 min to activate the Tag enzyme. Subsequently, the amplification process consisted of denaturation at 94°C for 15 s, annealing at 60°C for 60 s, and a total of 40 cycles. The cycle threshold (CT) value was determined using fluorescence quantitative PCR. Firstly, the difference between the CT value of the target gene (iNOS) and that of the internal control gene Beta-actin in each sample was calculated as ΔCT = CT(iNOS)-CT(beta-actin). Then, ΔΔCT was obtained by subtracting the ΔCT value of samples from the experimental group with that from the normal control group. Finally, calculation using 2^(−ΔΔCT) indicated the fold change in expression of iNOS gene in relation to its expression in samples from normal control group. This analysis aimed to evaluate expression levels of hub genes associated with coronavirus-associated viral sepsis.

The data analysis was conducted using R software (version 4.2.3). Measurement data that followed a normal distribution were presented as x ¯ ±s. An independent sample T-test was employed to compare the two groups, while non-normally distributed data were expressed as M(P25, P75), and the Mann–Whitney U test was utilized for group comparisons. A significance level of p < 0.05 was considered statistically significant. Receiver Operating Characteristic (ROC) curve analysis was performed to evaluate the diagnostic efficiency of ferroptosis-related hub genes through sensitivity and specificity metrics, with the Area Under the Curve (AUC) serving as the core evaluation metric. Interpretation criteria were defined as follows: an AUC of 0.5 indicates no discriminative capacity (equivalent to random chance); 0.6–0.8 represents moderate discriminative power; 0.8–0.9 reflects high diagnostic accuracy; and values exceeding 0.9 demonstrate exceptional predictive performance (14).

The Limma package was utilized to conduct differential analysis on the two groups of GSE164805 dataset, with a significance threshold set at p < 0.05 and |logFC| > 1. As a result, a total of 15,059 differentially expressed genes were identified. Among them, there were 7,458 up-regulated genes (represented by red dots in Figure 1A) and 7,601 down-regulated genes (represented by blue dots in Figure 1A). The findings were visualized using volcano plots (Figure 1A) and heat maps (Figure 1B).

The FerrDb V2 database provided a total of 564 ferroptosis genes. The intersection between the aforementioned coronavirus-associated viral sepsis DEGs and ferroptosis genes resulted in 14,870 differentially expressed genes associated with coronavirus-associated viral sepsis. Among them, 189 genes were identified as differentially expressed genes related to ferroptosis in coronavirus-associated viral sepsis (Supplementary Table 1). Additionally, a total of 375 genes were classified as ferroptosis genes (Figure 2).

The differentially expressed genes (DEGs) associated with ferroptosis resulting from coronavirus-associated viral sepsis were subjected to analysis using Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG). In the GO-BP analysis (Figure 3A), the primary pathway identified was oxidative stress. The GO-CC analysis (Figure 3B) revealed a significant enrichment of extracellular membrane-related processes. The results of the GO-MF enrichment analysis (Figure 3C) primarily indicated activation of pathways related to oxidative stress. KEGG analysis (Figure 3D) highlighted key signaling pathways including ferroptosis, HIF-1, PPAR, and NOD-like receptor signaling.

The Cibersort R software package was utilized to analyze immune cell infiltration in samples within the GSE164805 dataset. The box plot illustrating differences in immune cell populations revealed a higher presence of innate immune cells (such as dendritic cells, mast cells, and macrophages) in coronavirus-associated viral sepsis compared to the control group. Conversely, CD8⁺ T cells and NK cells exhibited lower levels of infiltration, as depicted in Figure 4. Furthermore, Figure 5 highlights a strong correlation between two hub genes, IL1-β and HMOX1, with the infiltration of CD8⁺ T cells and NK cells.

Among the 189 differentially expressed genes (DEGs) associated with ferroptosis in the context of coronavirus-associated viral sepsis, four of them are non-coding genes. These DEGs were imported into the STRING online database using default confidence settings (minimum interaction requirement score 0.4). Consequently, we obtained a protein–protein interaction (PPI) network comprising 185 nodes representing encoded proteins and 965 edges denoting interactions between these proteins. However, it is important to note that among these nodes, there were 11 proteins that remained unconnected to other molecules, that is, did not contribute to the formation of any molecular network (Supplementary Figure 1). To identify hub genes within this PPI network, we employed the Cytoscape CytoHubba plug-in and selected IL1-β and HMOX1 as the top two hub genes based on their highest scores according to the MCC algorithm (Figure 6).

The two screened ferroptosis hub genes IL1-β and HMOX1 were analyzed in the validation set (GSE199816). The results showed that compared to the bacterial sepsis group, the expression of ferroptosis hub genes IL1-β and HMOX1 was decreased in the coronavirus-associated viral sepsis group (Figure 7A). Further analysis of diagnostic efficacy showed that the AUROC of IL1-β and HMOX1 for the diagnosis of coronavirus-associated viral sepsis was 0.858 and 0.685, respectively (Figures 7B,C).

The expression levels of IL1-β and HMOX1 were assessed in the test group (n = 20) and the control group (n = 20) using RT-qPCR. The results demonstrated a significant down-regulation of IL1-β (0.390 ± 0.068) and HMOX1 (0.629 ± 0.117) expressions in patients with coronavirus-associated viral sepsis compared to the control group (1.101 ± 0.107), as depicted in Figures 8A,B, respectively, with statistical significance observed for both comparisons (both p < 0.05). Further analysis revealed that the AUROC values for diagnosing coronavirus-associated viral sepsis based on these two hub genes were 0.892 and 0.765, respectively, as shown in Figures 9A,B.