Hostname: page-component-848d4c4894-5nwft Total loading time: 0 Render date: 2024-05-15T20:31:35.963Z Has data issue: false hasContentIssue false

Identification of Potential Biomarkers of Septic Shock Based on Pathway and Transcriptome Analyses of Immune-Related Genes

Published online by Cambridge University Press:  01 January 2024

Jie Wang
Affiliation:
Department of Critical Care Medicine, The First Affiliated Hospital of Gannan Medical University, Ganzhou, Jiangxi 341000, China
Jie Cai
Affiliation:
Department of Critical Care Medicine, HUST Union Shenzhen Hospital (Nanshan Hospital), Shenzhen, Guangdong 518052, China
Linlin Yue
Affiliation:
Department of Critical Care Medicine, The First Affiliated Hospital of Gannan Medical University, Ganzhou, Jiangxi 341000, China
Xixi Zhou
Affiliation:
Department of Critical Care Medicine, The First Affiliated Hospital of Gannan Medical University, Ganzhou, Jiangxi 341000, China
Chunlin Hu
Affiliation:
Department of Emergency Medicine, The First Affiliated Hospital of Sun Yat-sen University, Guangzhou, 510080 Guangdong, China
Hongquan Zhu*
Affiliation:
Department of Critical Care Medicine, The First Affiliated Hospital of Gannan Medical University, Ganzhou, Jiangxi 341000, China
*
Correspondence should be addressed to Hongquan Zhu; zhq9712@gmu.edu.cn
Rights & Permissions [Opens in a new window]

Abstract

Immunoregulation is crucial to septic shock (SS) but has not been clearly explained. Our aim was to explore potential biomarkers for SS by pathway and transcriptional analyses of immune-related genes to improve early detection. GSE57065 and GSE95233 microarray data were used to screen differentially expressed genes (DEGs) in SS. Gene Ontology and KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway enrichment analyses of DEGs were performed, and correlations between immune cell and pathway enrichment scores were analyzed. The predictive value of candidate genes was evaluated by receiver operating characteristic (ROC) curves. GSE66099, GSE4607, and GSE13904 datasets were used for external validation. Blood samples from six patients and six controls were collected for validation by qRT-PCR and western blotting. In total, 550 DEGs in SS were identified; these genes were involved in the immune response, inflammation, and infection. Immune-related pathways and levels of infiltration of CD4 + TCM, CD8 + T cells, and preadipocytes differed between SS cases and controls. Seventeen genes were identified as potential biomarkers of SS (areas under ROC curves >0.9). The downregulation of CD8A, CD247, CD3G, LCK, and HLA-DRA in SS was experimentally confirmed. We identified several immune-related biomarkers in SS that may improve early identification of disease risk.

Type
Research Article
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Copyright
Copyright © 2023 Jie Wang et al.

1. Introduction

Septic shock (SS) is defined as an infection-related circulatory dysfunction and metabolic disorder and is clinically associated with myocardial dysfunction and decreased ejection fraction [Reference Kattan, Castro, Vera and Hernández1, Reference Spaggiari, Passini and Crestani2]. Patient signs and symptoms include fever, rapid heart rate, shortness of breath, weakness and sweating, hypoxia, and altered mental status [Reference Basodan, Al Mehmadi and Al Mehmadi3]. Kidney failure, malignancy, diabetes, chronic lung disease, congestive heart failure, and immunosuppression may also increase the risk of SS [Reference De La Rica, Gilsanz and Maseda4]. Currently, SS is the leading cause of death among hospitalized patients, with a fatality rate of up to 40% [Reference Cao, Wang and Wang5, Reference Font, Thyagarajan and Khanna6]. Although tissue perfusion can be quickly regained in patients with SS after a positive fluid resuscitation and symptomatic treatment with vasoactive drugs and anti-infection agents, the risk of recurrence after hospital discharge is high [Reference Ammar, Ammar and Wieruszewski7]. Furthermore, more than half of the survivors have impaired physical or neurocognitive function, mood disorders, and poor quality of life [Reference Angus and van der Poll8, Reference Thompson, Venkatesh and Finfer9]. To reduce the mortality rate and improve the quality of life of patients after SS resuscitation, clinical studies are ongoing; however, early untreatable multiorgan failure makes the recovery process difficult. Therefore, in addition to improving programmed management, increasing the specificity of early recognition can also effectively improve survival in SS.

Innate immunity and adaptive immunity play a key role in the response to SS. During the induced inflammatory response, monocytes, macrophages, and neutrophils are activated, exacerbating vascular damage by producing cytokines, proteases, kinases, and reactive oxygen species [Reference Mahapatra and Heffner10]. Immunosuppressive responses are active in patients with SS [Reference Beltrán-García, Osca-Verdegal and Jávega11], and the apoptosis of B cells and follicular dendritic cells is involved in this process [Reference Carson, Cavassani, Dou and Kunkel12]. T lymphocytes contribute markedly to the immune system, and T-cell abnormalities have been found in patients with SS [Reference Monserrat, de Pablo and Reyes13]. Furthermore, immune dysfunction has been shown to impair the ability to clear primary infections and increase secondary infections [Reference Rimmelé, Payen and Cantaluppi14, Reference Patil, Guo, Luan and Sherwood15]. Tolsma et al. found that neutropenia and specific immunodeficiency are independently associated with an increased risk of death in patients with SS [Reference Tolsma, Schwebel and Azoulay16]. However, the mechanism by which immune dysfunction contributes to the progression of SS remains unclear despite the important implications for the development of biomarkers.

In this study, we analyzed the immune-related pathways involved in SS progression and identified candidate biomarkers. First, we screened differentially expressed genes (DEGs) between patients with SS and healthy controls based on gene expression data from public databases, followed by enrichment analyses of DEGs. An immune cell deconvolution analysis and pathway-immune cell correlation analysis were performed to identify key immune-related genes. Finally, the expression of key genes was validated by quantitative real-time polymerase chain reaction (qRT-PCR) and western blotting. A flowchart of the study is displayed in Figure S1. Our findings are expected to contribute to the development of personalized immunotherapy regimens for patients with SS.

2. Materials and Methods

2.1. Data Acquisition

Two microarray datasets, GSE57065 and GSE95233, were downloaded from the Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo/) [Reference Barrett, Troup and Wilhite17]. GSE57065 comprised data for 28 patients whose blood samples were collected within 30 min, 24 h, and 48 h after SS and 25 healthy controls. The GSE95233 included 22 healthy controls and 51 patients with SS, sampled twice at admission and once at D2 or D3. All samples from the GSE57065 and GSE95233 datasets were detected using the GPL570 [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array. GSE4607 and GSE13904 were used for the external expression validation of candidate genes between healthy controls and SS. GSE66099, which includes data for 30 systemic inflammatory response syndrome (SIRS) and 181 SS samples, was used to compare the expression differences between SS samples and the samples of the diseases with similar or related phenotypes to SS.

2.2. Data Preprocessing and Gene Annotation

Based on the two datasets, the probe expression matrix after normalization and log2 transformation was downloaded. Thereafter, the annotation file from the detection platform was obtained to match the gene symbol by the probe number inside; probes that did not match a gene symbol were removed. If different probes were mapped to the same gene symbol, the mean value of the probes was used as the final expression value for the gene.

2.3. Screening for DEGs between SS Samples and Controls

The limma package version 3.10.3 in R [Reference Smyth18] (http://www.bioconductor.org/packages/2.9/bioc/html/limma.html) was used to analyze differences in gene expression between SS samples and healthy controls. The Benjamini & Hochberg (BH) method [Reference Thissen, Steinberg and Kuang19] was used to calculate the adjusted p value. Genes with an adjusted p value <0.05 and |logfold change (FC)| > 1 were selected as DEGs.

2.4. Functional Enrichment Analyses of DEGs

DAVID version 6.8 [Reference Huang, Sherman and Lempicki20] (https://david-d.ncifcrf.gov/) was applied for Gene Ontology (GO) [Reference Ashburner, Ball and Blake21] enrichment analysis of DEGs according to three main categories, biological process (BP), cell component (CC), molecular function (MF), as well as Kyoto Encyclopedia of Genes and Genomes (KEGG) [Reference Kanehisa and Goto22] pathway enrichment analyses. Results with p < 0.05 and gene count >2 were considered statistically significant.

2.5. Protein–Protein Interaction (PPI) Network Construction

The intersecting DEGs between the GSE57065 and GSE95233 datasets were obtained to predict the PPI using STRING version 10.0 [Reference Franceschini, Szklarczyk and Frankild23] (http://www.string-db.org/). During this analysis, the species was set to Homo, and the highest confidence score was set to 0.9. Cytoscape version 3.4.0 [Reference Shannon, Markiel and Ozier24] (http://chianti.ucsd.edu/cytoscape-3.4.0/) was used to visualize the PPI network, and CytoNCA plug-in version 2.1.6 [Reference Tang, Li, Wang, Pan and Wu25] (http://apps.cytoscape.org/apps/cytonca) was used to analyze the node degree with parameter setting “without weight.”

2.6. Gene Set Enrichment Analysis (GSEA)

The R package clusterProfiler version 3.16.0 [Reference Yu, Wang, Han and He26] (http://bioconductor.org/packages/release/bioc/html/clusterProfiler.html) was used to perform GSEA. The KEGG pathway gene set in the molecular signature database [Reference Liberzon, Subramanian, Pinchback, Thorvaldsdóttir, Tamayo and Mesirov27] (MSigDB, http://software.broadinstitute.org/gsea/msigdb/index.jsp) was used as the background gene set. DEGs between SS samples and controls with adjusted p < 0.05 and |logFC| > 0.263 were selected as input genes and were sorted in descending order by logFC. BH was finally used to adjust p values, and pathways with adjusted p values less than 0.05 were considered statistically significant.

2.7. Differential Analysis of Functions and Pathways

By considering c2.cp.kegg.v7.1.symbols.gmt in MSigDB as the background gene set, the R package GSVA version 1.36.2 [Reference Hänzelmann, Castelo and Guinney28] (http://bioconductor.org/packages/release/bioc/html/GSVA.html) was used to calculate the enrichment scores for KEGG pathways in each sample. A score matrix was obtained. Significant KEGG pathway differences were identified using the limma package with threshold values of adjusted p value <0.05 and |logFC| > 1.

2.8. Immune and Stromal Cell Deconvolution Analysis

xCell [Reference Aran, Hu and Butte29] (https://xcell.ucsf.edu/) was used to estimate the enrichment scores for 64 types of immune and stromal cells with a threshold of p < 0.05. The Wilcoxon test was then used to compare enrichment scores between SS samples and controls, and significance was set at p < 0.05. The shared differentially enriched immune cells between GSE57065 and GSE95233 were obtained as candidates for further analyses.

2.9. Analysis of Correlations between Pathways and Immune Cells

By considering the intersection of results obtained by GSEA and GSVA, key pathways were selected. Spearman correlation coefficients for relationships between the pathway enrichment score and immune cell deconvolution score in each sample were computed. The intersecting relations with adjusted p value <0.05 and |r| > 0.4 for GSE57065 and GSE95233 were selected as significant cell-pathway interactions.

2.10. Evaluation of Key Genes

The intersection of pathway-related genes involved in cell-pathway interactions and genes in the PPI network with degrees over 10 was identified as key genes. The R package ggstatsplot version 0.5.0 was used to construct violin plots of key genes in different groups, while the R package plotROC version 2.2.1 was used to generate the diagnostic receiver operator characteristic (ROC) curve for each gene.

2.11. Blood Sample Collection

Whole blood samples and white blood cell samples from 12 participants (six patients with SS and six healthy controls) were collected for qRT-PCR and western blotting, respectively. Patients were diagnosed according to the criteria defined by the European Society of Intensive Care Medicine/Society of Critical Care Medicine. Patients with SS required a vasopressor to maintain a mean arterial pressure greater than 65 mmHg and elevated serum lactate greater than 2 mmol/L, despite adequate fluid resuscitation [Reference Font, Thyagarajan and Khanna6, Reference Esposito, De Simone, Boccia, De Caro and Pagliano30]. Patients with cancer, diabetes, autoimmune diseases, or a history of viral infection were excluded from the study. Healthy volunteers were enrolled from the physical examination center, and subjects with a history of major disease or infection were excluded. All subjects were informed of the aims and procedures of this study, and signed informed consent was obtained accordingly. This study was approved by the Ethics Committee of the First Affiliated Hospital of Gannan Medical University (No : LLSC-2021101302) and conformed with the Declaration of Helsinki.

2.12. qRT-PCR Analysis

The core genes in the PPI network and genes with larger fold change values were selected for qRT-PCR. Total RNA was extracted from whole blood samples using TRIzol reagent (9109; TaKaRa, Kusatsu). After reverse transcription, cDNA was quantified by real-time PCR with the ABI 7900HT FAST (Applied Biosystems) and Power SYBR Green PCR Master Mix Kit (A25742; Thermo, Waltham, MA, USA). GAPDH was used as an internal reference. The forward and reverse primer sequences are shown in Table S1. The reaction conditions were as follows: 50.0°C for 2 min, 95.0°C for 10 min, and 40 cycles at 95.0°C for 15 s and 60.0°C for 60 s. The relative expression level was detected in triplicate and was normalized using the 2−ΔΔct method.

2.13. Western Blot Analysis

Genes with the largest degree of connections in the PPI network were selected for western blot analysis. Antibodies against HLA-DRA (A11787) and LCK (A2177) were obtained from ABclonal (Wuhan, China), while the GAPDH antibody (60004-1-lg) was purchased from Proteintech (Rosemont, IL, USA). Western blot analysis was performed using standard procedures with anti-HLA-DRA (1 : 1000) and anti-LCK (1 : 1000). After 2 h of incubation with secondary antibodies (anti-rabbit IgG: 111-035-045, Jackson; anti-mouse IgG: 115-035-003, Jackson), the protein bands were developed via chemiluminescence using the Millipore ECL system (Billerica, MA, USA). Thereafter, the bands were scanned and recorded using Tanon Image (Tanon, Shanghai, China).

2.14. Statistical Analysis

Tanon Image was used to analyze the gray level of the western blot data. The qRT-PCR and western blot results were analyzed using GraphPad Prism 5 (GraphPad Software, San Diego, CA, USA) and visualized using histograms. Differences in the mRNA and protein expression levels of candidate genes between the SS and control groups were analyzed by t-tests. Statistics p < 0.05, p < 0.01, and p < 0.001 represent significant, highly significant, and extremely significant, respectively.

3. Results

3.1. Screening of DEGs between SS Samples and Healthy Controls

Based on the expression matrix of GSE57065 and GSE95233, principal component analysis (PCA) was performed (Figures 1(a) and 1(b)), which revealed a sharp distinction between the SS and control groups in both datasets. There were no significant outliers, and all samples could be used for further analyses. In total, 763 and 933 DEGs were obtained from GSE57065 and GSE95233, respectively, as shown in Figures 1(c) and 1(d). DEGs with |logFC| > 2 were selected to construct a heatmap (Figures 1(e) and 1(f)); these genes were identified as significantly differentially expressed between the SS and control groups.

FIGURE 1: Screening of DEGs between SS samples and controls. PCA diagrams of samples in GSE57065 (a) and GSE95233 (b) suggest that SS samples and control samples could be clearly distinguished, with no significant outliers. A total of 763 and 933 DEGs from GSE57065 (c) and GSE95233 (d) were detected. Red triangles and blue squares represent upregulated and downregulated DEGs, respectively. The top 10 upregulated and downregulated DEGs based on logFC values are labeled. Heatmaps showing the expression differences of DEGs with |logFC| > 2 between SS samples and controls in GSE57065 (e) and GSE95233 (f).

3.2. Function and Pathway Enrichment Analyses of DEGs

GO and KEGG pathway enrichment analyses of the DEGs in the two datasets were performed. DEGs in GSE57065 were mainly enriched for 124 GO-BP terms, 34 GO-CC terms, 30 GO-MF terms, and 32 KEGG pathways. DEGs in GSE95233 were mainly enriched for 127 GO-BP terms, 37 GO-CC terms, 41 GO-MF terms, and 34 KEGG pathways. Based on the ranking of p values, the top 10 GO functions for GSE57065 and GSE95233 are displayed in Figures 2(a) and 2(b). Among them, the MHC class II protein complex binding of GO-MF, T-cell receptor complex of GO-CC, and T-cell activation of GO-BP showed the highest fold enrichment in both datasets. The enriched KEGG pathways were broadly consistent across the two datasets (Figures 2(c) and 2(d)), thereby suggesting good homogeneity and reliability among GSE57065 and GSE95233.

FIGURE 2: Functional and pathway enrichment analyses of DEGs. Top 10 GO functions of DEGs in GSE57065 (a) and GSE95233 (b) in the BP, CC, and MF categories. Enriched KEGG pathways for GSE57065 (c) and GSE95233 (d). The x-axis indicates fold enrichment, and the y-axis indicates GO functions and KEGG pathways. The larger the bubble, the greater the fold enrichment. The smaller the bubble, the smaller the p value.

3.3. Construction of a PPI Network

The intersection of DEGs in GSE57065 and GSE95233 included 550 genes (Figure 3(a)). These DEGs were used to construct a PPI using STRING. A total of 1,104 interactions involving 243 proteins were obtained (Figure 3(b)). Nodes with more connections had larger contributions to the PPI network.

FIGURE 3: PPI network construction. (a) Venn diagram showing 550 intersecting DEGs in the GSE57065 and GSE95233 datasets. (b) PPI network involving 243 DEGs. Red and green circles represent upregulated and downregulated DEGs, respectively. The larger the node, the greater the degree. Gray lines represent protein-protein interactions. The darker the line, the larger the value of |logFC|.

3.4. GSEA

With GSEA, 5 upregulated and 13 downregulated pathways were obtained in GSE57065, whereas 7 upregulated and 15 downregulated pathways were obtained in GSE95233. According to the normalized enrichment scores (NES), the top six upregulated pathways, including Alzheimer’s disease, complement and coagulation cascades, oxidative phosphorylation, and Parkinson’s disease, and the top six downregulated pathways, such as graft versus host disease, primary immunodeficiency, and T-cell receptor signaling pathway, were selected in both datasets and are displayed in Figure 4.

FIGURE 4: GSEA. Top 6 upregulated and downregulated pathways in GSE57065 and GSE95233 according to NES ranking.

3.5. Differential Analysis of Enriched Pathways

The enrichment scores for KEGG pathways for each sample were calculated using the GSVA algorithm. Based on the differential analysis of enriched pathways, both GSE57065 and GSE95233 had 15 significantly downregulated KEGG pathways. As shown in heatmaps (Figures 5(a) and 5(b)), enrichment scores for these pathways differed significantly between the SS and control groups.

FIGURE 5: Differential analysis of enriched pathways in the two datasets. Heatmaps show KEGG pathways with significant differences in enrichment scores between SS and control samples in GSE57065 (a) and GSE95233 (b).

3.6. Deconvolution Analysis of Immune and Stromal Cells

To evaluate differences in immune cell enrichment between the two groups, we obtained the enrichment scores for 64 types of immune cells and stromal cells in each sample and compared SS and healthy controls by a Wilcoxon test (Figures 6(a) and 6(b)). In total, 11 and 14 types of cells had significantly different enrichment scores between the two groups in GSE57065 and GSE95233, respectively. By considering the intersection (Figure 6(c)), six types of cells, including megakaryocytes, CD8 + T cells, CD4 + TCM, preadipocytes, osteoblasts, and epithelial cells, differed significantly between SS samples and controls in both datasets.

FIGURE 6: Differential analysis of immune and stromal cell enrichment. Violin plots show the immune and stromal cells with significant differences in enrichment scores between SS samples and controls in GSE57065 (a) and GSE95233 (b). (c) Venn diagram of the intersection of cells (megakaryocytes, CD8 + T-cells, CD4 + TCM, preadipocytes, osteoblasts, and epithelial cells) with significant differences in relative abundance between SS samples and controls in both datasets.

3.7. Correlations between Pathways and Immune Cells

The intersection of significant pathways obtained by GSEA and GSVA included 12 KEGG pathways. Spearman correlation coefficients were obtained to evaluate the relationships among these 12 key pathways and frequencies of six significant cell types in each sample. Accordingly, 41 and 27 cell-pathway relationships were obtained in GSE57065 and GSE95233, respectively. Finally, 21 overlapped cell-pathway relationships with significant positive correlations were identified. These 21 cell-pathway relationships contained 10 key pathways and 3 key cells, including CD8 + T cells, CD4 + T cells, and preadipocytes; their relationships in GSE57065 and GSE95233 are summarized in Table 1.

TABLE 1: Correlations between KEGG pathway enrichment score and deconvolution scores for key cells in GSE57065 and GSE95233.

KEGG, Kyoto Encyclopedia of Genes and Genomes; CI, confidence interval. r indicates the Spearman correlation coefficient. p < 0.05 indicates statistical significance.

3.8. Predictive Performance of Key Genes

The intersection of genes involved in 10 key pathways and genes in the PPI network (with degrees >10), and 17 key genes (CD8A, HLA-DPA1, HLA-DPB1, HLA-DQA1, HLA-DQB1, HLA-DRA, ZAP70, MAPK14, CD247, CD3D, CD3E, CD3G, LCK, PRKCQ, ITK, LAT, and FYN) were obtained, as shown in Table 2. The differences in these 17 genes between SS samples and controls in GSE57065 and GSE95233 were analyzed, and ROC curves were compared to explore their diagnostic value in SS. All 17 genes were significantly differentially expressed in the two datasets (p < 0.05). Among them, 16 genes were downregulated, and MAPK14 was upregulated in SS. Furthermore, the ROC curves based on the two datasets suggested that these candidate genes had superior diagnostic values for SS, with area under the curve (AUC) values greater than 0.9. The violin plots and ROC curves for LCK and HLA-DRA, which had the highest degrees of connections in the PPI network, are shown in Figures 7(a) and 7(b).

TABLE 2: Key genes in enriched KEGG pathways.

KEGG, Kyoto Encyclopedia of Genes and Genomes.

FIGURE 7: Predictive performance of key genes in the cell-pathway relations. (a) Violin plots show that LCK and HLA-DRA, which had the highest in the PPI network, were significantly differentially expressed between SS samples and controls in GSE57065 and GSE95233. (b) ROC curves suggest that LCK and HLA-DRA have excellent abilities to diagnose SS with AUCs >0.9.

3.9. Verification with External Datasets

To verify the expression differences in the 17 candidate genes between healthy controls and SS samples, GSE4607 and GSE13904 datasets were used for external verification. In the GSE4607 dataset (Figure 8(a)), all genes showed significant differences in expression between the two groups, of which only MAPK14 was upregulated in SS. In the GSE13904 dataset (Figure 8(b)), MAPK14 was still overexpressed in SS, and the other 15 genes were significantly downregulated in the SS samples, except HLA−DQA1. The GSE66099 dataset was subsequently included to verify whether these genes have expression specificity in SS that enables their distinction from other types of inflammation. Expression levels of CD3E, CD3G, FYN, HLA-DPA1, HLA-DPB1, and HLA-DRA differed significantly between SIRS and SS (Figure 9), indicating that these six genes may contribute to SS inflammation.

FIGURE 8: External expression verification of 17 candidate genes using the GSE4607 (a) and GSE13904 (b) datasets.

FIGURE 9: Differences in the expression of candidate genes between SIRS and SS samples from GSE66099.

3.10. Experimental Expression Validation

A total of 12 whole blood samples (6 cases and 6 controls) were collected to validate the expression of key genes. Among the 17 key genes, those with the highest degrees of connections in the PPI network (LCK and HLA-DRA) and genes with the highest fold change values in the GSE57065 or GSE95233 datasets (CD8A, CD247, and CD3G) were selected for validation by qRT-PCR. As depicted in Figure 10(a), the mRNA expression levels of these five genes were significantly lower in SS samples than in the control group. As determined by western blotting (Figure 10(b)), the protein expression levels of LCK and HLA-DRA were significantly lower in patients with SS than in healthy controls.

FIGURE 10: Validation of expression patterns of key genes. (a) qRT-PCR results suggest that the mRNA expression levels of CD8A, CD247, CD3G, LCK, and HLA-DRA were lower in SS samples than in control samples. (b) Western blotting results revealed that the expression levels of LCK and HLA-DRA were downregulated in SS at the protein level. *p < 0.05; **p < 0.01.

4. Discussion

Sepsis and SS are life-threatening diseases caused by a dysregulated immune response to infection and may lead to tissue and organ damage and even death [Reference Rahmel31, Reference Mitchell, Ryan, Gillion, Wells and Muthiah32]. Currently, there is no effective treatment for SS. Accordingly, the disease burden can only be reduced by early detection, resuscitation, and the prompt administration of appropriate antibiotics [Reference Thompson, Venkatesh and Finfer9]. Unlike sepsis, SS leads to uncontrolled and intensified inflammatory responses, but the exact timing of triggering this process is elusive, which underscores the importance of identifying gene expression differences between SS and normal controls, rather than sepsis, for diagnosis in patients early in disease progression [Reference Cazalis, Lepape and Venet33]. Therefore, this study compared mRNA expression profiles between SS and normal samples and then identified 17 immune-related genes involved in the progression of SS by bioinformatics analyses. These genes also showed the potential to identify the risk of SS development in the validation cohorts. Five of these 17 genes (including CD8A, CD247, CD3G, LCK, and HLA-DRA) were experimentally validated for expression, considering their central roles in the PPI network and their consistent downregulation in SS samples were finally confirmed. The novel biomarkers proposed in this study are important to improve early identification and the management of acute episodes and to reduce septicemic deaths and disability. There are several similar articles that reported the diagnostic biomarkers for SS [Reference Fan, Han and Li34, Reference Hu, Cheng, Zhong, Chen and Zhang35], but our study differs in that we focus more on the gene set involved in SS-related immune regulation. Therefore, we also proposed three key immune cells including CD4 + T cells, CD8 + T cells, and preadipocytes that may be regulated by immune-related candidate genes and involved in disease progression in SS. These potential immune regulatory mechanisms are important clues for understanding the role of candidate genes in disease diagnosis and for developing new drugs to prevent SS.

In this study, DEGs between SS samples and controls were mainly enriched in biological functions and pathways related to immunological and inflammatory responses. Using the GSVA algorithm, significant differences in immune-related functions and pathways were also found between them, including the T-cell receptor signaling pathway, autoimmune disease, and primary immunodeficiency, among others. A related bioinformatics analysis supported our findings and demonstrated that DEGs in SS were involved in immune response, chemokine-mediated signaling, neutrophil chemotaxis, and chemokine activity [Reference Mohammed, Cui, Mas and Kamaleswaran36]. Additionally, immunodeficiency is commonly observed in patients with severe sepsis and SS and is associated with an increased risk of short-term mortality [Reference Tolsma, Schwebel and Azoulay16].

Considering the role of immunodeficiency in SS development, we carried out deconvolution and correlation analyses to determine the effect of immune cells on SS. Based on our results, three key immune cell types, i.e., CD4 + T cells, CD8 + T cells, and preadipocytes, differed significantly between SS samples and controls. In terms of adaptive immunity, sepsis-induced apoptosis leads to lymphocytopenia in patients with SS, and this process involves all types of T cells, including T regulatory cells, CD4 + T cells, CD8 + T cells, and natural killer cells, which are conducive to immunosuppression [Reference Rimmelé, Payen and Cantaluppi14]. Immunosuppression is a compensatory anti-inflammatory response that explains the short-term death of SS patients, while survivors may experience a prolonged state of immunosuppression, which could be reactivated by pathogenic infection [Reference Beltrán-García, Osca-Verdegal and Romá-Mateo37]. During the immunosuppression in SS, the loss of T-cell function is associated with reduced resistance to secondary infections in patients with SS [Reference Darcy, Minigo and Piera38]. Furthermore, decreased expression of cytotoxic molecules weakens the lytic activity of CD8 + T cells [Reference Jin, Liang and Liu39]. In this study, the enrichment scores for CD4 + T cells and CD8 + T cells were found to be significantly lower in SS cases than in controls. Roger et al. further supported our findings and proposed that rates of CD4 + T cell and CD8 + T cell apoptosis were higher in patients with SS than in controls [Reference Roger, Hyvernat and Breittmayer40]. The above evidence indicated that DEGs identified in this study may trigger immunosuppression and lead to SS by down-regulating CD4 + T cell and CD8 + T cell levels. In addition, we also proposed a relationship between preadipocytes and immune deficiency and T-cell receptor-related pathways. Several factors secreted by preadipocytes have pro-inflammatory and anti-inflammatory effects and can contribute to diseases associated with immune system dysfunction [Reference MacLaren, Cui and Cianflone41]. Some immunodeficiency virus protease inhibitors inhibit preadipocyte differentiation and promote adipocyte death [Reference Dowell, Flexner, Kwiterovich and Lane42]. In this study, we found several genes involved in the T-cell receptor signaling pathway, of which FYN is expressed in human preadipocytes and is induced after the initiation of differentiation [Reference Feng, Zhang and Shan43]. Therefore, we hypothesized that preadipocytes participate in the activation of immune-related pathways during SS development by inducing the expression of FYN; however, further mechanistic investigations are still needed.

Furthermore, 17 key genes were identified in immune-related cell-pathway pairs, and core genes from the PPI network were selected for expression validation. Relevant results suggested that the mRNA expression levels of CD8A, CD247, and CD3G were downregulated in SS samples, while LCK and HLA-DRA were decreased at both the mRNA and protein levels. As a member of the Src family of kinases, LCK is involved in changes in the activity of CD4 + T cells and CD8+ T cells during T-cell development [Reference Horkova, Drobek and Mueller44]. LCK plays a crucial role in T-cell differentiation, survival, and activation [Reference Zamoyska, Basson, Filby, Legname, Lovatt and Seddon45]; however, the contribution of LCK to the development of SS has not been explored. With regard to HLA-DRA, studies have found that the reduced expression of HLA-DR mRNA is correlated with increased mortality after SS [Reference Cazalis, Friggeri and Cavé46]. Winkler et al. reported that HLA-DR expression is decreased in sepsis [Reference Winkler, Rissiek and Priefler47], and the reduced HLA-DR expression may be a characteristic feature of septic monocytes [Reference Rimmelé, Payen and Cantaluppi14]. Furthermore, the dynamic changes in HLA-DRA gene expression and helper T cell subsets in patients with sepsis are indicative of immunosuppression [Reference Xu, Li and Xiao48]. Combined with the results of this study, we speculated that the loss of LCK and HLA-DRA expression may lead to the failure of T-cell differentiation and dynamic changes of helper T cell subsets, thus resulting in immunosuppression and the onset of SS.

The main limitation of this study is that patients with autoimmune diseases were excluded from the independent clinical validation cohort, which to some extent affects the extrapolation of the results. Furthermore, we only performed a preliminary exploration of the potential roles of these genes in disease. The regulatory mechanisms by which these candidate genes contribute to immunosuppression during SS are not clearly established. In future studies, bioinformatic analyses will be used to predict the upstream regulatory mechanisms, and experimental approaches will also be carried out to confirm the regulatory effects of the candidate genes. Additionally, due to limited sample size, we did not further screen the most robust biomarkers out from the 17 identified key genes by LASSO regression and/or multivariate Cox analyses. In the future, a most valuable predictive signature should be developed based on more convincing methods and larger sample size.

In conclusion, we found that DEGs between SS cases and controls were mainly enriched in immune- and inflammation-related functions and pathways. In addition, CD4 + T cells, CD8 + T cells, and preadipocytes were proposed as key immune cells to involve in the SS progression. These immune cells were also associated with 17 key immune-related genes, among which the downregulation of CD8A, CD247, CD3G, LCK, and HLA-DRA in SS samples was further experimentally validated. Our findings reveal several novel biomarkers for the early identification of SS.

Data Availability

Microarray datasets including GSE57065, GSE95233, GSE4607, GSE13904, and GSE66099 used in this study can be downloaded from the GEO database at http://www.ncbi.nlm.nih.gov/geo/. The raw data of qRT-PCR and western blot are available at https://doi.org/10.4121/19074482.

Conflicts of Interest

The authors declare that there are no conflicts of interest.

Authors’ Contributions

JW contributed conception and design of the research; JC and HZ acquired the data; LY and CH performed analysis and interpretation of data; XZ performed statistical analysis; HZ drafted the manuscript; CH and JW performed revision of the manuscript for important intellectual content. All authors read and approved the final manuscript.

Acknowledgments

This work was supported by the Science and Technology Program of the First Affiliated Hospital of Gannan Medical University (grant no. YJYB202131).

Supplementary Materials

Figure S1. Study flowchart. Table S1. The sequences of primers used for qRT-PCR. (Supplementary Materials)

References

Kattan, E., Castro, R., Vera, M., and Hernández, G., “Optimal target in septic shock resuscitation,Annals of Translational Medicine, vol. 8, no. 12, p. 789, 2020.CrossRefGoogle ScholarPubMed
Spaggiari, V., Passini, E., Crestani, S. et al., “Neonatal septic shock, a focus on first line interventions,Acta BioMedica Atenei Parmensis, vol. 93, no. 3, Article ID e2022141, 2022.Google ScholarPubMed
Basodan, N., Al Mehmadi, A. E., Al Mehmadi, A. E. et al., “Septic shock: management and outcomes,Cureus, vol. 14, no. 12, Article ID e32158, 2022.Google ScholarPubMed
De La Rica, A. S., Gilsanz, F., and Maseda, E., “Epidemiologic trends of sepsis in western countries,Annals of Translational Medicine, vol. 4, no. 17, p. 325, 2016.CrossRefGoogle Scholar
Cao, J. D., Wang, Z. C., Wang, Y. L. et al., “Risk factors for progression of Urolith Associated with Obstructive Urosepsis to severe sepsis or septic shock,BMC Urology, vol. 22, no. 1, p. 46, 2022.CrossRefGoogle ScholarPubMed
Font, M. D., Thyagarajan, B., and Khanna, A. K., “Sepsis and Septic Shock-basics of diagnosis pathophysiology and clinical decision making,Medical Clinics of North America, vol. 104, no. 4, pp. 573585, 2020.CrossRefGoogle ScholarPubMed
Ammar, M. A., Ammar, A. A., Wieruszewski, P. M. et al., “Timing of vasoactive agents and corticosteroid initiation in septic shock,Annals of Intensive Care, vol. 12, no. 1, p. 47, 2022.CrossRefGoogle ScholarPubMed
Angus, D. C. and van der Poll, T., “Severe sepsis and septic shock,New England Journal of Medicine, vol. 369, no. 9, pp. 840851, 2013.CrossRefGoogle ScholarPubMed
Thompson, K., Venkatesh, B., and Finfer, S., “Sepsis and septic shock: current approaches to management: sepsis and septic shock,Internal Medicine Journal, vol. 49, no. 2, pp. 160170, 2019.CrossRefGoogle ScholarPubMed
Mahapatra, S. and Heffner, A. C., Septic Shock. StatPearls. Treasure Island (FL), StatPearls Publishing LLC, Petersburg, FL, USA, 2022.Google Scholar
Beltrán-García, J., Osca-Verdegal, R., Jávega, B. et al., Characterization of early peripheral immune responses in patients with sepsis and septic shock,Biomedicines, vol. 10, no. 3, 2022.CrossRefGoogle ScholarPubMed
Carson, W. F., Cavassani, K. A., Dou, Y., and Kunkel, S. L., “Epigenetic regulation of immune cell functions during post-septic immunosuppression,Epigenetics, vol. 6, no. 3, pp. 273283, 2011.CrossRefGoogle ScholarPubMed
Monserrat, J., de Pablo, R., Reyes, E. et al., “Clinical relevance of the severe abnormalities of the T cell compartment in septic shock patients,Critical Care, vol. 13, no. 1, p. R26, 2009.CrossRefGoogle Scholar
Rimmelé, T., Payen, D., Cantaluppi, V. et al., “Immune cell phenotype and function in sepsis,Shock, vol. 45, no. 3, pp. 282291, 2016.CrossRefGoogle ScholarPubMed
Patil, N. K., Guo, Y., Luan, L., and Sherwood, E., “Targeting immune cell checkpoints during sepsis,International Journal of Molecular Sciences, vol. 18, no. 11, p. 2413, 2017.CrossRefGoogle ScholarPubMed
Tolsma, V., Schwebel, C., Azoulay, E. et al., “Sepsis severe or septic shock: outcome according to immune status and immunodeficiency profile,Chest, vol. 146, no. 5, pp. 12051213, 2014.CrossRefGoogle ScholarPubMed
Barrett, T., Troup, D. B., Wilhite, S. E. et al., “NCBI GEO: mining tens of millions of expression profiles--database and tools update,Nucleic Acids Research, vol. 35, pp. D760D765, 2007.CrossRefGoogle ScholarPubMed
Smyth, G. K., “Limma: linear models for microarray data,” in Bioinformatics and Computational Biology Solutions Using R and Bioconductor, Springer, Berlin, Germany, 2013.Google Scholar
Thissen, D., Steinberg, L., and Kuang, D., “Quick and easy implementation of the benjamini-hochberg procedure for controlling the false positive rate in multiple comparisons,Journal of Educational and Behavioral Statistics, vol. 27, no. 1, pp. 7783, 2002.CrossRefGoogle Scholar
Huang, D. W., Sherman, B. T., and Lempicki, R. A., “Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources,Nature Protocols, vol. 4, no. 1, pp. 4457, 2009.CrossRefGoogle ScholarPubMed
Ashburner, M., Ball, C. A., Blake, J. A. et al., “Gene Ontology: tool for the unification of biology,Nature Genetics, vol. 25, no. 1, pp. 2529, 2000.CrossRefGoogle ScholarPubMed
Kanehisa, M. and Goto, S., “KEGG: kyoto encyclopedia of genes and genomes,Nucleic Acids Research, vol. 28, no. 1, pp. 2730, 2000.CrossRefGoogle ScholarPubMed
Franceschini, A., Szklarczyk, D., Frankild, S. et al., “STRING v9.1: protein-protein interaction networks, with increased coverage and integration,Nucleic Acids Research, vol. 41, no. D1, pp. D808D815, 2012.CrossRefGoogle ScholarPubMed
Shannon, P., Markiel, A., Ozier, O. et al., “Cytoscape: a software environment for integrated models of biomolecular interaction networks,Genome Research, vol. 13, no. 11, pp. 24982504, 2003.CrossRefGoogle ScholarPubMed
Tang, Y., Li, M., Wang, J., Pan, Y., and Wu, F. X., “CytoNCA: a cytoscape plugin for centrality analysis and evaluation of protein interaction networks,Biosystems, vol. 127, pp. 6772, 2015.CrossRefGoogle ScholarPubMed
Yu, G., Wang, L. G., Han, Y., and He, Q. Y., “clusterProfiler: an R package for comparing biological themes among gene clusters,OMICS: A Journal of Integrative Biology, vol. 16, no. 5, pp. 284287, 2012.CrossRefGoogle Scholar
Liberzon, A., Subramanian, A., Pinchback, R., Thorvaldsdóttir, H., Tamayo, P., and Mesirov, J. P., “Molecular signatures database (MSigDB) 3.0,Bioinformatics, vol. 27, no. 12, pp. 17391740, 2011.CrossRefGoogle ScholarPubMed
Hänzelmann, S., Castelo, R., and Guinney, J., “GSVA: gene set variation analysis for microarray and RNA-seq data,BMC Bioinformatics, vol. 14, no. 1, p. 7, 2013.CrossRefGoogle ScholarPubMed
Aran, D., Hu, Z., and Butte, A. J., “xCell: digitally portraying the tissue cellular heterogeneity landscape,Genome Biology, vol. 18, no. 1, p. 220, 2017.CrossRefGoogle ScholarPubMed
Esposito, S., De Simone, G., Boccia, G., De Caro, F., and Pagliano, P., “Sepsis and septic shock: new definitions, new diagnostic and therapeutic approaches,Journal of Global Antimicrobial Resistance, vol. 10, pp. 204212, 2017.CrossRefGoogle ScholarPubMed
Rahmel, T., “[SSC international guideline 2016-management of sepsis and septic shock],Anasthesiol Intensivmed Notfallmed Schmerzther, vol. 53, no. 2, pp. 142148, 2018.Google ScholarPubMed
Mitchell, A. B., Ryan, T. E., Gillion, A. R., Wells, L. D., and Muthiah, M. P., “Vitamin C and thiamine for sepsis and septic shock,The American Journal of Medicine, vol. 133, no. 5, pp. 635638, 2020.CrossRefGoogle ScholarPubMed
Cazalis, M. A., Lepape, A., Venet, F. et al., “Early and dynamic changes in gene expression in septic shock patients: a genome-wide approach,Intensive care medicine experimental, vol. 2, no. 1, p. 20, 2014.CrossRefGoogle ScholarPubMed
Fan, Y., Han, Q., Li, J. et al., “Revealing potential diagnostic gene biomarkers of septic shock based on machine learning analysis,BMC Infectious Diseases, vol. 22, no. 1, 2022.CrossRefGoogle ScholarPubMed
Hu, Y., Cheng, L., Zhong, W., Chen, M., and Zhang, Q., “Bioinformatics analysis of gene expression profiles for risk prediction in patients with septic shock,Medical Science Monitor, vol. 25, pp. 95639571, 2019.CrossRefGoogle ScholarPubMed
Mohammed, A., Cui, Y., Mas, V. R., and Kamaleswaran, R., “Differential gene expression analysis reveals novel genes and pathways in pediatric septic shock patients,Scientific Reports, vol. 9, no. 1, Article ID 11270, 2019.CrossRefGoogle ScholarPubMed
Beltrán-García, J., Osca-Verdegal, R., Romá-Mateo, C. et al., “Epigenetic biomarkers for human sepsis and septic shock: insights from immunosuppression,Epigenomics, vol. 12, no. 7, pp. 617646, 2020.CrossRefGoogle ScholarPubMed
Darcy, C. J., Minigo, G., Piera, K. A. et al., “Neutrophils with myeloid derived suppressor function deplete arginine and constrain T cell function in septic shock patients,Critical Care, vol. 18, no. 4, p. R163, 2014.CrossRefGoogle ScholarPubMed
Jin, B., Liang, Y., Liu, Y. et al., “Notch signaling pathway regulates T cell dysfunction in septic patients,International Immunopharmacology, vol. 76, Article ID 105907, 2019.CrossRefGoogle ScholarPubMed
Roger, P. M., Hyvernat, H., Breittmayer, J. P. et al., “Enhanced T-cell apoptosis in human septic shock is associated with alteration of the costimulatory pathway,European Journal of Clinical Microbiology & Infectious Diseases, vol. 28, no. 6, pp. 575584, 2009.CrossRefGoogle ScholarPubMed
MacLaren, R., Cui, W., and Cianflone, K., “Adipokines and the immune system: an adipocentric view,Advances in Experimental Medicine and Biology, vol. 632, pp. 121, 2008.CrossRefGoogle ScholarPubMed
Dowell, P., Flexner, C., Kwiterovich, P. O., and Lane, M. D., “Suppression of preadipocyte differentiation and promotion of adipocyte death by HIV protease inhibitors,Journal of Biological Chemistry, vol. 275, no. 52, Article ID 41325, 2000.CrossRefGoogle ScholarPubMed
Feng, J., Zhang, X., Shan, C. et al., “Src family kinases involved in the differentiation of human preadipocytes,Molecular and Cellular Endocrinology, vol. 533, Article ID 111323, 2021.CrossRefGoogle ScholarPubMed
Horkova, V., Drobek, A., Mueller, D. et al., “Dynamics of the coreceptor-LCK interactions during T cell development shape the self-reactivity of peripheral CD4 and CD8 T cells,Cell Reports, vol. 30, no. 5, pp. 15041514.e7, 2020.CrossRefGoogle Scholar
Zamoyska, R., Basson, A., Filby, A., Legname, G., Lovatt, M., and Seddon, B., “The influence of the src-family kinases, Lck and Fyn, on T cell differentiation, survival and activation,Immunological Reviews, vol. 191, no. 1, pp. 107118, 2003.CrossRefGoogle Scholar
Cazalis, M. A., Friggeri, A., Cavé, L. et al., “Decreased HLA-DR antigen-associated invariant chain (CD74) mRNA expression predicts mortality after septic shock,Critical Care, vol. 17, no. 6, p. R287, 2013.CrossRefGoogle ScholarPubMed
Winkler, M. S., Rissiek, A., Priefler, M. et al., “Human leucocyte antigen (HLA-DR) gene expression is reduced in sepsis and correlates with impaired TNFα response: a diagnostic tool for immunosuppression?,PLoS One, vol. 12, no. 8, Article ID e0182427, 2017.CrossRefGoogle ScholarPubMed
Xu, J., Li, J., Xiao, K. et al., “Dynamic changes in human HLA-DRA gene expression and Th cell subsets in sepsis: indications of immunosuppression and associated outcomes,Scandinavian Journal of Immunology, vol. 91, no. 1, Article ID e12813, 2020.CrossRefGoogle ScholarPubMed
Figure 0

FIGURE 1: Screening of DEGs between SS samples and controls. PCA diagrams of samples in GSE57065 (a) and GSE95233 (b) suggest that SS samples and control samples could be clearly distinguished, with no significant outliers. A total of 763 and 933 DEGs from GSE57065 (c) and GSE95233 (d) were detected. Red triangles and blue squares represent upregulated and downregulated DEGs, respectively. The top 10 upregulated and downregulated DEGs based on logFC values are labeled. Heatmaps showing the expression differences of DEGs with |logFC| > 2 between SS samples and controls in GSE57065 (e) and GSE95233 (f).

Figure 1

FIGURE 2: Functional and pathway enrichment analyses of DEGs. Top 10 GO functions of DEGs in GSE57065 (a) and GSE95233 (b) in the BP, CC, and MF categories. Enriched KEGG pathways for GSE57065 (c) and GSE95233 (d). The x-axis indicates fold enrichment, and the y-axis indicates GO functions and KEGG pathways. The larger the bubble, the greater the fold enrichment. The smaller the bubble, the smaller the p value.

Figure 2

FIGURE 3: PPI network construction. (a) Venn diagram showing 550 intersecting DEGs in the GSE57065 and GSE95233 datasets. (b) PPI network involving 243 DEGs. Red and green circles represent upregulated and downregulated DEGs, respectively. The larger the node, the greater the degree. Gray lines represent protein-protein interactions. The darker the line, the larger the value of |logFC|.

Figure 3

FIGURE 4: GSEA. Top 6 upregulated and downregulated pathways in GSE57065 and GSE95233 according to NES ranking.

Figure 4

FIGURE 5: Differential analysis of enriched pathways in the two datasets. Heatmaps show KEGG pathways with significant differences in enrichment scores between SS and control samples in GSE57065 (a) and GSE95233 (b).

Figure 5

FIGURE 6: Differential analysis of immune and stromal cell enrichment. Violin plots show the immune and stromal cells with significant differences in enrichment scores between SS samples and controls in GSE57065 (a) and GSE95233 (b). (c) Venn diagram of the intersection of cells (megakaryocytes, CD8 + T-cells, CD4 + TCM, preadipocytes, osteoblasts, and epithelial cells) with significant differences in relative abundance between SS samples and controls in both datasets.

Figure 6

TABLE 1: Correlations between KEGG pathway enrichment score and deconvolution scores for key cells in GSE57065 and GSE95233.

Figure 7

TABLE 2: Key genes in enriched KEGG pathways.

Figure 8

FIGURE 7: Predictive performance of key genes in the cell-pathway relations. (a) Violin plots show that LCK and HLA-DRA, which had the highest in the PPI network, were significantly differentially expressed between SS samples and controls in GSE57065 and GSE95233. (b) ROC curves suggest that LCK and HLA-DRA have excellent abilities to diagnose SS with AUCs >0.9.

Figure 9

FIGURE 8: External expression verification of 17 candidate genes using the GSE4607 (a) and GSE13904 (b) datasets.

Figure 10

FIGURE 9: Differences in the expression of candidate genes between SIRS and SS samples from GSE66099.

Figure 11

FIGURE 10: Validation of expression patterns of key genes. (a) qRT-PCR results suggest that the mRNA expression levels of CD8A, CD247, CD3G, LCK, and HLA-DRA were lower in SS samples than in control samples. (b) Western blotting results revealed that the expression levels of LCK and HLA-DRA were downregulated in SS at the protein level. *p < 0.05; **p < 0.01.

Supplementary material: File

Wang et al. supplementary material
Download undefined(File)
File 180.4 KB