Time-course transcriptome analysis of lungs from mice infected with inhaled aerosolized Stenotrophomonas maltophilia

Introduction

Stenotrophomonas maltophilia (SMA) is an important gram-negative multidrug-resistant bacterial pathogen that is widely distributed in the environment and has an effect on human health around the world (1). SMA is an opportunistic pathogen that causes a variety of serious infections in humans, including respiratory infections, bacteremia, endocarditis, gastrointestinal tract infections, liver infection, urinary tract infections, nervous system and spinal cord infections, eye infections, soft tissue and bone infections, and medical implant infections. Pneumonia and bacteremia are the most common clinical manifestations of this infection, with reported overall mortality rates ranging from 30% to 51% in patients with bacteremia (2), nearly 50% in patients with hospital-acquired pneumonia in intensive care units (3), and 100% in patients with hemorrhagic pneumonia (4). Despite the high mortality rate associated with SMA pneumonia, there is limited experience with the clinical management of pneumonia caused by this bacterium (5-7). In our previous research, outer membrane proteins (Omps) from SMA were isolated and potential vaccine candidates of Omps against SMA were identified (8). As a result, a better understanding of SMA pneumonia is needed.

The pathogenesis of SMA infection involves numerous virulence factors. Outer membrane vesicles (OMVs) secreted by SMA elicit early inflammatory responses in mouse lungs, including congestion and neutrophil infiltration, as well as elevated expression of proinflammatory cytokine and chemokine genes, and stimulate host defense responses (1,9). StmPr1, the major extracellular protease, degrades fibronectin, fibronectin receptor integrin a5b1, collagen, E-cadherin, tight junction proteins ZO-1, and occludin; disrupts the epithelial cell barrier; recruits neutrophils; inactivates host immune defense components and is a major player in lung cell injury (10). Flagella are directly involved in the adhesion to mouse tracheal mucus, and binding to mucus inhibits the movement of bacteria within the respiratory tract and protects them from the host’s immune response (11). Moreover, flagellin has been shown to increase the production of proinflammatory cytokines [interleukin (IL)-1β and tumor necrosis factor alpha (TNF-α)], increase the production of anti-inflammatory cytokines (IL-10), elevate the levels of myeloperoxidase and nitric oxide, significantly increase the engulfment activity of alveolar macrophages and the number of infiltrating leukocytes, and provide nonspecific protection against SMA challenge (12). Overall, these findings suggest that SMA infection may activate different gene expression patterns in various types of host cells, with varying effects on specific cell types. The pathogenesis of SMA infection and SMA-host interactions are complicated.

Several methods of inoculation have been reported for establishing SMA pneumonia models, including nasal inoculation, intratracheal instillation, and aerosol challenge in mice (13-15). In general, nasal inoculation may lead to unstable model establishment due to the inability to accurately control the dose of infection and the highly variable deposition of microorganisms in the lungs. Intratracheal instillation does not readily simulate naturally occurring aerosol-borne infections. Inhalational aerosol pneumonia models mostly involve the route of oral-nasal contact and systemic exposure, but accurate quantification of the administered dose may be limited. Unlike previous methods of inoculation, aerosol endotracheal inhalation uses a hand-held liquid aerosol pulmonary delivery device to directly place a known quantity of the organism into the lungs of mice (16). It is characterized by straightforward operation, high success rates, easy reproducibility and low animal mortality.

The lung is an important organ for SMA infection as well as a major immune organ in vertebrate animals. To investigate the mechanisms of SMA-host interactions, a detailed description of SMA pneumonia is needed. Due to the decreasing cost of next-generation sequencing (NGS), RNA sequencing (RNA-seq) has become a popular method for analyzing transcriptomes in recent years. Numerous studies have used RNA-seq to investigate host-pathogen interactions, host response mechanisms to genetic and infectious diseases (17). Similarly, we performed time-course RNA-seq of lung tissues from mice with SMA pulmonary infection via aerosolized intratracheal inhalation to investigate transcription profile changes in SMA-infected lungs and to identify genes or other valuable research targets. This research advances our understanding of the pathogenesis of SMA pneumonia. We present this article in accordance with the ARRIVE reporting checklist (available at https://jtd.amegroups.com/article/view/10.21037/jtd-23-1138/rc).

Methods

Bacterial strain and mice

The SMA strain K279a, which has been previously described (18), was stored in our lab. Bacteria were cultured overnight in lysogeny broth (LB) in a 37 ℃ shaker. On the second day, it was transferred to 50 mL LB medium at 1:100 for expansion to an optical density (OD) at 600 nm of 1.0 and centrifuged at 3,000 reactive centrifugal force (rcf) for 5 min. The bacterial pellets were washed twice with sterile phosphate-buffered saline (PBS) and adjusted to reach an OD at 600 nm of approximately 1.0 [approximately 5×10⁸ colony forming units (CFU)/mL)]. The bacterial pellets obtained from 1.0 mL volumes of bacterial culture media were resuspended with 50 µL of sterile PBS, and were used to challenge the mice. Male C57BL/6J mice were purchased from SPF (Beijing, China) Biotechnology Co., Ltd. (Beijing, China) and maintained under specific pathogen-free conditions at the Academy of Military Medical Sciences (AMMS) animal facilities for 1 week before use. Six- to seven-week-old mice weighing 18–20 grams were used in the experiment. A protocol was prepared before the study without registration. Animal experiments were performed under a project license (No. IACUC-DWZX-2022-064) granted by Animal Care and Use Committee of the AMMS, in compliance with institutional guidelines for the care and use of animals.

Experimental mouse model of acute lethal SMA pneumonia via aerosolized intratracheal inoculation

Mice were inoculated via the aerosolized intratracheal route with either PBS or SMA as described previously (16). In brief, isoflurane was administered through a small animal anesthesia machine (Yuyan Scientific Instrument Co., Ltd., Shanghai, China) to anesthetize the mice. Supine mice were secured after anesthesia. When the tracheal opening was fully exposed by a small animal laryngoscope (Huironghe Company, Beijing, China), 50 µL of the suspension was sprayed into the trachea (near the bifurcation of the trachea) by using an Intratracheal Aerosolizer Device (Huironghe Company, Beijing, China), followed by uniform dispersion of the solution throughout the lung. To establish an acute lethal pneumonia model, mice (n=10) were noninvasively intratracheally inoculated with 5×10⁸ CFU/mL of SMA in 50 µL of PBS. Post-infection mouse survival rates were monitored for 14 days. For transcriptome analysis and histopathological analysis, mice were humanely sacrificed at 4, 12, 24, 48, and 72 hours post-infection (hpi). So, the mice were divided into 6 groups and 4 mice each group. The mouse lungs were removed and split into two pieces for total RNA extraction or histopathological analysis. An inoculation of 50 µL of PBS was used as a control at 0 hpi.

Histopathology analysis

Paraformaldehyde-fixed and paraffin-embedded lung sections were used. In accordance with standard protocols, 3–5 µm sections were cut and stained with hematoxylin & eosin (H&E), and then the sealing slides were observed by light microscopy (200×). By using H&E staining and a lung injury scoring system, the histology of the lungs was evaluated. Based on histological evaluations, injury scores were 0= none, 1= slight, 2= mild, 3= moderate, and 4= severe. The categorization of the lesions was related to the degree of pathological lesions as follows: inflammatory cell infiltration, interstitial inflammation, congestion, and edema. One-way analysis of variance (ANOVA) was used to analyze pathology scores.

RNA quantification, library construction and sequencing

The mouse lungs, each in 5 mL of TRIzol (Invitrogen, ThermoFisher Scientific Inc., Carlsbad, CA, USA), were stored at −80 ℃ prior to processing and shipped on dry ice. The quantity and integrity of total RNA was analyzed using a Bioanalyzer 2100 system (Agilent Technologies, CA, USA). mRNA was enriched from 3 µg of total RNA using poly-T oligo-attached magnetic beads and fragmented by using divalent cations under elevated temperature. Moloney murine leukemia virus (M-MuLV) reverse transcriptase, random hexamer primers, DNA polymerase I and deoxynucleotide triphosphate (dNTP) were used to create double-stranded complementary DNA (cDNA), which were then further purified with AMPure XP beads. The final cDNA libraries were created by further purifying and enriching the cDNA fragments using polymerase chain reaction (PCR) (15 cycles) and sequenced by using an Illumina NovaSeq 6000. Construction of the library and sequencing were carried out by Novogene Company (Beijing, China).

Analysis of RNA-seq data

The high-throughput sequencer’s measurements of image data were converted into sequence data (reads) by CASAVA (Illumina, San Diego, CA, USA) base recognition. By eliminating adapter and N-base-containing reads as well as low-quality reads from the raw data, clean data (clean reads) were obtained. The Q20, Q30 and GC contents of the clean data were calculated. Direct downloads of the reference genome and gene model annotation files from the genome (mm10) website were performed. Hisat2 (v2.0.5) was used to create an index of the reference genome and to match clean paired-end reads to the reference genome. Read counts that were mapped to each gene were counted using featureCounts (v1.5.0-p3, SourceForge, San Diego, CA, USA). The length of each gene and the number of reads mapped to it were used to compute the fragments per kilobase of transcript per million mapped reads (FPKM) of each gene. Principal component analysis (PCA) was used to analyze the distribution of samples to judge the quality of the data.

Identification of differentially expressed genes (DEGs) and functional classification of DEGs

The DESeq2 R package (1.20.0) (http://www.bioconductor.org/) was used to analyze differential expression between the two groups. P values were adjusted using Benjamini and Hochberg’s method to reduce the false discovery rate. Using an adjusted P value (Padj) £0.05 and |log₂fold change| ³1 as cutoffs in this study, DEGs were identified and visualized using a volcano plot. DEGs were also functionally classified. Gene Ontology (GO) functional enrichment was conducted using the GO database (https://geneontology.org/), and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment was conducted using the KEGG database (http://www.genome.jp/kegg/).

Time courses analysis of gene expression

To cluster the genes identified in the time courses, we used the R package Mfuzz (http://www.bioconductor.org/) (19), which performs fuzzy c-means clustering. To remove the influence of expression magnitude and focus on the expression dynamics, DESeq2-normalized data were standardized in Mfuzz to have a mean value of zero and a standard deviation of one. A total of 9 clusters were used, and the fuzzifier coefficient, M, was set to 1.71.

Coexpression network establishment by weighted gene coexpression network analysis (WGCNA)

WGCNA is a computational technique for describing the mechanisms of gene relationships among various samples. Based on the coherence of gene sets and the association between gene sets and symptoms, it can be utilized to identify gene sets that are highly synergistically altered and suggest possible biomarkers or therapeutic targets (20). Using the WGCNA package in R, a weighted correlation network analysis was performed for the 22 lung samples to identify modules for all genes with various expression patterns. First, the data were filtered by calculating the coefficient of variation to remove genes with low coefficients of variation, which was 0.5 for this metric. The gene network was assumed to obey a scale-free network distribution, the gene coexpression correlation matrix was defined, and the coefficient of dissimilarity of different nodes was calculated. Using hierarchical clustering, genes with highly similar coexpression relationships were grouped according to their topological overlaps. Modules with the highest coefficients were selected as targets.

Immune cell infiltration score

The Immune Cell Abundance Identifier for mouse (ImmuCellAI-mouse, http://bioinfo.life.hust.edu.cn/ImmuCellAI-mouse/#!/) is a tool that estimates the relative proportions of 36 immune cell species based on mouse gene expression data (21). Using a hierarchical prediction strategy, ImmuCellAI-mouse was used to predict the abundance of immune cells and immune cell subtypes separately by simulating the flow cytometry analysis process. Layer 1 comprises seven different types of immune cells: T cells, B cells, natural killer (NK) cells, monocytes, granulocytes, macrophages, and dendritic cells (DCs). Layer 2 cells were subtypes of cells in the first layer, such as T-cell subtypes [CD4 T, NKT, CD8 T, and γδ T cells], B-cell subtypes (B1, memory B, marginal zone B, germinal center B, follicular B, and plasma B cells), subtypes of granulocytes (neutrophils, basophils, mast cells and eosinophils), subtypes of macrophages (M1 and M2 macrophages), and subtypes of DCs [monocyte-derived DCs (moDCs), plasmacytoid DCs (pDCs), type 1 conventional DCs (cDC1s) and type 2 conventional DCs (cDC2s)]. Finally, cells in layer 3 are CD4 T and CD8 T cell subtypes (including CD4 naïve, CD4 memory, T helper, regulatory T (Treg), CD8 naïve, CD8 effector memory, CD8 central memory, cytotoxic and exhausted cells). The prediction process used to analyze each sequential layer was performed separately.

Validation of RNA-seq data by quantitative real-time PCR (qRT-PCR)

To validate the RNA-seq data, 9 DEGs were randomly chosen for qRT-PCR validation of their expression patterns at five time points. Lung tissue from mice infected with SMA at 4, 12, 24, 48, and 72 hpi was dissected and snap-frozen in liquid nitrogen. RNA was extracted and measured using spectrophotometry. The TRIzol reagent was obtained from Sangon Biotech Co., Ltd. (Shanghai, China), as well as the AMV First Strand cDNA Synthesis Kit. PowerUp^TM SYBR^TM Green master mix (Thermo Scientific, USA) was used for qRT-PCR. The housekeeping reference gene Gapdh was used to compute relative quantification. Primers are listed in Table 1.

Statistical analysis

Comparisons between two groups were carried out using two-tailed unpaired t-tests. Nonparametric analysis was performed on data that lacked a normal distribution using the Mann-Whitney U test. One-way ANOVA was used to analyze pathology scores. Using the Kaplan-Meier method and log rank test, survival experiments were analyzed. Results with P<0.05 were defined as statistically significant. Statistics were analyzed using SPSS 16 (IBM, Arnmonk, NY, USA) and GraphPad Prism 8.0 (GraphPad, La Jolla, CA, USA).

Results

An acute lethal pneumonia model and histopathological changes in the lungs via aerosolized intratracheal inoculation of SMA

Mice infected with SMA via aerosolized intratracheal inoculation were monitored daily for survival events for 14 days. A control inoculation of 50 µL of PBS was used. Two mice whose death occurred within the initial 3-hour period were deemed as unexpected deaths and hence eliminated from the research analysis. Very few animals (3 mice) died within the first 4 hpi, probably due to the combined effects of asphyxia and anesthesia. These animals were excluded from the experimental group. Most mice in the model group died within 4 days, and mortality in this group differed from that in the control group (P<0.01) (Figure 1A). As a result, we selected six time points for subsequent experiments and analysis: 0, 4, 12, 24, 48, and 72 hpi (Figure 1B). Histopathological changes assessed using H&E staining revealed neutrophil infiltration in lung tissue at 4 hpi, followed by a significant increase in the infiltration of neutrophils and monocytes at 12–72 hpi. At 72 hpi, there was massive destruction of alveolar structures, thickening of alveolar septa, and hemorrhage. The H&E staining results were similar at 48 hpi (not presented) and 72 hpi. Pathological staining of the lungs of mice that died from the infection showed severe pulmonary hemorrhage and a large level of neutrophil infiltration. This result is consistent with findings reported previously (22). The severity of the lung lesions increased over time, according to histological scoring (Figure 1C). During histopathological verification, mice infected with SMA were shown to exhibit excessive acute inflammation and severe lung damage. As a result, the mouse model of acute lethal SMA pneumonia was successfully established, and all subsequent experiments were conducted using this animal model.

Figure 1 Survival curves and histopathological analysis of lung tissue of mice infected with SMA via aerosol intratracheal inoculation. (A) Survival curves of mice after pulmonary infection. Mice were challenged with either PBS or 5×10⁸ CFU/mL of the SMA strain (K279a) (n=10 per group). Survival differences between infected and control mice were statistically significant (**, P<0.01; log-rank test). In a repeat experiment, similar results were obtained. (B) Histopathological analysis and H&E staining of the lungs (200×). Histopathological changes assessed using H&E staining revealed neutrophil infiltration in lung tissue at 4 hpi, followed by a significant increase in the infiltration of neutrophils and monocytes at 12–72 hpi. At 72 hpi, there was massive destruction of alveolar structures, thickening of alveolar septa, and hemorrhage. (C) Pathological scores of lung sections in mice were determined at 0, 4, 12, 24, 48, and 72 hpi. *, P<0.05; ***, P<0.001; ****, P<0.0001 compared with 0 hpi. SMA, Stenotrophomonas maltophilia; hpi, hours post-infection; PBS, phosphate-buffered saline; CFU, colony forming units; H&E, hematoxylin & eosin.

Overview of mRNA sequencing data and transcriptome analysis

To identify the expression levels of mRNAs in the lungs of mice infected by SMA at different times, 22 cDNA libraries were constructed and sequenced (Table 2). All 22 libraries were shown to be appropriate for future research. PCA is commonly used to assess between-group differences and sample replicates within groups. Our PCA of the gene expression value (FPKM) for all samples revealed that samples were dispersed across groups and clustered together within groups. Principal component 1 (PC1) was 37.79%, and principal component 2 (PC2) was 26.62% (Figure 2A). DEGs were identified across five time points after SMA infection using the control group (0 hpi) as a reference. A total of 3,165 upregulated DEGs and 2,231 downregulated DEGs were detected at 4 hpi, 2,798 upregulated DEGs and 2,774 downregulated DEGs were detected at 12 hpi, 2,423 upregulated DEGs and 2,398 downregulated DEGs were detected at 24 hpi, 2,601 upregulated DEGs and 2,877 downregulated DEGs were detected at 48 hpi, and 2,512 upregulated DEGs and 2,468 downregulated DEGs were detected at 72 hpi (Figure 2B,2C). A total of 1,078 common DEGs upregulated and 635 downregulated were identified for all time points (Figure 2D,2E). We focused on analyzing gene expression profiles in the lungs within the first 4 hours to increase understanding of the rapid and intense inflammatory response of the host. The top 20 genes with the most significant levels of increased expression at 4 hpi were selected. The Padj and log₂fold change at each time point post-infection of these top 20 DEGs are shown in Table 3.

Figure 2 Overview of time-course transcriptome analysis based on the FPKM levels of all groups. (A) Principal component analysis was performed on all samples based on gene expression (FPKM). (B) Volcano map of DEGs at different time points compared with the control group. (C) Different time points of DEG upregulation and downregulation. (D) Venn diagram showing up-regulated DEGs identified in different groups. (E) Venn diagram showing down-regulated DEGs identified in different groups. Each transcriptome analysis included 3–4 biological replicates. PC, principal component; PBS, phosphate-buffered saline; hpi, hours post-infection; SMA, Stenotrophomonas maltophilia; Padj, adjusted P value; FPKM, fragments per kilobase of transcript per million mapped reads; DEGs, differentially expressed genes.

Time courses gene clustering

To investigate the dynamics of gene expression, we examined time-associated gene clusters of all genes from different time points using Mfuzz 2.58.0 (Figure 3A). In each cluster, the GO enrichment analysis is shown (Figure 3B). We found that the identified processes corresponded to the molecular pathophysiology of disease progression.

Figure 3 Clustering analysis of time series genes during SMA pneumonia. (A) Mfuzz clustering revealed nine distinct temporal patterns of gene expression. Red indicates that the gene variation was the most constrained to the center of the cluster, followed by blue and finally green. The Mfuzz membership value, which is used to determine whether the gene conforms to the trend of the gene cluster, ranges between 0.128 and 0.950. (B) The heatmap depicts the importance of GO terms in the biological processes that describe each of the nine clusters. The number (4,694, 3,189, 4,081, 3,349, 2,770, 3,008, 1,343, 2,820 and 5,397) of enriched genes in each cluster is shown at the bottom. hpi, hours post-infection; rRNA, ribosomal RNA; ERAD, endoplasmic reticulum-associated protein degradation; Padj, adjusted P value; ncRNA, non-coding RNA; SMA, Stenotrophomonas maltophilia; GO, Gene Ontology.

Gene expression in Cluster 1 and Cluster 2 rapidly peaked at 4 hpi and then declined. GO analysis showed that 4,694 genes in Cluster 1 were enriched for biological processes related to circulatory and muscular functions, including heart process, muscle system process, heart contraction, muscle contraction and cardiac muscle contraction. A total of 3,189 genes in Cluster 2 were enriched for biological processes related to the inflammatory signaling pathway and cytokine production, such as positive regulation of cytokine production, I-κB kinase (IKK)/NF-κB signaling, regulation of innate immune response and cytokine-mediated signaling pathway.

Gene expression levels in Clusters 5 and 7 began to rise immediately after infection. According to GO analysis, the 2,770 and 1,343 genes in Clusters 5 and 7, respectively, were primarily enriched for biological processes related to DNA replication and repair and protein repair, such as DNA repair, DNA replication, DNA-dependent DNA replication, chromosome segregation, double-strand break repair, ubiquitin-dependent endoplasmic reticulum-associated protein degradation (ERAD) pathway, ERAD pathway, response to topologically incorrect protein, cellular response to topologically incorrect protein, and response to unfolded protein.

After infection, the gene expression levels in Clusters 4 and 9 decreased. The 3,349 and 5,397 genes in Clusters 4 and 9, respectively, were mainly enriched for biological processes related to cell adhesion, cell junction, and smooth muscle cell contraction, such as cell-cell adhesion via plasma-membrane adhesion molecules, cell junction organization, smoothened signaling pathway, and smoothened signaling pathway regulation.

Cluster 6 gene expression peaked at 12 hours post-incubation. These 3,008 genes were enriched for biological processes related to transcriptional activity, such as transcriptional activator activity, RNA polymerase II transcription regulatory region sequence-specific DNA binding, transcription factor activity, and RNA polymerase II proximal promoter sequence-specific DNA binding, according to GO analysis.

Cluster 8 gene expression levels peaked at 24 hpi. GO analysis revealed that these 2,820 genes were enriched in biological processes related to noncoding RNA metabolic processes, response to endoplasmic reticulum stress (ERS), autophagy, processes utilizing autophagic mechanisms and the ERAD pathway.

WGCNA and enrichment analysis for WGCNA module genes

WGCNA hierarchical clustering was used to create dendrograms of 20 distinct modules of highly coexpressed genes; each module is colored differently, and genes with no cluster were found in the gray module (Figure 4A). The correlation between the infection process (hpi) and the module genes was then calculated (Figure 4B). Enrichment analysis was performed on modules with correlation coefficients greater than 0.9 to determine their importance in our study.

Figure 4 Gene modules in WGCNA. (A) WGCNA-identified coexpression modules are shown in a hierarchical clustering tree. Different modules are identified by different colors. A gene is represented by each leaf of the cluster tree. (B) A heatmap showing module-trait relationships. Red indicates a positive correlation with the stage of infection, and blue indicates a negative correlation. ME, module eigengene; PBS, phosphate-buffered saline; hpi, hours post-infection; SMA, Stenotrophomonas maltophilia; WGCNA, weighted gene coexpression network analysis.

According to our module-trait association analysis, the blue module had a positive relationship with 4 hpi. According to GO analysis, this module was primarily enriched in muscle system processes, such as striated muscle contraction, myofibril assembly, muscle contraction, and cardiac muscle contraction. The main enriched pathways, according to KEGG analysis, were related to the IL-17 signaling pathway, the citrate cycle (TCA cycle), cardiac muscle contraction, tight junction, complement, and coagulation cascades, and the TNF signaling pathway (Figure 5A).

Figure 5 Fold enrichment of top-level overrepresented GO and KEGG terms within modules that are highly positively correlated with different infection stages. (A) Blue module; (B) light green module. GO, Gene Ontology; Padj, adjusted P value; KEGG, Kyoto Encyclopedia of Genes and Genomes; TNF, tumor necrosis factor; IL-17, interleukin 17.

At 24 hpi, the light green module was critical. GO analysis revealed that this module was primarily enriched in ERS, including ERS response, endoplasmic reticulum unfolded protein response, cellular response to topologically incorrect protein, and response to unfolded protein. The main enriched pathways, according to KEGG analysis, were related to folate biosynthesis, amino sugar and nucleotide sugar metabolism, cell adhesion molecules, and endocytosis (Figure 5B).

Analysis of immune cell abundance

Using transcriptomic data, ImmuCellAI-mouse was used to assess the abundance of 36 immune cell types in lung tissue at different infection stages. ImmuCellAI-mouse was used to perform a hierarchical strategy to divide immune cell types into three layers. B cell, NK cell, and T cell are three of the lymphoid lineage cells in layer 1, and macrophage, DC, monocyte, and granulocyte are four of the myeloid lineage cells. Cells in layer 3 are subtypes of CD4 T and CD8 T cells, while cells in layer 2 are subtypes of cells in layer 1 (23). Figure 6A depicts abundant immune cell populations of various types in each sample. As shown in Figure 6B, regarding the levels of seven major immune cell types in layer 1, granulocytes, B cells, NK cells, and DCs showed an overall decreasing trend with time, while the levels of macrophages and monocytes showed an overall increasing trend. In layer 2, the abundance of 10 immune cell subtypes was examined further, and a significant increase in the levels of innate immune cells, including neutrophils, monocytes and macrophages, was found (Figure 6C). This suggests that innate immunity takes precedence after infection.

Figure 6 Analysis of the abundance of immune cells in different immune subtypes by ImmuCell-mouse. (A) The presence of 36 different types of immune cells in all samples. Abundance of immune cells in layer 1 (B) and layer 2 (C) at different infection stages. For each cell type, the x-axis represents the median cell abundance of all samples at each time point. The abundance represented the abundance of immune cells, which was predicted using Immune Cell Abundance Identifier (ImmuCellAI-mouse) based on transcriptomic data. *, P<0.05; **, P<0.01; ***, P<0.001; ****, P<0.0001, compared with PBS_0hpi. PBS, phosphate-buffered saline; hpi, hours post-infection; SMA, Stenotrophomonas maltophilia; NK, natural killer; NKT, natural killer T; γδ T, gamma delta T cell; cDC1, type 1 conventional DC; DC, dendritic cell; cDC2, type 2 conventional DC; MoDC, monocyte-derived DC; pDC, plasmacytoid DC; Tm, memory T cell; Treg, regulatory T cell; Tc, cytotoxic T cell; Tcm, central memory T cell; Tex, exhausted T cell; ns, no significance; ImmuCellAI-mouse, Immune Cell Abundance Identifier for mouse.

Validation of the transcriptome data

To validate the RNA-seq data, we randomly selected 9 DEGs, namely, Csf3, Fos, Fosl1, Ifng, Mmp9, Myo1f, Ncf1, Ptafr, and Tnfaip3, for qRT-PCR analysis at five different time points. The Gapdh gene was used as an internal standard. According to the qRT-PCR results, the 9 DEGs showed upregulated or downregulated expression patterns that were similar to the RNA-seq results (Figure 7). RNA-seq and qRT-PCR detected similar variations in the levels of DEGs.

Figure 7 qRT-PCR validation of RNA-seq data. Nine DEGs were chosen for validation, and the data showed a similar trend with RNA-seq data. The y-axis depicts the fold changes of various infection stages relative to the starting point, with positive values indicating upregulation and negative values indicating downregulation. Three biological replicates were used to obtain each data point. RNA-seq, RNA sequencing; FC, fold change; qRT-PCR, quantitative real-time polymerase chain reaction; hpi, hours post-infection; DEGs, differentially expressed genes.

Discussion

In this study, C57BL/6 mice were challenged with 5×10⁸ CFU/mL of SMA by aerosolized intratracheal inhalation, lung sections from challenged mice were subjected to H&E staining, and lung time-course transcriptome expression was analyzed by RNA-seq. A mouse model of acute lethal SMA pneumonia was successfully validated by aerosolized intratracheal inoculation. Using this model, we investigated, for the first time, changes in the time-course transcriptome profile of infected mouse lungs at different time points to better understand SMA lung infection and to identify genes or targets with research potential.

To avoid the drawback of using a single dataset per infection time, we used the temporal clustering of samples to group infection times that shared similar expression patterns. Multiple bioinformatics analysis methods, including Mfuzz time clustering, WGCNA, and ImmuCellAI-mouse, were used in this study to gain insight into the gene expression pattern during SMA pulmonary infection.

Transcriptome analysis revealed a rapid increase in the expression of genes encoding proinflammatory mediators (chemokines, cytokines, and transcription factors) early in the infection (4 hpi), indicating the initiation of an acute inflammatory response. In the Mfuzz analysis, gene expression levels for Cluster 2 peaked at 4 hpi, and GO analysis revealed increased gene enrichment for biological processes related to inflammatory signaling pathways, regulation of innate immune response, and cytokine production. The blue module, which was found to be strongly associated with the 4 hpi timepoint in WGCNA, was linked to inflammatory signaling pathways, such as the IL-17 signaling pathway, complement and coagulation cascades, and the TNF signaling pathway. Similarly, the abundance of innate immune cells determined during ImmuCellAI-mouse analysis increased at 4 hpi. IKK is primarily involved in inflammatory and innate immune responses, mediating the rapid response of recruited immune cells to inflammation and injury via NF-κB signaling (24). The NF-κB signaling pathways play a crucial role in inflammatory cytokine production as part of the innate immune system (25). GO analysis revealed that the blue module strongly associated with the 4 hpi timepoint was primarily enriched in muscle system processes; genes in Clusters 1 and 3 were mainly enriched in biological processes related to circulation and muscle function, as well as mitochondrial respiratory chain function. Energy is needed to fight infection and for both anti-inflammatory and proinflammatory immune responses (26). For energy to be supplied to the body, cardiac metabolism is essential (27). Energy is produced primarily through the mitochondrial respiratory chain. In the Mfuzz analysis, GO analysis of the expression levels of genes that peaked at 24 hpi revealed that these genes were enriched in biological processes associated with processes that respond to ERS, autophagy and the use of autophagic mechanisms. ERS and autophagy played an important role in infectious pneumonia. GO analysis linked the light green module, which was found to be strongly associated with 24 hpi in WGCNA, to ERS. Infections that cause inflammation are closely linked to ERS (28). In infected cells, autophagy occurs downstream of ERS signals (29). GO analysis showed that genes for biological processes enriched in Clusters 4 and 9, in which gene expression levels were reduced, were associated with cell adhesion, cell junctions and smooth muscle cell contraction. These findings suggest that intercellular adhesion may be compromised, which may be related to increased permeability of the alveolar-capillary barrier and the development of pulmonary edema.

In the immune cell infiltration analysis, the levels of neutrophils, monocytes, and macrophages were consistently increased after SMA infection. In response to bacterial infection, innate immune cells, such as neutrophils, monocytes, and macrophages become activated during sepsis, resulting in an inflammatory response and recruitment of immune cells into tissues (30). When neutrophils are overrecruited and activated following severe infection, the exaggerated inflammatory response can lead to severe tissue damage and result in organ failure (31). This was reflected in our pathological findings, which showed that after infection, damage to lung tissue structures occurred after lung immune cells were significantly activated and recruited. However, the number of NK cells decreased from 24 hpi. In contrast to the response to SMA infection, the levels of NK cells remain unchanged or increase during the majority of bacterial lung infections (32-35). When pathogens infiltrate the lungs, NK cells play a critical role in the immune response, and they are closely associated with fighting bacterial infections. The response to the host may be bidirectional (36,37). As a result, we hypothesize that SMA lung infection reduces the levels of NK cells, resulting in reduced inflammatory cell infiltration. The role of NK cells in lethal SMA pneumonia, which could lead to new therapeutic targets, warrants further investigation and clarification.

Several of the top 20 genes upregulated during infection are linked to gram-negative bacteria-induced inflammatory responses, while some are linked to innate immune responses. Cxcl10 encodes interferon gamma (IFN-γ) inducible protein 10 (IP-10), a neutrophil chemotactic factor (38). CD14 binds to monomeric lipopolysaccharides and transports them to the LY96/TLR4 complex, where they are processed and used to mediate the innate immune response to bacterial lipopolysaccharides (39). Gbp5, an activator of Nlrp3 inflammasome assembly, is involved in innate immunity and inflammation (40). CXCR2 binds to IL-8, CXCL3, growth regulated protein (GRO)/melanoma growth stimulatory activity (MGSA) and neutrophil-activating peptide 2 (NAP-2) with high affinity and plays an important role in neutrophil recruitment (41). EGR1, a nucleic zinc finger transcription factor, regulates inflammatory and immune responses (42). In a model of skin infection, phenol-soluble modulins (PSMs), produced and secreted by Staphylococcus aureus, were found to directly activate Egr1 and rapidly recruit leukocytes to the site of infection (43). The Sell gene encodes l-selectin, which is responsible for recruiting circulating leukocytes to endothelial inflammation sites (44). By hydrolyzing guanosine triphosphate (GTP) to guanosine monophosphate (GMP), Gbp2 increases oxidative killing and delivers antimicrobial peptides to autophagolysosomes (45). Given the continuing elevated expression levels and multiple functions of these genes in SMA-infected lungs, as well as the fact that they have not previously been reported in SMA pneumonia, they are likely to be suitable therapeutic targets for intervention in SMA-induced lung injury and warrant further investigation.

Conclusions

Our work is the first study to investigate the pulmonary transcriptional response to SMA infection. The findings suggest that SMA causes acute lethal inflammation in the lung that worsens over time, causing lung structure damage. NK cells are expected to be a new therapeutic target for SMA pneumonia after attention is drawn to their possible role in lethal SMA pneumonia. Cxcl10, Cd14, Gbp5, Cxcr2, Tnip1, Zc3h12a, Egr1, Sell, and Gbp2, which have consistently high or unreported levels of expression in SMA pneumonia, may be important targets for future research. And, the roles of infiltrated NK cells and mechanism of recruitment of NK cells in protecting host is worth of further study. Finally, our findings validate the findings of previous studies on SMA pneumonia while also providing new insights into SMA pneumonia that hold promise for the development of therapeutic approaches to reduce SMA lung inflammation.

Acknowledgments

Funding: None.

Footnote

Reporting Checklist: The authors have completed the ARRIVE reporting checklist. Available at https://jtd.amegroups.com/article/view/10.21037/jtd-23-1138/rc

Data Sharing Statement: Available at https://jtd.amegroups.com/article/view/10.21037/jtd-23-1138/dss

Peer Review File: Available at https://jtd.amegroups.com/article/view/10.21037/jtd-23-1138/prf

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://jtd.amegroups.com/article/view/10.21037/jtd-23-1138/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. Animal experiments were performed under a project license (No. IACUC-DWZX-2022-064) granted by Animal Care and Use Committee of the Academy of Military Medical Sciences, in compliance with institutional guidelines for the care and use of animals.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Brooke JS. Advances in the Microbiology of Stenotrophomonas maltophilia. Clin Microbiol Rev 2021;34:e0003019. [Crossref] [PubMed]
2. Jian J, Xie Z, Chen L. Risk Factors for Mortality in Hospitalized Patients with Stenotrophomonas maltophilia Bacteremia. Infect Drug Resist 2022;15:3881-6. [Crossref] [PubMed]
3. Wang N, Tang C, Wang L. Risk Factors for Acquired Stenotrophomonas maltophilia Pneumonia in Intensive Care Unit: A Systematic Review and Meta-Analysis. Front Med (Lausanne) 2021;8:808391. [Crossref] [PubMed]
4. Imoto W, Yamada K, Yamairi K, et al. Clinical Characteristics of Rapidly Progressive Fatal Hemorrhagic Pneumonia Caused by Stenotrophomonas maltophilia. Intern Med 2020;59:193-8. [Crossref] [PubMed]
5. Pek Z, Cabanilla MG, Ahmed S. Treatment refractory Stenotrophomonas maltophilia bacteraemia and pneumonia in a COVID-19-positive patient. BMJ Case Rep 2021;14:e242670. [Crossref] [PubMed]
6. Imoto W, Yamada K, Kuwabara G, et al. In which cases of pneumonia should we consider treatments for Stenotrophomonas maltophilia? J Hosp Infect 2021;111:169-75. [Crossref] [PubMed]
7. Puech B, Canivet C, Teysseyre L, et al. Effect of antibiotic therapy on the prognosis of ventilator-associated pneumonia caused by Stenotrophomonas maltophilia. Ann Intensive Care 2021;11:160. [Crossref] [PubMed]
8. Xu G, Tang X, Shang X, et al. Identification of immunogenic outer membrane proteins and evaluation of their protective efficacy against Stenotrophomonas maltophilia. BMC Infect Dis 2018;18:347. [Crossref] [PubMed]
9. Menetrey Q, Sorlin P, Jumas-Bilak E, et al. Achromobacter xylosoxidans and Stenotrophomonas maltophilia: Emerging Pathogens Well-Armed for Life in the Cystic Fibrosis Patients' Lung. Genes (Basel) 2021;12:610. [Crossref] [PubMed]
10. Molloy K, Cagney G, Dillon ET, et al. Impaired Airway Epithelial Barrier Integrity in Response to Stenotrophomonas maltophilia Proteases, Novel Insights Using Cystic Fibrosis Bronchial Epithelial Cell Secretomics. Front Immunol 2020;11:198. [Crossref] [PubMed]
11. Wu CM, Huang HH, Li LH, et al. Molecular Characterization of Three Tandemly Located Flagellin Genes of Stenotrophomonas maltophilia. Int J Mol Sci 2022;23:3863. [Crossref] [PubMed]
12. Zgair AK, Chhibber S. Stenotrophomonas maltophilia flagellin restricts bacterial colonization in BALB/c mouse lung in vivo. FEMS Immunol Med Microbiol 2012;66:191-200. [Crossref] [PubMed]
13. Fratoni AJ, Nicolau DP, Kuti JL. Minocycline pharmacodynamics against Stenotrophomonas maltophilia in the neutropenic murine infection model: implications for susceptibility breakpoints. J Antimicrob Chemother 2022;77:1052-60. [Crossref] [PubMed]
14. Imoto W, Kaneko Y, Yamada K, et al. A mouse model of rapidly progressive fatal haemorrhagic pneumonia caused by Stenotrophomonas maltophilia. J Glob Antimicrob Resist 2020;23:450-5. [Crossref] [PubMed]
15. Han P, Zhang W, Pu M, et al. Characterization of the Bacteriophage BUCT603 and Therapeutic Potential Evaluation Against Drug-Resistant Stenotrophomonas maltophilia in a Mouse Model. Front Microbiol 2022;13:906961. [Crossref] [PubMed]
16. Hu X, Yu Y, Feng J, et al. Pathologic changes and immune responses against Coxiella burnetii in mice following infection via non-invasive intratracheal inoculation. PLoS One 2019;14:e0225671. [Crossref] [PubMed]
17. Price A, Caciula A, Guo C, et al. DEvis: an R package for aggregation and visualization of differential expression data. BMC Bioinformatics 2019;20:110. [Crossref] [PubMed]
18. Hisatomi A, Shiwa Y, Fujita N, et al. Identification and structural characterisation of a catecholate-type siderophore produced by Stenotrophomonas maltophilia K279a. Microbiology (Reading) 2021; Erratum in: Microbiology (Reading) 2021 doi: 10.1099/mic.0.001081. [Crossref]
19. Li Q, Chen JX, Wu Y, et al. The mechanism of FZXJJZ decoction suppresses colorectal liver metastasis via the VDR/TGF-β/Snail1 signaling pathways based on network pharmacology-TCGA data-transcriptomics analysis. J Ethnopharmacol 2022;287:114904. [Crossref] [PubMed]
20. Liu K, Chen S, Lu R. Identification of important genes related to ferroptosis and hypoxia in acute myocardial infarction based on WGCNA. Bioengineered 2021;12:7950-63. [Crossref] [PubMed]
21. Miao YR, Zhang Q, Lei Q, et al. ImmuCellAI: A Unique Method for Comprehensive T-Cell Subsets Abundance Prediction and its Application in Cancer Immunotherapy. Adv Sci (Weinh) 2020;7:1902880. [Crossref] [PubMed]
22. Di Bonaventura G, Pompilio A, Zappacosta R, et al. Role of excessive inflammatory response to Stenotrophomonas maltophilia lung infection in DBA/2 mice and implications for cystic fibrosis. Infect Immun 2010;78:2466-76. [Crossref] [PubMed]
23. Miao YR, Xia M, Luo M, et al. ImmuCellAI-mouse: a tool for comprehensive prediction of mouse immune cell abundance and immune microenvironment depiction. Bioinformatics 2022;38:785-91. [Crossref] [PubMed]
24. Barnabei L, Laplantine E, Mbongo W, et al. NF-κB: At the Borders of Autoimmunity and Inflammation. Front Immunol 2021;12:716469. [Crossref] [PubMed]
25. Caetano-Silva ME, Rund LA, Vailati-Riboni M, et al. Copper-Binding Peptides Attenuate Microglia Inflammation through Suppression of NF-kB Pathway. Mol Nutr Food Res 2021;65:e2100153. [Crossref] [PubMed]
26. Yuan Q, Zeng ZL, Yang S, et al. Mitochondrial Stress in Metabolic Inflammation: Modest Benefits and Full Losses. Oxid Med Cell Longev 2022;2022:8803404. [Crossref] [PubMed]
27. Hu X, Ou-Yang Q, Wang L, et al. AdipoRon prevents l-thyroxine or isoproterenol-induced cardiac hypertrophy through regulating the AMPK-related pathway. Acta Biochim Biophys Sin (Shanghai) 2019;51:20-30. [Crossref] [PubMed]
28. Sweet LA, Kuss-Duerkop SK, Keestra-Gounder AM. IRE1α-Driven Inflammation Promotes Clearance of Citrobacter rodentium Infection. Infect Immun 2022;90:e0048121. [Crossref] [PubMed]
29. Cao X, Wan H, Wan H. Urolithin A induces protective autophagy to alleviate inflammation, oxidative stress, and endoplasmic reticulum stress in pediatric pneumonia. Allergol Immunopathol (Madr) 2022;50:147-53. [Crossref] [PubMed]
30. Elemam NM, Ramakrishnan RK, Hundt JE, et al. Innate Lymphoid Cells and Natural Killer Cells in Bacterial Infections: Function, Dysregulation, and Therapeutic Targets. Front Cell Infect Microbiol 2021;11:733564. [Crossref] [PubMed]
31. Burn GL, Foti A, Marsman G, et al. The Neutrophil. Immunity 2021;54:1377-91. [Crossref] [PubMed]
32. Schuijs MJ, Png S, Richard AC, et al. ILC2-driven innate immune checkpoint mechanism antagonizes NK cell antimetastatic function in the lung. Nat Immunol 2020;21:998-1009. [Crossref] [PubMed]
33. Pan Z, Xu ML. T-cell and NK-cell lymphomas in the lung. Semin Diagn Pathol 2020;37:273-82. [Crossref] [PubMed]
34. Esaulova E, Das S, Singh DK, et al. The immune landscape in tuberculosis reveals populations linked to disease and latency. Cell Host Microbe 2021;29:165-178.e8. [Crossref] [PubMed]
35. Witkowski M, Tizian C, Ferreira-Gomes M, et al. Untimely TGFβ responses in COVID-19 limit antiviral functions of NK cells. Nature 2021;600:295-301. [Crossref] [PubMed]
36. Frank K, Paust S. Dynamic Natural Killer Cell and T Cell Responses to Influenza Infection. Front Cell Infect Microbiol 2020;10:425. [Crossref] [PubMed]
37. Bergantini L, d'Alessandro M, Cameli P, et al. NK and T Cell Immunological Signatures in Hospitalized Patients with COVID-19. Cells 2021;10:3182. [Crossref] [PubMed]
38. Li M, Chen Y, Li H, et al. Serum CXCL10/IP-10 may be a potential biomarker for severe Mycoplasma pneumoniae pneumonia in children. BMC Infect Dis 2021;21:909. [Crossref] [PubMed]
39. Ciesielska A, Krawczyk M, Sas-Nowosielska H, et al. CD14 recycling modulates LPS-induced inflammatory responses of murine macrophages. Traffic 2022;23:310-30. [Crossref] [PubMed]
40. Ding K, Li X, Ren X, et al. GBP5 promotes liver injury and inflammation by inducing hepatocyte apoptosis. FASEB J 2022;36:e22119. [Crossref] [PubMed]
41. Pavlidis P, Tsakmaki A, Pantazi E, et al. Interleukin-22 regulates neutrophil recruitment in ulcerative colitis and is associated with resistance to ustekinumab therapy. Nat Commun 2022;13:5820. [Crossref] [PubMed]
42. Banerji R, Saroj SD. Early growth response 1 (EGR1) activation in initial stages of host-pathogen interactions. Mol Biol Rep 2021;48:2935-43. [Crossref] [PubMed]
43. Nguyen TH, Cheung GYC, Rigby KM, et al. Rapid pathogen-specific recruitment of immune effector cells in the skin by secreted toxins. Nat Microbiol 2022;7:62-72. [Crossref] [PubMed]
44. Mocanu CA, Fuior EV, Voicu G, et al. P-selectin targeted RAGE-shRNA lipoplexes alleviate atherosclerosis-associated inflammation. J Control Release 2021;338:754-72. [Crossref] [PubMed]
45. Zhang Y, Liao Y, Hang Q, et al. GBP2 acts as a member of the interferon signalling pathway in lupus nephritis. BMC Immunol 2022;23:44. [Crossref] [PubMed]