Single-cell and spatial transcriptomics reveal a stress-induced EMT-like epithelial subset driving immune activation in silica-injured lung

The mechanism that lung epithelial cells regulate immune responses during chronic injury still remains unclear. Here, we combined single-cell RNA sequencing with spatial transcriptomics to track epithelial dynamics in silica (SiO₂)-exposed mouse lungs. By day 56, SiO₂ induced significant epithelial proliferation, followed with a distinct C0 subset emerging as the dominant population. C0 cells co-expressed surfactant genes (Sftpc, Scgb3a2), mesenchymal markers (Vim, Mmp12), and pro-inflammatory cytokines (Ccl6, S100a8/a9), reflecting a hybrid phenotype. Spatial and cell-cell interaction analyses showed C0 cells engaging macrophages and neutrophils through SPP1-CD44, APP-CD74, and GRN-MARCO signaling, driving immune recruitment and activation. Pseudotime and CytoTRACE analyses indicated that C0 cells represent a late-stage, low-stemness state with epithelial-mesenchymal transition (EMT)-like features. Taken together, these findings reveal a novel, stress-induced epithelial subset that amplifies immune crosstalk and tissue remodeling, offering new perspectives on silica-induced lung injury.

1 Introduction

Lung injury and inflammation resulting from environmental or occupational exposure to silica (SiO₂) remain significant challenges in clinical practice and scientific research, particularly in industrial and mining settings. Silica dust, a pervasive environmental pollutant, is a well-known cause of lung diseases such as pneumoconiosis (1), which often leads to pulmonary fibrosis, a disease characterized by excess collagen buildup and permanent lung scarring. Globally, tens of millions of workers remain exposed to silica dust across various industries, including mining, construction, and engineered stone manufacturing. For instance, China reported approximately 9,000 new cases annually due to extensive silica exposure in the mining and construction industries, while India and Brazil have over 10 million and 6 million workers respectively who are regularly exposed to silica dust (2). Moreover, occupational silicosis is responsible for more than 12,900 annual deaths and approximately 0.65 million disability-adjusted life years worldwide (2).

At the cellular and molecular levels, SiO₂ induces apoptosis in lung epithelial cells and stimulates epithelial cells to release exosomes that promote collagen deposition and activate lung fibroblasts (3, 4). Another effect of SiO₂ is that it induces a proinflammatory response in pulmonary immune cells. For instance, exposure to SiO₂ suppresses pulmonary T-cell responses and induces a shift in macrophage polarization from a resting (M0) state to a proinflammatory (M1) phenotype (5–7). Macrophages, through the phagocytosis of SiO₂ particles, trigger a cascade of proinflammatory and profibrotic responses. More recently, lung epithelial cells have gained attention for their role in innate and adaptive immunity, mediated by the secretion of cytokines (8). Despite these insights, the precise mechanisms underlying epithelial cell behavior and their interplay with immune cells following SiO₂ exposure remain poorly understood, underscoring the need for further research.

Central to unraveling these mechanisms is the role of epithelial cells in the lung’s response to SiO₂-induced stress. The alveolar epithelium, a critical structural and functional barrier, comprises two main cell types: type I alveolar cells (AT1), which mediate gas exchange across most of the alveolar surface, and type II alveolar epithelial cells (AT2), which produce surfactant to prevent alveolar collapse and support repair by proliferating and differentiating into AT1 cells following injury (9, 10). Repeated epithelial injury from chronic SiO₂ exposure disrupts these repair processes, potentially resulting in aberrant cellular responses such as epithelial-to-mesenchymal transition (EMT), which further promotes fibrosis. Notably, epithelial cells exposed to stressors like SiO₂ express specific gene programs that profoundly influence their interactions with immune cells, balancing inflammatory and repair processes. For example, in the lung injury model, exogenous SCGB3A2 was shown to accelerate tissue repair and mitigate fibrosis by modulating TGF-β signaling (11, 12). Similarly, other evidence indicates that SPP1, highly expressed in the lung tissue of pulmonary fibrosis patients, promotes M2 macrophage polarization and induces EMT-like changes in epithelial cells (13, 14). However, how epithelial cells regulate immune responses via these signaling pathways remains unclear and warrants further investigation.

Recent technological breakthroughs, notably single-cell RNA sequencing and spatial transcriptomics, have provided powerful tools to dissect the complex cellular dynamics in lung tissue after SiO₂ exposure. In this study, we aim to elucidate the cellular and molecular mechanisms underlying SiO₂-induced lung injury, with a particular focus on epithelial cells and their interactions with immune cells. Leveraging advanced approaches such as single-cell RNA sequencing and spatial transcriptomics, we sought to track epithelial cell dynamics and gene expression alterations following chronic SiO₂ exposure. We hypothesized that chronic SiO₂-induced epithelial injury could lead to the emergence of novel epithelial subsets that undergo partial EMT, while orchestrating inflammation and fibrosis through direct interactions with immune cells.

These findings shed light on a new epithelial-immune mechanism in SiO₂-induced lung injury, highlighting the pivotal role of C0 cells in orchestrating immune responses and contributing to disease progression. Given the persistent global health burden of silicosis, particularly in highly industrialized and developing regions, insights into these mechanisms could pave the way for novel therapeutic strategies targeting epithelial–immune cell interactions in fibrotic lung diseases.

2 Materials and methods

2.1 Preprocessing of single-cell transcriptomes

Single-cell RNA sequencing data were obtained in matrix format from the Gene Expression Omnibus (GEO) under accession number GSE183682. Lung tissues were harvested from mice of normal saline (NS) and SiO₂ groups at 7 days and 56 days (mouse and modeling methods as described in Chen et al.), immediately frozen in OCT on dry ice, and stored at -80°C before further processing (15). Quality control was performed using Seurat and DoubletFinder to identify and remove multiplets and low-quality cells. The Seurat package was also used to filter out foreign cells, where a gene was considered expressed if detected in more than 3 cells, and each cell was required to have at least 200 expressed genes. Following these quality control procedures, a total of 36,122 high-quality single cells were retained for downstream analysis. Gene expression data were subsequently normalized using Seurat’s “LogNormalize” method to reduce variability in gene expression counts.

2.2 Dimensionality reduction, clustering, and annotation

To identify highly variable genes (HVGs), Seurat’s “FindVariableFeatures” function was applied with default parameters. Principal component analysis (PCA) was then performed on the HVGs, and the top 30 principal components were retained for subsequent analysis. Uniform manifold approximation and projection (UMAP) was used to reduce the dimensionality of the data and to visualize the cellular clusters in two dimensions. Clustering of the cells was performed using Seurat’s “FindClusters” function. To correct for potential batch effects, the Harmony R package was employed to integrate datasets across different conditions.

2.3 Intracellular crosstalk analysis

Cell-cell interactions were analyzed using the CellChat R package, which relies on known ligand-receptor pairs to explore intercellular communication. The normalized counts of the merged samples were pre-processed using the identifyOverExpressedGenes, identifyOverExpressedInteractions, and projectData functions from the CellChatDB.human database. To identify significant ligand-receptor pairs, key CellChat functions, including computeCommunProb, computeCommunProbPathway, and aggregationNet, were applied. Interaction networks were visualized using CellChat’s netVisual_bubble function to assess the strength and significance of the interactions.

2.4 Differential expression and enrichment analyses

The “FindMarkers” function in Seurat was used to identify differentially expressed genes (DEG) in specific cell subgroups, with the following cutoff criteria: an average log fold change (Avg_logFC) greater than 0.25 and an adjusted p-value (p_val_adj) less than 0.05. To provide further biological insights into the DEGs, functional enrichment analysis was conducted using the ClusterProfiler R package to perform Gene Ontology (GO) analyses.

2.5 Trajectory analysis

Trajectory analysis was conducted using the Monocle 2.0 R package to investigate the differentiation processes of epithelial cells. The “differentialGeneTest” function was used to identify genes associated with pseudotime.

2.6 CellTrek analysis

To further explore the spatial transcriptomic data in conjunction with single-cell RNA-seq, the feature-barcode matrix was downloaded from GEO (accession number GSE183683). The “traint” function in the CellTrek R package was applied to co-embed the spatial transcriptomics (ST) data with the single-cell RNA-seq data. Subsequently, the celltrek function, using default parameters, was employed to project single cells onto the ST coordinates.

2.7 EMT score calculation

To score EMT related gene sets of each single cell in the scRNA-Seq data, the ‘AddModuleScore’ method from the Seurat package was utilized. Scores were calculated by expressions of EMT related genes noted in result section.

2.8 Animal models and ethical compliance

Male C57BL/6 mice (8–12 weeks; body weight range 20–24 g) were purchased from Cavion Experimental Animal Co., LTD., Changzhou City, Jiangsu Province, China (animal license number SCXY (Su) 2011-0003). Animals were raised under the condition of no specific pathogen (SPF), with a temperature of 22-24 °C, relative humidity of 50 ± 5%, and a day-night cycle of 12 hours of light and 12 hours of darkness (lights at 08:00). After a 7–14 day adaptation period, all mice were randomly divided into two groups of 12 mice each: the NS group (NS) and the silica exposure group (SiO₂).

The mice were then anesthetized by inhaling isoflurane and intranasal instillation of silica (20ug/ul, 80ul) diluted in saline was used (SiO₂). The animals in the NS group were given the same dose of normal saline nasal drops, and the mice were sacrificed at the corresponding time points, and then the lung tissues were collected and stored for further experiments.

All procedures followed the Guide for the Care and Use of Laboratory Animals (NIH Publication No. 8023, revised 1978) and were approved by the Institutional Animal Care and Ethics Committee of Anhui University of Science and Technology.

2.9 Immunohistochemistry and immunofluorescence staining

For immunohistochemistry (IHC) staining, paraffin-embedded lung tissue sections were first cleared of paraffin with xylene and rehydrated in a series of ethanol dilutions. To expose antigens, sections were heated in Citrate Antigen Retrieval Solution (Beyotime, P0090) for 30 minutes. Endogenous peroxidase activity was blocked by treating the slides with 3% hydrogen peroxide for 10 minutes. Permeabilization was conducted with 0.3% Triton^® X-100 (BioFroxx, 1139) at room temperature for 10 minutes. After a 1-hour block with 3% BSA-PBST, primary antibodies were applied and left on the sections overnight at 4°C. The following day, HRP-linked secondary antibodies were introduced for 1 hour at room temperature, and staining was developed using 3,3’-diaminobenzidine (DAB). Hematoxylin was used to stain nuclei, after which sections were dehydrated, coverslipped, and examined under a light microscope.

For immunofluorescence staining, paraffin sections of mice lungs were dewaxed, rehydrated in a series of ethanol dilutions, immersion in hydrogen peroxide, repairing antigen, permeabilization, and then blocking for 1 hour in 10% goat serum at room temperature. Subsequently, the primary antibody was incubated overnight in a 4°C refrigerator after staining, rewarmed at room temperature, and then stained with the secondary antibody at room temperature for 1 hour. Finally, nuclei were stained with DAPI (1:1000) for 10 minutes. All tissue sections were mounted on glass slides and covered with an anti-fade mounting medium. Images were captured using a confocal microscope (Olympus, FV-3000). The antibodies used in this experiment are listed in Table 1.

Table 1

Table 1. Details of antibodies used in IHC and IF experiments.

2.10 Statistical analysis

Data are presented as mean ± standard error of the mean (SEM). Homogeneity of variances across groups was examined using the Brown-Forsythe and Bartlett’s tests. Comparisons involving two groups were performed using an independent-samples t-test. For comparisons involving three or more groups, two-way ANOVA was conducted, followed by Bonferroni’s multiple comparisons test. GraphPad Prism 9.0.0 software was used to perform the statistical analyses, where p< 0.05 was set as the level of significance for determining any significant differences.

3 Results

3.1 Marked epithelial cell proliferation at Day 56 upon SiO₂ treatment

To investigate the effects of SiO₂ exposure on lung tissues, we first performed dimensionality reduction using UMAP to visualize the distribution of different cell types across various experimental groups (NS_7d, SiO₂_7d, NS_56d, SiO₂_56d) using data obtained from the GEO under accession number GSE183682 (Figure 1a). We utilized Seurat for preprocessing the single-cell RNA sequencing data, including quality control and normalization. UMAP plots revealed distinct cell clusters, the SiO₂-treated group at 56 days exhibited markedly greater changes compared to other groups. In particular, significant alterations were observed not only within immune cell populations but also prominently among epithelial cells (Figure 1b). This marked epithelial cell response prompted us to further investigate their potential roles in mediating SiO₂-induced lung injury.

Figure 1

Figure 1. Single-cell RNA sequencing analysis of lung tissue in SiO₂-induced injury. (a) UMAP plot showing the distribution of different cell types across various treatment groups (NS_7d, SiO₂_7d, NS_56d, SiO₂_56d). (b) The proportion of different cell types in each treatment group. Epithelial cell percentages are as follows: NS_7d (1.23%), SiO₂_7d (1.25%), NS_56d (8.00%), SiO₂_56d (10.00%). (c) GO enrichment analysis of differentially expressed genes in epithelial cells. (d-g) UMAP plots (left) and spatial distribution of epithelial cells (right) for each group at 7 and 56 days post-treatment: NS_7d (d), SiO₂_7d (e), NS_56d (f), and SiO₂_56d (g). (h, i) Immunofluorescence staining and quantification of Sftpc (green) and CD206 (red) in NS and SiO₂ group lung tissues, with DAPI (blue) nuclear staining. Significance levels: P **** p ≤ 0.0001; P *** p ≤ 0.001; P ** p ≤ 0.01; P * p≤ 0.05; NS P>0.05.

The increased number of epithelial cells in the SiO₂-treated group at 56 days suggests a possible involvement of cellular proliferation and repair processes induced by SiO₂ exposure. GO enrichment analysis of differentially expressed genes in epithelial cells from the disease groups revealed significant upregulation of genes involved in epithelial cell migration, proliferation and responses to external stimuli (Figure 1c). These findings indicate that SiO₂ exposure enhances epithelial cell expansion, likely through heightened proliferative and stress-responsive activities.

We then combined spatial transcriptomics with single-cell RNA sequencing, which revealed distinct patterns in the distribution and gene expression of various cell types across the tissue sections (Supplementary Figure 1). In the 7-day treatment groups, both NS and SiO₂ exposures exhibited sparse distribution of epithelial cells across the tissue. By day 56, however, epithelial cells markedly increased in abundance and formed distinct clusters, especially within SiO₂-exposed tissues (Figures 1d-g). Additionally, spatial analysis revealed that macrophages were more abundant in the SiO₂-treated tissues at 56 days, which coincided with the increase in epithelial cells (Supplementary Figure 1).

To further validate these observations, we performed immunofluorescence staining on lung tissue sections from the NS_56d and SiO₂_56d groups (Figures 1h-i). In the NS_56d group, Sftpc⁺ epithelial cells (green) were evenly distributed with moderate intensity, accompanied by a sparse presence of CD206⁺ M2 macrophages (red). In contrast, the SiO₂_56d group exhibited a significant increase in Sftpc⁺ epithelial cells and CD206⁺ macrophages. Quantitative analysis confirmed that both Sftpc and CD206 signals were significantly elevated in the SiO₂_56d group compared to NS_56d. Interestingly, Sftpc-positive cells were found to be localized around CD206-positive macrophages in the SiO₂-treated group, suggesting a coordinated expansion of epithelial and macrophage populations in response to SiO₂-induced injury (Figure 1h).

Taken together, these results revealed that SiO₂ exposure leads to significantly enhanced epithelial cell proliferation and a heightened immune response. These findings imply that SiO₂ exposure may induce a complex tissue response, involving both immune activation and epithelial cell proliferation.

3.2 Spatial transcriptomics confirms late-stage expansion of C0 cells with elevated Sftpc and Ccl6

Next, we aimed to investigate the relationship between different subclusters. The gene expression data show that the epithelial cells were classified into 5 distinct groups, with SiO₂ exposure significantly altering their composition (Figures 2a-c). The UMAP analysis (Figure 2b) shows that at 7 days, the NS and SiO₂ groups have similar distributions of epithelial subpopulation. However, we noticed C0 cells a marked shift in the SiO₂_56d group to dominate. This was further confirmed by the proportional distribution of epithelial subpopulations (Figure 2c), where C0 cells were significantly more abundant in the SiO₂_56d group compared to other groups.

Figure 2

Figure 2. Characterization of epithelial subpopulations in SiO₂-induced lung injury. (a) UMAP plot showing the distribution of epithelial cell subpopulations across all treatment groups. (b) UMAP plot focusing on the distribution of epithelial subpopulations (C0-C5) across NS_7d, SiO₂_7d, NS_56d, and SiO₂_56d groups. (c) The proportion of each epithelial subpopulation (C0-C5) in different treatment groups. (d) The expression of differentially expressed genes across epithelial subpopulations (C0-C5). (e) The expression levels of Sftpc, Scgb3a2, Ccl6, and Krt18 across epithelial subpopulations (C0-C5).

Violin plots and heatmap analyses revealed distinct gene expression signatures for each cell cluster (Figures 2d, e). Cluster C1 cells, marked by high expression of Foxj1 and Tppp3, were identified as ciliated airway epithelial cells. C2 cells exhibited elevated Scgb1a1, a canonical marker of secretory club cells. C3 cells exhibited co-expressing Krt17 (progenitor/stress marker) and Ager (AT1 marker), suggesting that this population may represent transitional epithelial cells committed toward the AT1 lineage (16). C3 cells also expressed Cldn18, Cldn4, and Bcam, indicating roles in tight junction formation and barrier maintenance within the respiratory epithelium. C4 cells were identified as Lamp3⁺ and Sftpc⁺ AT2 cells. C4 cells also expressed innate immune response genes, including Lyz1, Lyz6, Cxcl15, and Sfta2. Lastly, C5 cells, another secretory population, displayed high levels of Scgb3a2 and Scgb3a1, along with Calml3, Chad, and Aldh3a1, suggesting functions in epithelial cell development and response to wounding.

In contrast, C0 cells expressed epithelial markers such as Scgb3a2 and Sftpc but were distinguished by EMT-related genes, including Spp1, Mmp12, Ctss, and Vim. Unlike other clusters, C0 cells lacked Krt18, a marker of epithelial differentiation expressed across all other subpopulations, suggesting a mesenchymal transition (17–19). Further analysis revealed high expression of inflammatory activation markers S100a8 and S100a9, alongside immune- and inflammation-related genes such as Ccl6, Cxcl2, and Il1b. These genes, linked to chemotaxis and immune cell recruitment (5, 20–22), collectively indicate that C0 cells actively contribute to immune modulation and inflammatory responses in the SiO₂_56d group.

To map these dynamics over time, we employed spatial transcriptomics. At 7 days post-exposure, both the NS and SiO₂ groups exhibited sparse C0 cell populations (Figures 3a, b). By 56 days, however, C0 cells were significantly more abundant in the SiO₂ group (Figures 3c, d). To assess C0 cell dynamics, we examined the spatial expression of Sftpc and Ccl6 (Figure 3c). SiO₂_56d tissue showed elevated and localized Sftpc expression, while Ccl6 expression was markedly elevated and widespread (Figure 3d). Immunohistochemical staining further confirmed these findings (Figures 3e, f). In SiO₂_56d lungs, Spp1 expression was markedly elevated, co-occurring with a significant increase in Sftpc-positive cells compared to controls. Quantitative analysis showed significant increases in both markers’ area percentages, consistent with a partial EMT and inflammatory signatures in SiO₂-induced injury.

Figure 3

Figure 3. Spatial distribution of epithelial subpopulations in SiO₂-induced lung injury. (a-d) Spatial expression of Sftpc (left) and Ccl6 (right) in epithelial subpopulations (C0-C5) across different time points: (a) NS_7d, (b) SiO₂_7d, (c) NS_56d, and (d) SiO₂_56d. (e, f) Immunohistochemistry staining and quantification of Spp1 (brown) and Sftpc (brown) in NS and SiO₂-treated lung tissues. P **** p ≤ 0.0001; P *** p ≤ 0.001.

Together, these findings illuminate the dynamic epithelial response to SiO₂ exposure. The C0 cell population expresses epithelial markers such as Sftpc and Scgb3a2 but is distinguished by co-expression of EMT-related genes like Mmp12 and Vim, indicating an ongoing transition toward a mesenchymal phenotype. Notably, C0 cells also upregulate a suite of immune-related genes—including Ccl6, S100a8, and Cxcl2—driving immune cell recruitment and amplifying inflammatory responses. In contrast, other clusters maintain distinct epithelial identities while contributing to cellular defense. For instance, the Lamp3⁺ and Sftpc⁺ AT2 cells (C4) express innate immune response genes such as Lyz1 and Cxcl15, supporting antimicrobial defense and tissue protection.

3.3 Pseudotime analysis reveals C0 emergence from C1

Having established the expansion and phenotypic characteristics of C0 cells, we next sought to investigate their developmental origin. To achieve this, we employed CytoTRACE analysis to assess the differentiation potential of epithelial cell subpopulations. The CytoTRACE analysis plot reveals that C0 cells are at the final stage of differentiation, representing mature or post-differentiation cells, while C1 and other cells are at an earlier, progenitor-like stage (Figure 4a). The stemness score shows that C1 cells exhibit the highest stemness score, with C0 the lowest (Figure 4b). This suggests that C1 cells may be involved in early-stage differentiation or progenitor functions, while C0 cells likely participate in later stages of differentiation or mature epithelial functions. This is further supported by pseudotime analysis of the epithelial cells, both with and without SiO₂ treatment (Figure 4c). For this analysis, we selected C1 cells as the pseudotime starting point, revealing that C0 cells represent the most differentiated state at the end of the trajectory. These findings align with our earlier observation that Krt18 was absent in C0 cells, further suggesting that they occupy a more advanced and differentiated state.

Figure 4

Figure 4. Differentiation and trajectory of epithelial subpopulations in SiO₂-induced lung injury. (a) CytoTRACE analysis predicting the differentiation potential of epithelial subpopulations (C0-C5). (b) The predicted ordering and differentiation potential across different epithelial subpopulations. (c) Pseudotime trajectory analysis of epithelial subpopulations at different time points (NS_7d, SiO₂_7d, NS_56d, SiO₂_56d) showing the progression of differentiation from C1 to C0 cells. (d) Spatial expression of epithelial subpopulations (C0, C1) in NS_56d and SiO₂_56d lung tissues.

Next, pseudotime trajectory analysis and differential expression analysis revealed a clear distinction between C0 and C1 cells, as visualized in the pseudotime heatmap (Supplementary Figure 2) and volcano plot (Supplementary Figure 4). C0 cells were characterized by high expression of immune response genes, such as Ccl6, S100a8, S100a9, and Il 1b, which are involved in inflammation and immune signaling. Additionally, C0 cells exhibited elevated levels of EMT markers such as Spp1 and Vim, suggesting their active involvement in tissue repair or remodeling. Over the pseudotime trajectory, we observed significant downregulation of genes associated with epithelial functions, including tight junction components (Cldn4, Cldn18), an undifferentiated epithelial marker (Krt7), canonical surfactant genes (Sftpa1, Sftpb, Sftpd), and secretory epithelial markers (Scgb1a1, Scgb3a1, Scgb3a2). These shifts indicate that C0 cells diverge from traditional epithelial roles, reinforcing their specialization in immune modulation and mesenchymal transition during SiO₂-induced injury.

To further explore the spatial distribution of C0 and C1 cells, we performed spatial transcriptomics on the NS_56d and SiO₂_56d groups. Intriguingly, the data showed that C0 and C1 cells are consistently juxtaposed in both treatment groups, regardless of SiO₂ exposure (Figures 3c, d). This spatial proximity, coupled with pseudotime trajectory analysis (Figure 4d), suggests a potential developmental relationship between these clusters. Specifically, C0 cells exhibit elevated expression of EMT markers like Vim, Spp1, and Mmp12, while C1 cells retain ciliated epithelial traits marked by Foxj1, hinting that C0 cells may differentiate from C1 cells along a trajectory involving mesenchymal transition. We propose that under SiO₂-induced stress, C0 cells emerge from C1 precursors and subsequently undergo partial EMT, as evidenced by their upregulation of immune genes (Ccl6, S100a8) and EMT markers, enabling their migration and expansion to support tissue remodeling and inflammation in response to chronic injury.

3.4 C0 cells exhibit increased immune interactions and upregulation of key ligand–receptor interactions in silica-injured lung

We next investigate the signaling pathways and intercellular communication profiles of C0 cells. Firstly, using the CellChat analysis, we identified a substantial increase in both the number and strength of cell-cell interactions in the SiO₂_56d group compared to the NS_56d group. C0 cells in the disease group exhibited significantly enhanced interaction strength and a greater number of interactions with immune cells, particularly macrophages and neutrophils, with a smaller portion directed towards cycling cells (Figures 5a–d).

Figure 5

Figure 5. Cell-cell interactions in SiO₂-induced lung injury. (a, b) The differential number of interactions between epithelial subpopulations and other cell types. (c, d) The differential interaction strength between epithelial subpopulations and other cell types. (e, f) Interaction network of C0 cells with macrophages and other immune cells in SiO₂-treated tissues. (g) Relative information flow between signaling interactions in the SiO₂-treated group. (h) The differential signaling interactions (SPP1, APP, GRN, and SCGB3A2) and their interactions with different cell types.

Further analysis of ligand-receptor interactions showed a marked increase in relative information flow for APP, SPP1, SCGB3A2, and GRN signaling in the SiO₂_56d group (Figures 5e–g). These pathways were categorized by their target cells: SPP1 and APP predominantly engaged macrophages and neutrophils, while SCGB3A2 and GRN mainly interacted with cycling cells (Figure 5h). Expression of SPP1 and APP correlated with pro-inflammatory responses, while SCGB3A2 and GRN aligned with anti-inflammatory and proliferative functions (1, 5, 23). These patterns suggest C0 contributes to both immune activation and tissue repair in SiO₂-induced injury.

Subsequent analysis of spatial transcriptomics data between cycling cells and C0 cells (Figures 6a-c) revealed that, as expected, Scgb3a2 and Grn were predominantly expressed in C0 cells, while their corresponding receptors, Sort1 and Marco, exhibited moderate degree of expression in cycling cells. Notably, both pairs of ligands and receptors were significantly increased in the disease group. Similarly, the ligands SPP1 and APP, along with the receptors CD44 and CD74, were observed to be expressed in C0 cells and macrophage populations, respectively (Figures 6d-f). Furthermore, ligand genes Spp1 and App were significantly upregulated in the SiO₂_56d samples compared to the NS_56d groups (Figures 6d-f), while their ligands, CD44 and CD74, showed a milder upregulation. Additionally, spatial transcriptomics data highlighted a significant rise in CD44 expression in neutrophils (Supplementary Figure 3). Collectively, these findings indicate that C0 cells amplify ligand-driven signaling, via SCGB3A2, SORT, SPP1 and APP, in SiO₂-induced injury, likely enhancing immune cell recruitment and tissue repair through interactions with macrophages and cycling cells.

Figure 6

Figure 6. Spatial expression of key immune and epithelial markers in SiO₂-induced lung injury. (a-c) Spatial expression of immune markers and epithelial markers in NS_56d and SiO₂_56d lung tissues, visualized with spatial transcriptomics. (d-f) Spatial expression of immune-related markers in NS_56d and SiO₂_56d lung tissues.

3.5 Partial EMT dynamics during C1 to C0 transition are enhanced by neutrophil-derived signaling

Previous analyses revealed shifts in epithelial cell populations, with C0 cells emerging as key mediators of immune responses and tissue remodeling after SiO₂ exposure. Interestingly, these cells exhibit mesenchymal traits, suggesting EMT. To explore this, we analyzed EMT-related gene expression along the pseudotime trajectory.

We observed distinct temporal patterns in EMT markers. Mesenchymal genes, including Vim, Fn1, Zeb1, and Dab2, were upregulated, while canonical epithelial markers Epcam and Cdh1 were downregulated as cells progressed toward the C0 phenotype (Figure 7a). These changes peaked in mature C0 cells, indicating a progressive but incomplete EMT. Since C0 cells also expressed differentiated epithelial markers like Sftpc and Scgb3a2, this suggests partial EMT rather than a full mesenchymal transition.

Figure 7

Figure 7. Partial EMT gene expression dynamics and neutrophil-epithelial interactions during pseudotime differentiation. (a) Temporal expression patterns of key EMT-related genes along pseudotime, showing significant upregulation of mesenchymal markers (Vim, Fn1, Zeb1, Dab2 etc.) and downregulation of epithelial markers (Epcam, Cdh1). (b) Quantitative EMT index plotted across pseudotime. (c) UMAP visualization of EMT index distribution. (d) CellChat analysis of ligand-receptor interactions.

To quantify EMT dynamics, we used the AddModuleScore function to calculate an EMT index, which increased along the pseudotime trajectory, confirming a partial EMT from C1 progenitors to C0 cells (Figure 7b). UMAP visualization showed elevated EMT scores in C0 cells, confirming their distinct mesenchymal-like phenotype (Figure 7c).

To investigate drivers of this process, we examined cell-cell interactions using CellChat. Neutrophils exhibited strong ligand-receptor communication with C0 cells via neutrophil-derived ThBS1 interacting with Cd47 and Cd36 receptors on C0 cells (Figure 7d). ThBS1 is known to promote EMT and epithelial plasticity in cancer and tissue regeneration (24, 25).

These findings propose a model where neutrophil-derived ThBS1 drives partial EMT in epithelial cells, promoting differentiation of C1 progenitors into inflammatory, mesenchymal-like C0 cells that retain epithelial markers. This immune-epithelial crosstalk underscores epithelial-mesenchymal plasticity as a dynamic process in SiO₂-induced lung remodeling.

5 Discussion

In this study, we utilized single-cell RNA sequencing and spatial transcriptomics to characterize the temporal and spatial dynamics of lung epithelial cells following SiO₂ exposure. Our analyses identified a distinct epithelial subset, C0 cells, which expanded significantly by day 56 post-exposure. These cells co-expressed alveolar cell markers (Sftpc, Scgb3a2) with mesenchymal transition and immune-related genes (Spp1, Mmp12), yet lacked Krt18, suggesting a shift from traditional epithelial differentiation. These findings position C0 cells as key contributors to immune responses and tissue remodeling in SiO₂-induced lung injury.

C0 cells expressed a robust set of chemokines and cytokines (Ccl6, S100a8, S100a9), alongside CellChat analyses showing enhanced interactions with macrophages and neutrophils. This expression profile, supported by spatial co-localization with CD206⁺ macrophages (Figure 1h), underscores C0’s role in driving the recruitment and infiltration of immune cells in response to SiO₂-induced injury. Pseudotime analysis further revealed that C0 cells likely emerge from Foxj1⁺ C1 precursors, downregulating epithelial markers (Krt7, Cldn4, Scgb1a1) while upregulating EMT genes (Spp1, Vim, Mmp12). These shifts, validated by immunohistochemistry showing elevated Spp1 and Sftpc in SiO₂_56d lungs (Figures 3e-f), suggest that C0 cells undergo partial EMT and enhance their immunomodulatory capacities. Such plasticity may also amplify inflammation and contribute to fibrotic repairment.

Furthermore, our detailed analysis of EMT-related gene expression dynamics along the pseudotime trajectory (Section 3.5) demonstrated upregulation of mesenchymal genes like Vim, Fn1, and Zeb1, as well as downregulation of epithelial markers such as Epcam and Cdh1. Intriguingly, cell-cell interaction analyses identified neutrophil-derived ThBS1 as a critical signaling molecule potentially driving this partial EMT via its receptors Cd47 and Cd36 expressed on C0 cells (Figure 7d). This finding aligns with Nickel et al.’s results that ThBS1 modulates TGF-β signaling to regulate EMT, further supporting Thbs1’s key role in driving epithelial plasticity across different biological contexts (24).

Signaling pathway analysis highlighted key ligand-receptor pairs mediating C0-immune interactions, including SPP1-CD44, APP-CD74, SCGB3A2-SORT1, and GRN-MARCO. In SiO₂_56d, Spp1 and App were markedly upregulated in C0 cells (Figures 6d-f), driving interactions with macrophages and neutrophils, while Scgb3a2 and Grn targeted cycling cells, consistent with their respective pro-inflammatory and repair roles. The SPP1-CD44 axis, previously linked to neutrophil activation in other contexts, may similarly promote immune infiltration here. These networks suggest C0 cells coordinate a dual response, balancing inflammation and repair, though functional assays are needed to confirm these effects.

The hybrid phenotype of C0 cells, combining epithelial and mesenchymal features, echoes transitional epithelial states reported in lung injury models. Choi et al. (26) described damage-associated transient progenitors derived from AT2 cells, which exhibit immune-responsive traits (26). Similarly, Kobayashi et al. (27) identified a pre-alveolar type-1 transitional state with a partial epithelial-mesenchymal phenotype during fibrosis (27). Like these states, C0 cells lack Krt18 and upregulate Spp1 and Vim, serving as an inflammation-associated epithelial subtype under SiO₂ exposure. This convergence suggests that C0 cells may represent a conserved response to chronic lung stress, distinct from the canonical AT2-to-AT1 differentiation pathway.

However, an unresolved question is why these C0 cells maintain or even upregulate epithelial markers such as Sftpc and Scgb3a2 despite their mesenchymal background. Our analyses offer potential insights into this paradox through identifying candidate upstream signaling mechanisms. First, enhanced ligand–receptor interactions involving SCGB3A2 and its receptor SORT1 (Figures 5g-h, Figures 6a-c) suggest that SCGB3A2 may exert autocrine or paracrine signaling to maintain epithelial gene expression within cells undergoing EMT-like transitions (28). Second, the observed neutrophil-derived ThBS1 signaling via receptors Cd47 and Cd36 expressed on C0 cells (Figure 7d) indicates that neutrophil-mediated signaling pathways might further drive this epithelial-mesenchymal plasticity. These findings highlight the complexity of signaling interactions influencing epithelial marker retention despite an ongoing mesenchymal transition, warranting further experimental validation.

Despite these insights, several limitations must be addressed. The lack of in vivo lineage tracing or C0-specific depletion limits our ability to definitively confirm their origin from C1 cells or their direct impact on fibrosis progression. Additionally, reliance on mouse models and omics data raises questions about C0’s relevance to human silicosis, and whether their emergence is reversible remains unexplored. Future studies should employ targeted interventions and extend analyses to human samples to address these gaps.

In conclusion, we identified C0 cells as a novel epithelial subset in SiO₂-induced lung injury, marked by a hybrid epithelial-immune-mesenchymal phenotype. Their expansion, EMT-like features, and immune interactions highlight epithelial plasticity as a driver of inflammation and remodeling in silicosis. Targeting C0-mediated pathways, such as SPP1-CD44, may offer therapeutic strategies to mitigate silica-induced lung diseases.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors.

Ethics statement

The animal study was approved by Institutional Animal Care and Ethics Committee of Anhui University of Science and Technology. The study was conducted in accordance with the local legislation and institutional requirements.

Author contributions

TW: Conceptualization, Data curation, Formal analysis, Investigation, Writing – original draft, Writing – review & editing. JH: Conceptualization, Data curation, Formal analysis, Methodology, Writing – original draft. BL: Conceptualization, Data curation, Formal analysis, Writing – original draft. AH: Writing – original draft. JW: Investigation, Visualization, Writing – original draft. HX: Investigation, Writing – review & editing. QW: Formal analysis, Writing – review & editing. WG: Investigation, Writing – review & editing. CC: Writing – review & editing. CY: Formal analysis, Writing – review & editing. XG: Writing – review & editing. TL: Funding acquisition, Writing – review & editing. MZ: Conceptualization, Funding acquisition, Project administration, Supervision, Writing – review & editing. YX: Conceptualization, Funding acquisition, Project administration, Supervision, Writing – review & editing. XT: Conceptualization, Funding acquisition, Project administration, Supervision, Writing – review & editing.

Funding

The author(s) declare that financial support was received for the research and/or publication of this article. This work was supported by Medical Special Cultivation Project of Anhui University of Science and Technology (No. YZ2023H1A003), Anhui Provincial Special Program for Clinical and Translational Medical Research (No. 202304295107020000), Clinical medicine major project of Anhui University of Science and Technology medical special cultivation project (No. YZ2023H1A002), Open Fund Project of Anhui Provincial Occupational Health and Safety Engineering Research Center (No. AYZJSGCLK202202003), Natural Science Research Project of Anhui Educational Committee (No. KJ2020A0337), Kunshan Special Science and Technology Program for Social Development (No. KSZ2346), and Huainan Science and Technology Bureau Natural Science Project (No. 2022175).

Acknowledgments

The data analysis, article-editing and revising process, and kind guidance from Prof. Jianhua Wang, Cancer Institute, Fudan University Shanghai Cancer Center.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2025.1609616/full#supplementary-material

References

1. Dong Y, Tan H, Wang L, and Liu Z. Progranulin promoted the proliferation, metastasis, and suppressed apoptosis via JAK2-STAT3/4 signaling pathway in papillary thyroid carcinoma. Cancer Cell Int. (2023) 23:191. doi: 10.1186/s12935-023-03033-2

PubMed Abstract | Crossref Full Text | Google Scholar

2. Handra C M, Gurzu I L, Chirila M, and Ghita I. Silicosis: new challenges from an old inflammatory and fibrotic disease. Front Biosci (Landmark Ed). (2023) 28:96. doi: 10.31083/j.fbl2805096

PubMed Abstract | Crossref Full Text | Google Scholar

3. Lee K I, Su C C, Fang K M, Wu C C, Wu C T, Chen Y W, et al. Ultrafine silicon dioxide nanoparticles cause lung epithelial cells apoptosis via oxidative stress-activated PI3K/Akt-mediated mitochondria- and endoplasmic reticulum stress-dependent signaling pathways. Sci Rep. (2020) 10:9928. doi: 10.1038/s41598-020-66644-z

PubMed Abstract | Crossref Full Text | Google Scholar

4. Li Y, Xu H, Wang Y, Zhu Y, Xu K, Yang Z, et al. Epithelium-derived exosomes promote silica nanoparticles-induced pulmonary fibroblast activation and collagen deposition via modulating fibrotic signaling pathways and their epigenetic regulations. J Nanobiotechnology. (2024) 22:331. doi: 10.1186/s12951-024-02609-y

PubMed Abstract | Crossref Full Text | Google Scholar

5. Chen Y, Hao X, Li M, Tian Z, and Cheng M. UGRP1-modulated MARCO (+) alveolar macrophages contribute to age-related lung fibrosis. Immun Ageing. (2023) 20:14. doi: 10.1186/s12979-023-00338-8

PubMed Abstract | Crossref Full Text | Google Scholar

6. Parimon T, Yao C, Stripp B R, Noble P W, and Chen P. Alveolar epithelial type II cells as drivers of lung fibrosis in idiopathic pulmonary fibrosis. Int J Mol Sci. (2020) 21(7). doi: 10.3390/ijms21072269

PubMed Abstract | Crossref Full Text | Google Scholar

7. Purev E, Bahmed K, and Kosmider B. Alveolar organoids in lung disease modeling. Biomolecules. (2024) 14(1). doi: 10.3390/biom14010115

PubMed Abstract | Crossref Full Text | Google Scholar

8. Leiva-Juárez M M, Kolls J K, and Evans S E. Lung epithelial cells: therapeutically inducible effectors of antimicrobial defense. Mucosal Immunol. (2018) 11:21–34. doi: 10.1038/mi.2017.71

PubMed Abstract | Crossref Full Text | Google Scholar

9. Chuquimia O D, Petursdottir D H, Periolo N, and Fernández C. Alveolar epithelial cells are critical in protection of the respiratory tract by secretion of factors able to modulate the activity of pulmonary macrophages and directly control bacterial growth. Infect Immun. (2013) 81:381–9. doi: 10.1128/IAI.00950-12

PubMed Abstract | Crossref Full Text | Google Scholar

10. Wang Y, Tang Z, Huang H, Li J, Wang Z, Yu Y, et al. Pulmonary alveolar type I cell population consists of two distinct subtypes that differ in cell fate. Proc Natl Acad Sci U S A. (2018) 115:2407–12. doi: 10.1073/pnas.1719474115

PubMed Abstract | Crossref Full Text | Google Scholar

11. Cai Y, Yoneda M, Tomita T, Kurotani R, Okamoto M, Kido T, et al. Transgenically-expressed secretoglobin 3A2 accelerates resolution of bleomycin-induced pulmonary fibrosis in mice. BMC Pulm Med. (2015) 15:72. doi: 10.1186/s12890-015-0065-4

PubMed Abstract | Crossref Full Text | Google Scholar

12. Cai Y and Kimura S. Secretoglobin 3A2 exhibits anti-fibrotic activity in bleomycin-induced pulmonary fibrosis model mice. PloS One. (2015) 10:e0142497. doi: 10.1371/journal.pone.0142497

PubMed Abstract | Crossref Full Text | Google Scholar

13. Yang X, Liu Z, Zhou J, Guo J, Han T, Liu Y, et al. SPP1 promotes the polarization of M2 macrophages through the Jak2/Stat3 signaling pathway and accelerates the progression of idiopathic pulmonary fibrosis. Int J Mol Med. (2024) 54(4). doi: 10.3892/ijmm.2024.5413

PubMed Abstract | Crossref Full Text | Google Scholar

14. Xu C, Sun L, Jiang C, Zhou H, Gu L, Liu Y, et al. SPP1, analyzed by bioinformatics methods, promotes the metastasis in colorectal cancer by activating EMT pathway. BioMed Pharmacother. (2017) 91:1167–77. doi: 10.1016/j.biopha.2017.05.056

PubMed Abstract | Crossref Full Text | Google Scholar

15. Chen M, Wang J, Yuan M, Long M, Sun Y, Wang S, et al. AT2 cell-derived IgA trapped by the extracellular matrix in silica-induced pulmonary fibrosis. Int Immunopharmacol. (2023) 122:110545. doi: 10.1016/j.intimp.2023.110545

PubMed Abstract | Crossref Full Text | Google Scholar

16. Chung MI and Hogan BLM. Ager-creER (T2): A new genetic tool for studying lung alveolar development, homeostasis, and repair. Am J Respir Cell Mol Biol. (2018) 59:706–12. doi: 10.1165/rcmb.2018-0125OC

PubMed Abstract | Crossref Full Text | Google Scholar

17. Battaglia R A, Delic S, Herrmann H, and Snider N T.V. Vimentin on the move: new developments in cell migration. F1000Res. (2018) 7. doi: 10.12688/f1000research

PubMed Abstract | Crossref Full Text | Google Scholar

18. Chen Y, Zhou Y, Bai Y, Jia K, Zhang H, Chen Q, et al. Macrophage-derived CTSS drives the age-dependent disruption of the blood-CSF barrier. Neuron. (2025) 113(7):1082–97.e8. doi: 10.1016/j.neuron.2025.01.023

PubMed Abstract | Crossref Full Text | Google Scholar

19. Velez-Delvalle C, Marsch-Moreno M, Castro-Muñozledo F, Galván-Mendoza I J, and Kuri-Harcuch W. Epithelial cell migration requires the interaction between the vimentin and keratin intermediate filaments. Sci Rep. (2016) 6:24389. doi: 10.1038/srep24389

PubMed Abstract | Crossref Full Text | Google Scholar

20. Feng X, Ji Y, Zhang C, Jin T, Li J, Guo J, et al. CCL6 promotes M2 polarization and inhibits macrophage autophagy by activating PI3-kinase/Akt signalling pathway during skin wound healing. Exp Dermatol. (2023) 32:403–12. doi: 10.1111/exd.14718

PubMed Abstract | Crossref Full Text | Google Scholar

21. Wang S, Song R, Wang Z, Jing Z, Wang S, Ma J, et al. S100A8/A9 in inflammation. Front Immunol. (2018) 9:1298. doi: 10.3389/fimmu.2018.01298

PubMed Abstract | Crossref Full Text | Google Scholar

22. Nakamura Y, Aihara R, Iwata H, Kuwayama T, and Shirasuna K. IL1B triggers inflammatory cytokine production in bovine oviduct epithelial cells and induces neutrophil accumulation via CCL2. Am J Reprod Immunol. (2021) 85:e13365. doi: 10.1111/aji.13365

PubMed Abstract | Crossref Full Text | Google Scholar

23. Hsu L C, Li C J, Lin L T, Pan L F, Wen Z H, Sheu J J, et al. Multi-omics reveals the role of osteopontin/secreted phosphoprotein 1 in regulating ovarian aging. J Pers Med. (2024) 14(1). doi: 10.3390/jpm14010078

PubMed Abstract | Crossref Full Text | Google Scholar

24. Nickel S, Vlaic S, Christ M, Schubert K, Henschler R, Tautenhahn F, et al. Mesenchymal stromal cells mitigate liver damage after extended resection in the pig by modulating thrombospondin-1/TGF-β. NPJ Regener Med. (2021) 6:84. doi: 10.1038/s41536-021-00194-4

PubMed Abstract | Crossref Full Text | Google Scholar

25. Jayachandran A, Anaka M, Prithviraj P, Hudson C, McKeown S J, Lo P H, et al. Thrombospondin 1 promotes an aggressive phenotype through epithelial-to-mesenchymal transition in human melanoma. Oncotarget. (2014) 5:5782–97. doi: 10.18632/oncotarget.v5i14

PubMed Abstract | Crossref Full Text | Google Scholar

26. Choi J, Park J E, Tsakogeorga G, Yanagita M, Koo B K, Han N, et al. Inflammatory signals induce AT2 cell-derived damage-associated transient progenitors that mediate alveolar regeneration. Cell Stem Cell. (2020) 27:366–82.e7. doi: 10.1016/j.stem.2020.06.020

PubMed Abstract | Crossref Full Text | Google Scholar

27. Kobayashi Y, Tata A, Konkimalla A, Katsura H, Lee R F, Ou J, et al. Persistence of a regeneration-associated, transitional alveolar epithelial cell state in pulmonary fibrosis. Nat Cell Biol. (2020) 22:934–46. doi: 10.1038/s41556-020-0542-8

PubMed Abstract | Crossref Full Text | Google Scholar

28. Kadur Lakshminarasimha Murthy P, Sontake V, Tata A, Kobayashi Y, Macadlo L, Okuda K, et al. Human distal lung maps and lineage hierarchies reveal a bipotent progenitor. Nature. (2022) 604:111–9. doi: 10.1038/s41586-022-04541-3

PubMed Abstract | Crossref Full Text | Google Scholar