Identification and Validation of Autophagy-Related Genes in Sepsis-Induced Acute Respiratory Distress Syndrome and Immune Infiltration.

Purpose: Autophagy-related genes (ARGs) play an important role in the pathophysiology processes of sepsis-induced acute respiratory distress syndrome (ARDS). However, expression profiles of ARGs have rarely been used to explore the relationship between autophagy and sepsis-induced ARDS. Therefore, we aim to identify and validate the potential ARGs of sepsis-induced ARDS through bioinformatics analysis and experiment validation.
Methods: We downloaded GSE32707 data from the Gene Expression Omnibus (GEO) database. The potential differentially expressed genes (DEGs) and differentially expressed ARGs (DEARGs) of sepsis-induced ARDS were screened by R software. Then, we performed functional enrichment analyses to explore the potential biological functions of DEARGs and constructed protein-protein interaction (PPI) networks. Subsequently, correlation analysis and receiver operating characteristic (ROC) curve were used for the DEARGs. In addition, we estimated the proportions of 22 immune cell subsets by using CIBERSORT algorithm. Finally, RNA expression of seven DEARGs were validated by qRT-PCR in blood samples from sepsis-induced ARDS and healthy controls.
Results: We identified 28 DEARGs, including 11 up-regulated genes and 17 down-regulated genes, which were primarily involved in autophagy and apoptosis. Seven genes (BAG3, CTSD, ERBB2, MYC, PEA15, RAB24 and SIRT1) with AUC >0.70 were considered possible to be sepsis-induced ARDS hub genes for ROC curve analysis. CIBERSORT results shown that sepsis-induced ARDS contained a higher proportion of naive CD4+ T cells, gamma delta T cells, monocytes, and neutrophils, and lower levels of CD8+ T cells, memory resting CD4+ T cells, follicular helper T cells were relatively lower. The results of qRT-PCR also demonstrated that the expression levels of BAG3, CTSD, ERBB2, MYC and SIRT1 in sepsis-induced ARDS patients and healthy controls had differences.
Conclusion: We identified an association between DEGs and immune infiltration in sepsis-induced ARDS and validated BAG3, CTSD, ERBB2, MYC and SIRT1 that may be have excellent diagnostic performance.

Introduction

Acute respiratory distress syndrome (ARDS) is an acute inflammatory lung injury characterized by rapidly progressing respiratory failure, pulmonary edema, diffuse alveolar damage, and inflammatory cell infiltration.1 Sepsis is a common risk factor for ARDS and patient prognosis is worsened when sepsis is present.2 Previous studies suggested that multiple plasma biomarkers, such as interleukin-6 (IL-6), interleukin-8 (IL-8), angiopoietin-2 (ANG-2), Kerbs von Lungren 6 antigen (KL-6), receptor for advanced glycation endproducts (RAGEs), and surfactant protein D (SP-D), had predictive value for ARDS prognosis.3–5 Accumulating evidence showed that several biological functions also are involved in the pathogenesis of sepsis-induced ARDS, including cell proliferation, apoptosis and autophagy.6–8 Among these biological functions, autophagy plays an essential roles in the development of sepsis-induced ARDS.

Autophagy (“self-eating”) is the process through which parts of the cell are degraded in the lysosome.9 In recent decades, there have been some signaling pathways, which have been reported to affect the biological functions of ARDS through autophagy. For instance, one previous study demonstrated that autophagy attenuates pulmonary inflammation, lung permeability and hypoxemia by downregulating NLR family pyrin domain containing 3 (NLRP3) and interleukin-1β (IL-1β) inflammatory vesicle activity in Lipopolysaccharide(LPS) and mechanical ventilation(MV)-induced acute lung injury(ALI) mouse models.10,11 In addition, autophagy activation participated in the pathophysiologic process of sepsis, and alleviated the cytokine excessive release and lung injury in sepsis.6 However, ARGs of sepsis-induced ARDS remain largely unknown and need to be further explored, which would provide us with potential biomarkers to treat ARDS.

In this study, we downloaded GSE32707 data from the Gene Expression Omnibus (GEO) database. Next, we identified the key DEARGs that are relevant for diagnosis, which may discriminate sepsis-induced ARDS patients from controls. Subsequently, we used CIBERSORT for the first time to analyze differences in the immune infiltration between sepsis-induced ARDS patients and controls in whole blood samples. Finally, the expression levels of DEARGs were further verified in sepsis-induced ARDS patients and healthy individuals.

Materials and Methods

Autophagy-Related Genes Datasets and Microarray Data

A total of 222 genes were obtained from The Human Autophagy Database (). We extracted information from datasets downloaded from the Gene Expression Omnibus (GEO, ). The GSE32707 gene expression profiles,12 obtained using the GPL10558 platform, contained 144 human blood samples, including 58 from patients with sepsis, 31 from patients with sepsis-induced ARDS, 21 from patients with systemic inflammatory response syndrome and 34 healthy controls. All patients were undergoing mechanical ventilation. The 21 systemic inflammatory response syndrome samples and 58 sepsis samples were excluded from our analysis. The normalize Between Arrays function in the limma package was applied to normalize all the gene expression profiles.13 All the analyses were performed using R software version 3.5.3. (). The workflow of the present study is shown in Figure 1.

Figure 1

Workflow chart of the study.

Differentially Expressed Analysis of Autophagy-Related Genes (DEARGs)

The normalized expression matrix of microarray data was downloaded from the GSE32707 dataset. Then the “limma” package of R software was used to identify the DEARGs. Genes with an adjusted P-value <0.05 and absolute fold-change value >1.5 were considered as differentially expressed genes. The volcano plot, heatmap and box plot were conducted using “ggplot2” and “heatmap” packages of R software.

GO and KEGG Pathway Enrichment Analysis of Autophagy-Related Genes

GO and KEGG pathway enrichment analysis were conducted in R software using the package “clusterProfiler”. The GO analysis consisted of cellular component (CC), biological process (BP) and molecular function (MF).

Correlation Analysis and PPI Analysis of the Differentially Expressed Autophagy-Related Genes

The correlation analysis of the DEARGs was identified using Spearman correlation in the “corrplot” package of R software. PPI analysis of the DEARGs was analyzed using STRING database () and Cytoscape software (version 3.5.1).

Receiver Operating Characteristic (ROC) Curve Analysis

The receiver operating characteristic curves (ROCs) were performed to calculate the area under the curve (AUC) on screened DEARGs which their diagnostic value were assessed using ROC curve analysis with the p ROC software package in R.14

Estimation of Immune Cell Subtype Distribution

The CIBERSORT algorithm was used to transform the normalized gene expression value into the component of peripheral blood immune cells with the LM22 reference signature matrix set at 1000 permutations.15 We retained the CIBERSORT output and calculated the percentage of immune cell subtypes in all samples. We applied CIBERSORT to characterize the immune cells in sepsis-induced ARDS and control samples. In R we used box plots to depict the infiltration levels of immune cells from different subpopulations between the two groups.

ARDS Patients and Healthy Individuals

A total of 6 sepsis-induced ARDS patients (cases) and 6 age and sex-matched healthy individuals (controls) were obtained from the First Affiliated Hospital of Zhengzhou University between October 2021 and November 2021. Sepsis-induced ARDS was defined according to the Berlin consensus definition.16 Sepsis was defined according to standard criteria.17 The 6 healthy individuals were recruited from the hospital’s health checkup center. This study was conducted in accordance with the Declaration of Helsinki and approved by the Medical Ethics Committee of the hospital. Written informed consent was obtained from all the participants. Venous blood was collected from all cases and controls who participated in the study.

RNA Extraction and Quantitative Real-Time Polymerase Chain Reaction (qRT-PCR)

The blood samples from each participant were processed to isolate peripheral blood mononuclear cells (PBMCs) using Ficoll solution (Solarbio Life Sciences, Beijing, China). Total RNA was extracted from PBMCs with RNA Extraction Kit (TRIzon Reagent, Shanghai, China). Reverse transcription was conducted using Prime Script RT Master Mix Kit (Takara, Dalian, China). The mRNA level was assessed using TB Green Premix Ex Taq Kit (Takara, Dalian, China) following the instructions. Relative expression levels were normalized to GAPDH. Primers are available in . The relative expression of mRNA was calculated by 2−ΔΔCt method with the normalization to GAPDH.

Results

Identification of DEARGs Between Normal and Sepsis-Induced ARDS

We analyzed the expression of 222 ARGs in 31 sepsis-induced ARDS patients and 34 healthy individuals, and 28 ARGs were identified using the criteria of adjusted P-value <0.05 and absolute fold-change value >1.5, including 11 up-regulated genes and 17 down-regulated genes (Table 1). Following the analysis of the GSE32707 dataset with R software, the 28 DEARGs between sepsis-induced ARDS and normal groups were presented in volcano plot and heatmap (Figure 2A and B). Moreover, box plots showed the expression patterns of 28 DEARGs between ARDS and normal samples (Figure 3A and B).

Table 1

The 28 DEAGs in ARDS Samples Compared to Healthy Sample

Gene Symbol	P value	Adj.P.val	logFC	Regulation
ATG9A	0.00213	0.021722	0.6034436	Up
RAB24	0.000966	0.013275	0.7408382	Up
SH3GLB1	0.00728	0.049258	0.74861416	Up
ATG16L2	0.0026	0.024716	0.84585503	Up
NLRC4	0.000475	0.008903	0.88533125	Up
MAP1LC3A	0.001	0.013542	0.89767508	Up
CTSD	0.00022	0.005822	0.9949901	Up
CTSB	0.00425	0.034403	1.06796432	Up
TNFSF10	0.00199	0.020789	1.22727772	Up
CFLAR	0.00263	0.024901	1.63170317	Up
SERPINA1	0.00292	0.02674	1.78327746	Up
MYC	0.00000599	0.001596	−1.2602977	Down
BNIP3	0.000235	0.006019	−1.12920532	Down
HSP90AB1	0.000625	0.010274	−1.09186858	Down
PARP1	0.000109	0.004251	−1.01156063	Down
HSPA8	0.00307	0.02759	−0.97181352	Down
PELP1	0.000257	0.006311	−0.86776787	Down
RHEB	0.000428	0.008357	−0.80968619	Down
ATIC	0.000407	0.008145	−0.79693892	Down
BIRC5	0.00523	0.039543	−0.77814748	Down
EEF2	0.00234	0.023072	−0.77500085	Down
ITGB1	0.00235	0.023151	−0.77400807	Down
SPNS1	0.00194	0.02048	−0.74292339	Down
PEA15	0.000441	0.008493	−0.72592288	Down
BAG3	0.0000199	0.002406	−0.67128778	Down
ITGA3	0.00105	0.013959	−0.60020641	Down
SIRT1	0.0000217	0.002444	−0.59783683	Down
ERBB2	0.0000184	0.002391	−0.58585733	Down

Abbreviations: ATG9A, autophagy related 9A; RAB24, RAB24, member RAS oncogene family; SH3GLB1, SH3 domain containing GRB2 like, endophilin B1; ATG16L2, autophagy related 16 like 2; NLRC4, NLR family CARD domain containing 4; MAP1LC3A, microtubule associated protein 1 light chain 3 alpha; CTSD, cathepsin D; CTSB, cathepsin B; TNFSF10, TNF superfamily member 10; CFLAR, CASP8 and FADD like apoptosis regulator; SERPINA1, serpin family A member 1; MYC, MYC proto-oncogene; BNIP3, BCL2 interacting protein 3; HSP90AB1, heat shock protein 90 alpha family class B member 1; PARP1, poly(ADP-ribose) polymerase 1; HSPA8, heat shock protein family A (Hsp70) member 8; PELP1, proline, glutamate and leucine rich protein 1; RHEB, Ras homolog, mTORC1 binding; ATIC, 5-aminoimidazole-4-carboxamide ribonucleotide formyltransferase/IMP cyclohydrolase; BIRC5, baculoviral IAP repeat containing 5; EEF2, eukaryotic translation elongation factor 2; ITGB1, integrin subunit beta 1; SPNS1, sphingolipid transporter 1 (putative); PEA15, proliferation and apoptosis adaptor protein 15; BAG3, BAG cochaperone 3; ITGA3, integrin subunit alpha 3; SIRT1, sirtuin 1; ERBB2, erb-b2 receptor tyrosine kinase 2; ARDS, acute respiratory distress syndrome; ARGs, autophagy-related genes; GEO, Gene Expression Omnibus; DEGs, differentially expressed genes; DEARGs, differentially expressed ARGs; PPI, protein–protein interaction; ROC, receiver operating characteristic; AUC, area under the curve; 95% CI, 95% confidence interval; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; BP, biological process; CC, cellular component; MF, molecular function; qRT-PCR, quantitative real-time polymerase chain reaction; CI, confidence interval; IL-6, interleukin-6; IL-8, interleukin-8; ANG-2, angiopoietin-2; KL-6, Kerbs von Lungren 6 antigen; RAGEs, receptor for advanced glycation endproducts; SP-D, surfactant protein D; NLRP3, NLR family pyrin domain containing 3; IL-1β, interleukin-1β; LPS, lipopolysaccharide; MV, mechanical ventilation; GAPDH, glyceraldehyde-3-phosphate dehydrogenase; ATG-12, autophagy related-12.

Figure 2

DEARGs in sepsis-induced ARDS and healthy samples. (A) Volcano plot of the 222 DEARGs. The red dots represent the significantly up-regulated genes and the blue dots indicate the significantly down-regulated genes. (B) Heatmap of the 28 DEARGs in sepsis-induced ARDS and healthy samples.

Figure 3

The boxplot of 28 DEARGs in sepsis-induced ARDS and healthy samples. (A) The boxplot of top 14 DEARGs in sepsis-induced ARDS and healthy samples. (B) The boxplot of last 14 DEARGs in sepsis-induced ARDS and healthy samples. P-values were calculated using a two-sided unpaired Student’s t-test. *P<0.05; **P<0.01.

GO and KEGG Enrichment Analysis of the DEARGs

We performed GO and KEGG enrichment analyses using R software to determine the potential function of DEARGs in the development of sepsis-induced ARDS. The BP analysis showed that these genes are associated with the autophagy, process utilizing autophagic, as well as with the macroautophagy in the sepsis-induced ARDS (Figure 4A). In terms of the CCs, these genes were involved in the autophagosome, vacuolar membrane, and ficolin-1-rich granule (Figure 4A). With regard to MF, these genes participated in certain key functions, such as cell adhesion molecule binding, cadherin binding, and protease binding (Figure 4A) (). In addition, the KEGG pathway enrichment analysis indicated that these genes were associated with Autophagy-animal, Apoptosis, PI3K-Akt signaling pathway. (Figure 4B, ).

Figure 4

Functional enrichment of the DEARGs. (A) BP, biological process; CC, cellular component; MF, molecular function. (B) KEGG pathway enrichment analysis.

Correlation Analysis and PPI Network of the DEARGs

To explore the expression correlation of these ARGs, correlation analysis was performed. The results revealed the relationship of the 28 DEARGs in the GSE32707 dataset (Figure 5). To determine the interactions between the DEARGs, PPI analysis was performed. The results demonstrated that these ARGs interacted with each other (Figure 6).

Figure 5

Spearman correlation analysis of the 28 DEARGs.

Figure 6

Protein–protein interactions (PPI) analysis the 28 DEARGs. Red and blue circles indicate upregulated and downregulated differentially expressed genes, respectively. Circle size indicates the node degree.

We performed ROC curve analysis using R Studio to further validate the reliability of these ARGs in the peripheral blood of patients. Seven genes (BAG3, CTSD, ERBB2, MYC, PEA15, RAB24, and SIRT1) with an AUC greater than 0.70 were considered possible to be hub genes, indicating that they may be able to diagnose patients with sepsis-related ARDS with excellent specificity and sensitivity (Figure 7).

Figure 7

ROC curves to evaluate the diagnostic value of the DEARGs.

Distribution Pattern of Immune Cell Subtype

The overall distribution of different immune subsets in all samples was shown in a histogram (Figure 8A). Different colors and heights represent various types and percentages of immune cells in the sample, and the sum of the proportion of various immune cells is 1. CD8+ T cells, naive CD4+ T cells, memory resting CD4+ T cells, follicular helper T cells, gamma delta T cells, monocytes, mast cells resting, eosinophils and neutrophils were the primary infiltrating cells. There are individual differences in the proportion of immune cells between the two groups. The different immune cell subsets were weakly-to-moderately correlated (Figure 8B). Compared with healthy controls, the proportion of naive CD4+ T cells, gamma delta T cells, monocytes, and neutrophils was higher in sepsis-related ARDS samples, while the fraction of CD8+ T cells, memory resting CD4+ T cells, follicular helper T cells, mast cells resting, and eosinophils was relatively lower (Figure 8C).

Figure 8

The profiles of immune cell subtype distribution pattern in GSE32707 cohort. (A) The bar plot visualizing the relative percent of 22 immune cell in each sample. (B) Correlation heatmap of all 22 immune cells. (C) Boxplot of all 22 immune cells differentially infiltrated fraction. P-values were calculated using a two-sided unpaired Student’s t-test. *P<0.05; **P<0.01.

Validation the DEARGs by Real-Time PCR Analysis

To validate the reliability of the seven DEARG expression levels that were further identified by qRT-PCR in our clinical samples. The qRT-PCR results indicated that CTSD was upregulated and other four genes (BAG3, ERBB2, MYC and SIRT1) were downregulated. Only the expression of PEA15 and ERBB2 were inconsistent with our integrated analysis (Figure 9).

Figure 9

RNA expression of seven DEARGs were measured in sepsis-induced ARDS and healthy samples. (A) RNA expression of BAG3, CTSD, ERBB2, MYC, PEA15, RAB24 and SIRT1were measured in blood samples using qRT-PCR. P-values were calculated using a two-sided unpaired Student’s s t-test. **P<0.01.

Discussion

Sepsis-induced ARDS has poorer recovery from lung injury, higher overall disease severity and higher mortality than non-sepsis-induced ARDS.18 Unfortunately, mechanisms of sepsis-induced ARDS still remain unclear and effective therapies are lacking. Increasing evidences suggested autophagy plays critical roles in the pathogenesis of sepsis-induced ARDS. Autophagy is a double-edged sword. In the early stages, sepsis-induced ARDS can activate various autophagy-related proteins, which reduce cytokine release, inflammatory response, and apoptosis in the lung.6,19,20 Another evidence suggested that autophagy, by virtue of excessive autophagosomes accumulation, may play a maladaptive role in the late stages sepsis-related ARDS.21,22 Therefore, there is an urgent need to improve the understanding of autophagy in pathogenesis of sepsis-induced ARDS by extensive validations.

Recently, a series of studies have explored that autophagy is associated with sepsis-induced ARDS. Colell et al results indicated that glyceraldehyde-3-phosphate dehydrogenase (GAPDH) protected against sepsis-related lung injury by promoting autophagy related 12 (ATG12)-mediated autophagy and ameliorated oxidative stress.23,24 In addition, a recent study found that MicroRNA-377-3p released from mesenchymal stem cell exosomes improves LPS-induced ALI by targeting RPTOR to induce autophagy.25 Another study reported that CD39+ regulatory T cells attenuate LPS-induced ALI via autophagy and ERK/FOS pathway.26 However, bioinformatics analysis of ARGs has not been well explored in sepsis-induced ARDS.

In this study, we identified 28 potential DEARGs of sepsis-induced ARDS by bioinformatics analysis. The potential biological functions of these DEARGs were performed by GO and KEGG enrichment analysis which mainly indicated a key role of autophagy in the pathophysiological mechanisms of sepsis-induced ARDS.

We further performed ROC curves to assess the diagnostic ability of DEARGs in sepsis-induced ARDS. The seven genes (BAG3, CTSD, ERBB2, MYC, PEA15, RAB24 and SIRT1) were indicated that they may be able to diagnose sepsis-induced ARDS patients with excellent specificity and sensitivity. Abundant research has confirmed that the immune system plays a central role in sepsis-induced ARDS. Therefore, we used CIBERSORT, a transcriptome deconvolution algorithm, to calculate the proportion of blood immune cell subsets which can better understand the role of leukocyte infiltration and the inflammatory response in the pathogenesis and development of sepsis-related ARDS. Besides, there were also obviously differences in the proportion of most of 22 immune cells between sepsis-related ARDS and healthy controls. These differences include an increase in the proportion of naive CD4+ T cells, gamma delta T cells, monocytes, and neutrophils, while the fraction of CD8+ T cells, memory resting CD4+ T cells, follicular helper T cells, mast cells resting, and eosinophils was relatively lower in sepsis-induced ARDS.

Finally, based on bioinformatics analysis results, the expression levels of seven DEARGs were further identified by qRT-PCR in our clinical samples. CTSD was upregulated and the other four genes (BAG3, ERBB2, MYC and SIRT1) were downregulated. Only the expression of PEA15 and ERBB2 were inconsistent with our integrated analysis. Some of these verified genes have been reported to be associated with sepsis-related ARDS. Liu et al revealed that SIRT1 can enhance cellular autophagy via regulating of ATG7 and alleviate LPS-induced ARDS.27 Dong et al found that LC3 and ATG16 expressions are elevated in activated autophagy, which maintain barrier and adherens functional integrity in the murine model of ALI. Besides, activated autophagy can alleviate cytokine release and further improve ALI evoked by LPS stimulation in mice.28 Nevertheless, the explicit mechanism whereby these verified genes remain largely unclear in sepsis-induced ARDS to need further exploration.

Our research inevitably presents some limitations. Firstly, our research represents the second mining and analysis of previously published data sets, so we did not know the severity of patients, treatment methods and relevant clinical data, which may have confounded factors in our analysis. Secondly, the small number of clinical samples were used for analysis, which further studies with a larger sepsis-induced ARDS cohort will be needed to further confirm and extend the conclusions of this study. Thirdly, we only verified the expression level of the differentially expressed autophagy-related genes in clinical samples, owing to the lack of sepsis-induced ARDS detailed information, it is difficult to draw a clear association between selected genes and the severity of sepsis-induced ARDS while using the same samples. Therefore, mechanisms of sepsis-induced ARDS are needed to further explore in vitro or in vivo experiments.

Conclusion

We identified an association between DEG and immune infiltration in sepsis-induced ARDS and validated BAG3, CTSD, ERBB2, MYC and SIRT1 that may be have excellent diagnostic performance.