Monocytes (Mos)/macrophages (Mϕs) orchestrate biological processes critical for efficient skin wound healing. However, current understanding of skin wound Mo/Mϕ heterogeneity is limited by traditional experimental approaches such as flow cytometry and immunohistochemistry. Therefore, we sought to more fully explore Mo/Mϕ heterogeneity and associated state transitions during the course of excisional skin wound healing in mice using single-cell RNA sequencing. The live CD45+CD11b+Ly6G cells were isolated from skin wounds of C57BL/6 mice on days 3, 6, and 10 postinjury and captured using the 10x Genomics Chromium platform. A total of 2813 high-quality cells were embedded into a uniform manifold approximation and projection space, and eight clusters of distinctive cell populations were identified. Cluster dissimilarity and differentially expressed gene analysis categorized those clusters into three groups: early-stage/proinflammatory, late-stage/prohealing, and Ag-presenting phenotypes. Signature gene and Gene Ontology analysis of each cluster provided clues about the different functions of the Mo/Mϕ subsets, including inflammation, chemotaxis, biosynthesis, angiogenesis, proliferation, and cell death. Quantitative PCR assays validated characteristics of early- versus late-stage Mos/Mϕs inferred from our single-cell RNA sequencing dataset. Additionally, cell trajectory analysis by pseudotime and RNA velocity and adoptive transfer experiments indicated state transitions between early- and late-state Mos/Mϕs as healing progressed. Finally, we show that the chemokine Ccl7, which was a signature gene for early-stage Mos/Mϕs, preferentially induced the accumulation of proinflammatory Ly6C+F4/80lo/− Mos/Mϕs in mouse skin wounds. In summary, our data demonstrate the complexity of Mo/Mϕ phenotypes, their dynamic behavior, and diverse functions during normal skin wound healing.

Skin wound healing involves a series of sequential but overlapping phases, including inflammation, proliferation, and remodeling, that restores skin barrier integrity, function, and homeostasis (1, 2). Wound healing is a complex and dynamic process that requires collaboration and crosstalk between a number of different cell types (1). After injury, circulating leukocytes rapidly infiltrate injured sites, playing important roles in diverse biological processes (13). Among these cells, monocytes (Mos)/macrophages (Mϕs) play critical roles in all phases of wound healing (2, 3). First, Mos/Mϕs help to initiate and amplify the inflammatory response in part by secreting cytokines and other mediators that influence activities of wound cells, contributing to host defense and debridement of the injury site. During the proliferation and remodeling phases, Mo/Mϕ-derived factors promote resolution of inflammation, as well as proliferation and migration of fibroblasts, endothelial cells, and epithelial cells, resulting in angiogenesis/vessel formation, granulation formation, and collagen deposition (2, 46).

Mos/Mϕs exhibit a high degree of plasticity and perform diverse functions under different pathophysiological conditions (4, 6). A widely used classification scheme holds that Mϕs polarize into two major categories: proinflammatory M1 (classically activated) and anti-inflammatory M2 (alternatively activated). However, this M1/M2 terminology was established based on responses simplified in in vitro conditions: IFN-γ with LPS/TNF induces the M1 phenotype, and IL-4/IL-13 induces the M2 phenotype (7, 8). However, this scheme has been challenged by in vivo studies (6, 9, 10), including studies of skin wound Mϕs, which adopt a wide spectrum of phenotypes depending on stage of healing, wound environment, lineage, differentiation, and mitochondrial metabolism (6, 11). Indeed, Mϕs can coexpress both M1 and M2 markers at various time points of the wound healing process (1214). Thus, the attempt to fit Mos/Mϕs into M1 or M2 categories based on the expression of a small number of markers does not describe well the full spectrum of Mo/Mϕ heterogeneity during wound healing and leads to biased interpretation of results.

The exploration of Mo/Mϕ heterogeneity has traditionally been limited by experimental approaches that rely on measurement of a small number of markers at the single-cell level via flow cytometry or immunohistochemistry (6, 12, 13, 15, 16). Other studies have explored Mϕ phenotypes at the bulk level by quantitative PCR (qPCR) or more recently by RNA sequencing (RNA-seq), but these techniques cannot distinguish between different subsets present in the same sample (1720). The recent introduction of single-cell RNA-seq (scRNA-seq) has enabled more comprehensive and less biased investigations into cell heterogeneity in an array of pathophysiological processes (21, 22). Mo/Mϕ heterogeneity has been explored in a variety of conditions, including atherosclerosis, acute lung inflammation, spinal cord injury, and cancer (2326). During skin wound healing, scRNA-seq has revealed plasticity of wound fibroblasts, keratinocytes, hair follicle stem cells, and epidermal basal cells, as well as diverse Mϕ populations that exhibit varying mitochondrial metabolism status (2731). However, much remains to be learned about the heterogeneity of skin wound Mos/Mϕs and potential state transitions throughout the course of healing.

Therefore, our study aimed to explore Mo/Mϕ heterogeneity and associated state transitions during the course of skin wound healing. Using Ly6GCD11b+CD45+ cells sorted from skin wounds collected on days 3, 6, and 10 postinjury, we identified eight clusters of cells that could be categorized into three groups: early-stage Mos/Mϕs that are predominantly proinflammatory, late-stage Mos/Mϕs that exhibit prohealing characteristics, and APCs, all of which were implicated in biological processes associated with the healing process. Moreover, cell trajectory analysis indicated state transitions between Mos/Mϕs at different stages of healing, and such state transitions were confirmed by adoptive transfer experiments. Finally, we demonstrate that chemokine Ccl7, which was a signature gene for the early-stage Mos/Mϕs, but had not been studied during skin wound healing, preferentially induced the accumulation of Ly6C+F4/80lo/− Mos/Mϕs in wounds during the inflammatory phase of healing.

Male C57BL6/J (CD45.2) and B6.SJL-Ptprca Pepcb/BoyJ (CD45.1) mice were purchased from The Jackson Laboratory (Bar Harbor, ME) and housed in an environmentally controlled animal facility for at least 2 wk before experiments. Adult mice (CD45.2, age 10–11 wk) were subjected to excisional wounding of the dorsal skin with an 8-mm biopsy punch (two wounds per mouse) as previously described (32). Tegaderm (3M, 1626W) was used to cover the wounds until sacrifice. Skin wounds were collected on days 3, 6, and 10 postinjury using a 12-mm biopsy punch to include both wound and surrounding tissue. After harvesting, tissue was either snap-frozen in liquid nitrogen for protein measurement or dissociated into a single-cell suspension by mechanical and enzymatic digestion as described (32), for subsequent analysis. Briefly, wound tissue was placed in a petri dish filled with 1 ml of enzyme mixture containing collagenase I, collagenase XI, and hyaluronidase (5 mg/ml, Sigma-Aldrich, St. Louis, MO). Forceps were used to flatten the sample, the epidermis was carefully scraped from the wound using a Lance scalpel blade (Cincinnati Surgical, Cincinnati, OH), and then the dermis was minced and digested in enzyme mixture for 1 h at 37°C. Finally, cells were filtered through a 40-µm cell strainer and resuspended in cold 1× Dulbecco’s PBS (DPBS).

All animal procedures were approved by the Animal Care and Use Committee of the University of Illinois at Chicago.

Skin cells were pooled from four wounds (two wounds from each of two mice) collected on days 3, 6, and 10 postinjury as described above. After washing, cells were resuspended in 1× DPBS and incubated with a Zombie Violet fixable viability kit (BioLegend, San Diego, CA) for cell viability analysis and then followed by Fc blocking with TruStain FcX Ab (BioLegend, anti-mouse CD16/32, clone S17011E). Next, surface markers were labeled with optimized dilutions of anti–Ly6G-PE-CF594 (BD Biosciences, San Jose, CA, clone 1A8), CD11b-APC/Fire 750 (BioLegend, clone M1/70), and CD45-FITC (BioLegend, clone S18009F). Mo/Mϕ cell populations were sorted as ZombieLy6GCD11b+CD45+ cells using a BD FACSAria III sorter (BD Biosciences).

After sorting, fresh skin Mos/Mϕs (ZombieLy6GCD11b+CD45+) were collected at a final concentration of ∼200 cells/µl in cold 1× DPBS with 10% FBS. Cell viability was evaluated in each sample immediately prior to capture (day 3, 88.7%; day 6, 71.3%; and day 10, 89.3%). Next, cells were loaded on the 10x Chromium Next GEM (gel beads in emulsion) Chip G with gel beads and master mix for cell capture and GEM generation targeting 2000 cells recovery per sample. Following GEM reverse transcription cleanup, cDNA amplification and 3′ gene expression library construction was performed according to the manufacturer’s instructions (10x Genomics, San Francisco, CA). Libraries were sequenced on HiSeq sequencing systems (Illumina, San Diego, CA) with paired-end reads aiming for 100,000 reads/cell.

Raw scRNA-seq data were demultiplexed and aligned to reference mouse genome GRCm38 (mm10) into single cells using Cell Ranger software (10x Genomics). Further data filtering and analysis were conducted using R package Seurat (version 3.1.5) (33). Quality control included omitting cells with >10% mitochondrial gene expression, unique features/gene expression <500 and >10,000, or <1000 unique molecular identifier counts; 2831 cells passed these quality control filters and were included in downstream analysis (Supplemental Fig. 1A–C). Data were then normalized, and the top 6000 most variably expressed genes were used for unsupervised principal component analysis (Supplemental Fig. 1D). The top 50 principal components were used for clustering and dimensional reduction via uniform manifold approximation and projection (UMAP) analysis (Supplemental Fig. 1E). Differentially expressed (DE) gene analysis across all clusters was performed using the Wilcoxon rank sum test to identify marker genes for each cluster. For cluster composition comparisons, we compared the number of cells in each clustering originating from each sample. For each sample and cluster combination, a 2 × 2 contingency table was constructed based on 1) the number of cells in the sample and the cluster, 2) other cells in the cluster, 3) other cells in the sample, and 4) all other cells in the dataset, and association statistics (odds ratio and p value) were computed using a Fisher’s exact test. To quantify high-level differences in overall gene expression levels between clusters, analysis of similarities (ANOSIM) R statistics were computed in the vegan package in R (https://cran.r-project.org/web/packages/vegan/vegan.pdf).

Gene Ontology (GO) enrichment analysis of the genes from DE analysis from each cluster was performed using the Database for Annotation, Visualization and Integrated Discovery (DAVID) 2021 (false discover rate < 5%). Functional annotation clustering was conducted using medium classification stringency on all GO terms in each cluster, and similar terms were summarized as representative GO terms. The overall enrichment score was calculated based on the Expression Analysis Systematic Explorer (EASE) score of each group member (34).

RNA velocity for cell transition and trajectory analysis was conducted using the scVelo package using the dynamical model (35). Partition-based graph abstraction (PAGA) was established based on the results of RNA velocity results (36). We also performed pseudotime analysis using R package Monocle 3 (37). Two “roots” were used as starting points based on the cells with low velocity to order the cells in pseudotime. In addition, we used Ly6c2hi cells (expression > median Ly6c2 expression) in cluster 3 from the day 3 sample (early stage) as root cells and constructed cell trajectory between early- and late-stage clusters.

To verify scRNA-seq results, expression of selected genes was evaluated in skin Ly6C+ versus Ly6C cells. First, skin cells isolated from skin wounds on days 3 and 10 postinjury were processed similarly as above with a Zombie Green fixable viability kit, TruStain FcX (anti-mouse CD16/32), anti–Ly6G-PE/Dazzle 594 (BioLegend, clone 1A8), CD11b-APC/Fire 750 (BioLegend, clone M1/70), CD45-allophycocyanin (BioLegend, clone I3/2.3), and Ly6C-BV421 (BD Biosciences, clone AL-21) Abs. Then, target cells were sorted as ZombieLy6GCD11b+CD45+Ly6Chi/− on a MoFlo Astrios cell sorter (Beckman Coulter, Brea, CA).

Next, freshly sorted Ly6C+ versus Ly6C cells from day 3 and 10 skin wounds were washed with cold 1× DPBS and resuspended in ∼5–10 μl of 1× DPBS. Cells were immediately lysed, and cDNA was prepared using SYBR Green Fast Advanced Cells-to-CT kit (Invitrogen, Waltham, MA) following the manufacturer’s instructions. Gene expression of housekeeping and selective signature genes from each cluster was assessed by real-time RT-PCR using a Power SYBR Green PCR system (Applied Biosystems, Waltham, MA). Primers were designed and obtained from Integrated DNA Technologies (Coralville, IA). Sequences of all primers are listed in Table I. We assessed expression of a number of commonly used housekeeping genes in our samples and chose the three most stable genes, that is, Ywhaz, Rpl27, and Rpl13a, for normalization of the data. The housekeeping gene expression level was calculated based on the geometric mean of three genes, that is, Ywhaz, Rpl27, and Rpl13a, to avoid variations between samples. Relative expression was calculated using the ΔΔCt method after normalizing for expression of the housekeeping gene.

To test the trajectory inferences from the RNA velocity analysis of the scRNA-seq data, an adoptive transfer experiment was performed. Bone marrow (BM) mononuclear cells from CD45.1 donor mice were first enriched by density gradient centrifugation. Cells were then prepared for FACS sorting as described above with a Zombie Aqua fixable viability kit, TruStain FcX (anti-mouse CD16/32), anti–Ly6G-BV510 (BioLegend, clone 1A8), CD11b-APC/Fire 750, CD115-PE (BioLegend, clone AFS98), c-Kit-PerCP/Cy5.5 (BioLegend, clone 2B8), and Ly6C-BV421. Naive BM Ly6Chi Mos were sorted as ZombieLy6Gc-KitCD11b+CD115+Ly6Chi on a MoFlo Astrios cell sorter (Beckman Coulter).

Donor cells were washed and resuspended in cold 1× DPBS. Next, 2 million donor cells (CD45.1) were delivered to skin wounds of recipient mice (CD45.2, n = 3 per group) by intradermal injections on day 2 or day 5 after wounding. At 18 h after cell transfer (on day 3 and day 6 postinjury), skin wounds from recipient mice were collected and single-cell suspensions were prepared as described above for flow cytometry analysis. Skin wound cells from recipient mice were incubated with a Zombie Aqua fixable viability kit, followed by TruStain FcX (anti-mouse CD16/32) for Fc blocking, and then CD11b-APC/Fire 750, F4/80-PE (BioLegend, clone 1A8), Ly6C-PerCP/Cy5.5 (BioLegend, clone 2B8), and CD45.1-BV421 (BioLegend, clone A20) Abs. Adoptively transferred donor cells (CD45.1 mice) were identified as ZombieCD11b+CD45.1+ cells. All samples were analyzed on an LSRFortessa with a high-throughput sampler (BD Biosciences). Data were analyzed using FlowJo (FlowJo, Ashland, OR).

To verify one of the key DE genes at the protein level, noninjured skin (n = 4) and wounds harvested on days 3 and 10 postinjury (n = 10 per group) were mechanically homogenized with cold 1× DPBS (10 µl/mg tissue) with 1% protease inhibitor mixture (Sigma-Aldrich). After centrifuging, supernatants from each sample were collected for Ccl7 protein measurement by a custom BioLegend LEGENDplex kit following the manufacturer’s instructions. Sample acquisition was performed using a CytoFLEX S (Beckman Coulter) cytometer. Data were analyzed using the cloud-based LEGENDplex data analysis software suite (BioLegend).

To assess the function of Ccl7 in wounds, 30 or 120 ng of recombinant mouse Ccl7 protein (carrier-free, BioLegend) was injected intradermally around the periphery of each wound at four evenly distributed sites daily starting immediately after wounding until day 2 postinjury (n = 6 mice per group, three independent experiments). Control mice received injections of an equal volume of sterile 1× DPBS (n = 6). Twenty-four hours after the last injection (day 3 postinjury), all wounds were harvested and cells were dissociated from wounds and analyzed by flow cytometry as described above with a Zombie Aqua fixable viability kit, TruStain FcX (anti-mouse CD16/32), Ly6G-BV605 (clone 1A8), CD11b-APC/Fire 750, F4/80-PE, and Ly6C-PerCP/Cy5.5 Abs (BioLegend or BD Biosciences). Data were analyzed using FlowJo (FlowJo).

A loss-of-function experiment was also performed using 20 µg of mouse Ccl7 neutralizing Ab (NAb; catalog no. AF-456-NA, R&D Systems, Minneapolis, MN) injected intradermally around the periphery of skin wounds similarly as for the recombinant protein experiment, immediately after wounding (day 0) and on day 2 postinjury (n = 5 mice per group, three independent experiments). Control mice received injections of an equal amount of IgG protein (catalog no. AB-108-C, R&D Systems). Two different doses of Ccl7 blocking Ab were tested in a preliminary dose-response experiment. A high dose of Ccl7 NAb (20 µg/wound) appeared to reduce Ly6C+F4/80lo/− Mo/Mϕ accumulation in wounds whereas a low dose (5 µg/wound) showed no effect; therefore, the higher dose was chosen for the full experiment. All wounds were collected on day 3 (24 h after the last injection) and cells were processed for flow cytometry analysis as for the recombinant protein experiments.

To assess the effects of Ccl7 on wound closure, digital images were taken of the wound surface. Wound size was measured on the days when injections were performed and the day of wound collection using Fiji ImageJ and expressed as a percentage of the original wound area, or (day x area/day 0 area) ×100.

Data are expressed as mean ± SEM. All data were tested for normality. Statistical significance of differences was evaluated by a one- or two-way ANOVA with post hoc tests. A p value of <0.05 was considered statistically significant.

scRNA-seq data in this study are available in Gene Expression Omnibus (www.ncbi.nlm.nih.gov/geo) under accession number GSE203244. All other relevant data are available from the corresponding author.

To better understand the heterogeneity of wound Mos/Mϕs that accumulate in skin wounds during the course of healing, we generated scRNA-seq data using sorted the live Ly6GCD11b+CD45+ Mos/Mϕs from wounds at an inflammatory stage (day 3 postinjury), a proliferative stage (day 6), and a remodeling stage (day 10); complete re-epithelialization typically occurs by day 10 in this model (16). After quality control filtering, 2831 cells pooled from all three time points (day 3, 894; day 6, 585; day 10, 1334 cells) were used for downstream analysis (Supplemental Fig. 1). Unsupervised dimensionality reduction by UMAP and graph-based clustering analysis identified eight Mo/Mϕ clusters (Fig. 1A). In this UMAP plot, clusters 3, 4, and 6 appeared to form a group and clusters 0, 2, and 5 seemed to form another group, whereas clusters 1 and 7 seemed to be distinct from other clusters. These impressions were confirmed by ANOSIM analysis, for which cluster 7 showed the highest dissimilarity with other clusters. Additionally, clusters 0, 2, and 5 showed low dissimilarity to each other and were hierarchically grouped together whereas clusters 1, 3, 4, and 6 formed another hierarchical group (Fig. 1B).

FIGURE 1.

scRNA-seq analysis reveals heterogeneity of Mos/Mϕs in skin wounds in mice. The live CD11b+Ly6G cells were isolated from excisional skin wounds of C57BL/6 mice on days 3, 6, and 10 postinjury and captured using the 10x Genomics Chromium platform. After sequencing, demultiplexing, and alignment to the reference mouse genome, a total of 2813 cells passed quality control, and the top 6000 most variably expressed genes were used for unsupervised principal component analysis. (A) The top 50 principal components were used for unsupervised dimensional reduction by UMAP, and graph-based clustering analysis identified eight Mo/Mϕ clusters. (B) Cluster dissimilarity analysis by ANOSIM indicated hierarchical grouping of Mo/Mϕ clusters. (C) Heatmap displaying the top 10 differentially expressed genes in each cluster. (D) Feature maps showing expression patterns of selected signature genes representing each cluster. (E) GO clusters representing the top five most enriched GO term groups in each cluster. The overall enrichment score was calculated based on the Expression Analysis Systematic Explorer (EASE) score of each group member.

FIGURE 1.

scRNA-seq analysis reveals heterogeneity of Mos/Mϕs in skin wounds in mice. The live CD11b+Ly6G cells were isolated from excisional skin wounds of C57BL/6 mice on days 3, 6, and 10 postinjury and captured using the 10x Genomics Chromium platform. After sequencing, demultiplexing, and alignment to the reference mouse genome, a total of 2813 cells passed quality control, and the top 6000 most variably expressed genes were used for unsupervised principal component analysis. (A) The top 50 principal components were used for unsupervised dimensional reduction by UMAP, and graph-based clustering analysis identified eight Mo/Mϕ clusters. (B) Cluster dissimilarity analysis by ANOSIM indicated hierarchical grouping of Mo/Mϕ clusters. (C) Heatmap displaying the top 10 differentially expressed genes in each cluster. (D) Feature maps showing expression patterns of selected signature genes representing each cluster. (E) GO clusters representing the top five most enriched GO term groups in each cluster. The overall enrichment score was calculated based on the Expression Analysis Systematic Explorer (EASE) score of each group member.

Close modal

Next, we began to explore the characteristics of each cluster using DE gene analysis. Expression of the top 10 DE genes in each cluster is shown in the heatmap in (Fig. 1C. In addition, feature plots of signature genes for each cluster and violin plots of the top 10 genes for each cluster are presented in (Fig. 1D and Supplemental Fig. 2, respectively. As shown in those plots, genes involved in inflammation and cell differentiation (Pf4/Cxcl4) (38, 39), lipid metabolism (Apoe) (40), lysosomal function (Ctsb, Ctsd, and Lgmn), and the complement system (C1qa and C1qc) were highly expressed in cluster 0 (38, 4144). Cluster 2 was enriched with genes involved in matrix remodeling and angiogenesis (Pdpn and Mmp12) (45, 46) and metal ion homeostasis (Mt1 and Mt2), as well as the prohealing/anti-inflammatory marker gene Arg-1 (47). Cluster 5 showed upregulation of complement system–associated genes (C1qa, C1qb, and C1qc) and a widely used M2 Mϕ marker, Cd163 (48). Clusters 3 and 4 shared similar features, as both were enriched with genes coding for proinflammatory cytokines and protein markers, particularly related to IFN signaling pathways, including Ccl7, Ly6c2, Irf7, and Ly6e for cluster 3 and Cxcl10, Cxcl3, Thbs1, and Vcan for cluster 4 (4955). Common Mo/Mϕ marker genes (Sell and Lyz2) were highly expressed in cluster 6, together with genes encoding both proinflammatory (Ifitm6) and anti-inflammatory proteins (Hp) as well as anti-oxidative regulators (Mgst1 and Msrb1) (9, 10, 5659). In cluster 1, Ag processing/presentation–associated genes were highly expressed (H2-Eb1, H2-Ab1, H2-Aa, and Cd209) (60). Finally, cluster 7 was the smallest and most unique subpopulation of cells and was dominated by expression of chemokines and their receptors (Ccr7, Ccl22, Ccl5, and Ccl17) as well as mature dendritic cell (DC) markers (Fscn1 and Cacnb3) (10, 61).

To further elucidate the potential functional heterogeneity of the Mo/Mϕ clusters, we performed GO analysis based on the gene signatures for each cluster. GO terms that shared similar biological meanings were grouped together in each GO cluster and the top five groups were presented and ranked by enrichment score in (Fig. 1E. As expected, there were several overlapping biological processes between some clusters. GO terms of cell motility and angiogenesis, which are two well-defined functions of wound Mos/Mϕs (28), were enriched in cluster 0 as well as cell death and endocytosis. Similar GO terms were enriched in clusters 2 and 5 (cell motility, angiogenesis, and/or cell death), supporting the idea that clusters 0, 2, and 5 shared healing-associated features with each other. Clusters 3 and 4 were highly enriched for GO terms associated with response to cytokines and immune responses, supporting the possibility that those two clusters are proinflammatory and help drive inflammation in wounds. Moreover, cluster 6 was appreciably enriched in GO terms associated with protein biosynthesis together with phagocytosis and endocytosis, suggesting that highly active cytokine-producing and phagocytic Mos/Mϕs populated this cluster. Interestingly, protein synthesis and Ag processing and presentation–related terms appeared among the top GO term groups in cluster 1, indicating that those Mos/Mϕs worked mainly as APCs in the healing process. In addition to cluster 1, Ag processing and presentation was also identified in cluster 7, which had upregulation of mature DC marker genes Fscn1 and Cacnb3. Therefore, cluster 7 contained a small population of APCs that may include Mϕs and/or DCs.

Table I.

Sequences of primers

GeneForward (5′→3′)Reverse (5′→3′)
Ywhaz TTGAGCAGAAGACGGAAGGT GAAGCATTGGGGATCAAGAA 
Rpl27 GCGATCCAAGATCAAGTCCTTTG TCAAAGCTGGGTCCCTGAACAC 
Rpl13a AGCCTACCAGAAAGTTTGCTTAC GCTTCTTCTTCCGATAGTGCATC 
Ccl7 CAGAAGGATCACCAGTAGTCGG ATAGCCTCCTCGACCCACTTCT 
Ly6e CGGGCTTTGGGAATGTCAAC GTGGGATACTGGCACGAAGT 
Ms4a4c CCAGCGAGGTGAGTAGCAGTA CTTAAGACTGAGCAAGGGGCA 
Tmsb10 CCACGAGTTGTAAGAAAATGGCAG GGGAAATCTTCCTAGGCTTTTAGGAG 
Arg1 AAGCCAGGGACTGACTACCTTAAA TGATGCCCCAGATGGTTTTC 
Pdpn GCCAACT-TGTAACTGCACTTGAT GAAGGCACAGTGGATGCTTAG 
Bnip3 GCTCCAAGAGTTCTCACTGTGAC GTTTTTCTCGCCAAAGCTGTGGC 
Cd163 CAGGACAGCCAGCTAAATGGA ATGAAGATGCAGCTGTGGTCAT 
Cxcl3 TGAGACCATCCAGAGCTTGACG CCTTGGGGGTTGAGGCAAACTT 
Isg15 CTGAAGAAGCAGATTGCCCAGAAG CGCTGCAGTTCTGTACCACTAGC 
Ms4a7 ATTGGAGTGCCTAGTGAAGCCACA GGTTGAAGTAGGAAGCTGACACCA 
C1qa GGATGGGGCTCCAGGAAATC CTGATATTGCCTGGATTGCC 
GeneForward (5′→3′)Reverse (5′→3′)
Ywhaz TTGAGCAGAAGACGGAAGGT GAAGCATTGGGGATCAAGAA 
Rpl27 GCGATCCAAGATCAAGTCCTTTG TCAAAGCTGGGTCCCTGAACAC 
Rpl13a AGCCTACCAGAAAGTTTGCTTAC GCTTCTTCTTCCGATAGTGCATC 
Ccl7 CAGAAGGATCACCAGTAGTCGG ATAGCCTCCTCGACCCACTTCT 
Ly6e CGGGCTTTGGGAATGTCAAC GTGGGATACTGGCACGAAGT 
Ms4a4c CCAGCGAGGTGAGTAGCAGTA CTTAAGACTGAGCAAGGGGCA 
Tmsb10 CCACGAGTTGTAAGAAAATGGCAG GGGAAATCTTCCTAGGCTTTTAGGAG 
Arg1 AAGCCAGGGACTGACTACCTTAAA TGATGCCCCAGATGGTTTTC 
Pdpn GCCAACT-TGTAACTGCACTTGAT GAAGGCACAGTGGATGCTTAG 
Bnip3 GCTCCAAGAGTTCTCACTGTGAC GTTTTTCTCGCCAAAGCTGTGGC 
Cd163 CAGGACAGCCAGCTAAATGGA ATGAAGATGCAGCTGTGGTCAT 
Cxcl3 TGAGACCATCCAGAGCTTGACG CCTTGGGGGTTGAGGCAAACTT 
Isg15 CTGAAGAAGCAGATTGCCCAGAAG CGCTGCAGTTCTGTACCACTAGC 
Ms4a7 ATTGGAGTGCCTAGTGAAGCCACA GGTTGAAGTAGGAAGCTGACACCA 
C1qa GGATGGGGCTCCAGGAAATC CTGATATTGCCTGGATTGCC 

Forward and reverse primers given in the table were used for real-time RT-PCR using the Power SYBR Green PCR system.

Taken together, our initial analysis of live Ly6GCD11b+CD45+ Mo/Mϕ gene expression on the single-cell level suggest that these cells comprise a heterogeneous population of cells with different specific functions in mouse skin wounds.

To investigate the impact of time after wounding on changes in Mo/Mϕ phenotypes during wound healing, we further analyzed how cells from day 3, 6, and 10 samples were distributed among clusters. As presented in (Fig. 2A, cells from day 3 wounds were predominantly found in clusters 3, 4, and 6 along with cluster 1. In contrast, day 6 and 10 cells were distributed in clusters 0, 2, and 5 as well as cluster 1. Next, the log2(odds ratio) was calculated to quantify the contribution of cells from each time point to cluster composition, avoiding the influence of total cell number that passed quality control in each sample. As shown in (Fig. 2B, clusters 3, 4, and 6 were mainly populated by cells from day 3 wounds whereas clusters 0, 2, and 5 were populated by cells from day 6 and 10 wounds. The two clusters with Ag-presenting characteristics, that is, cluster 1 and 7, had more even contributions from all time points than the other clusters. Together with cluster signature gene expression, we separated the eight clusters into three major groups: early stage (clusters 3, 4, and 6, proinflammatory), late stage (clusters 0, 2, and 5, prohealing), and APCs (clusters 1 and 7). Signature genes of each cluster were then displayed by average expression levels and percentage of cells in the cluster that express the gene (Fig. 2C). Again, genes associated with inflammation and oxidative stress were highly expressed by cells in early-stage clusters whereas complement system and tissue repair–associated genes dominated the late-stage clusters.

FIGURE 2.

Skin wound Mos/Mϕs cluster into distinct phenotypes in a healing progression-dependent manner. (A) Feature plots showing distribution of cells from day 3, 6, and 10 wounds over the eight clusters identified. (B) Heatmap of cell contribution from each time point to cluster composition based on log2(odds ratio). Note that day 3 cells are overrepresented in clusters 3, 4, and 6 and day 6 and 10 cells are overrepresented in clusters 0, 2, and 5. Clusters 1 and 7 show more even representation over all time points. (C) Dot plot showing expression of selected signature genes for each cluster arranged by early-stage, late-stage, and APC cluster groups. Size of dots indicates percentage of cells in each cell cluster expressing the marker gene; color intensity represents averaged scaled expression levels. Expression levels of selected signature genes for early- (D) and late-stage (E) clusters were measured by real-time qPCR in FACS-sorted skin Ly6C+ versus Ly6C Mos/Mϕs from wounds harvested on day 3 and day 10 postinjury. Data are mean ± SEM. *p < 0.05 between four groups by two-way ANOVA.

FIGURE 2.

Skin wound Mos/Mϕs cluster into distinct phenotypes in a healing progression-dependent manner. (A) Feature plots showing distribution of cells from day 3, 6, and 10 wounds over the eight clusters identified. (B) Heatmap of cell contribution from each time point to cluster composition based on log2(odds ratio). Note that day 3 cells are overrepresented in clusters 3, 4, and 6 and day 6 and 10 cells are overrepresented in clusters 0, 2, and 5. Clusters 1 and 7 show more even representation over all time points. (C) Dot plot showing expression of selected signature genes for each cluster arranged by early-stage, late-stage, and APC cluster groups. Size of dots indicates percentage of cells in each cell cluster expressing the marker gene; color intensity represents averaged scaled expression levels. Expression levels of selected signature genes for early- (D) and late-stage (E) clusters were measured by real-time qPCR in FACS-sorted skin Ly6C+ versus Ly6C Mos/Mϕs from wounds harvested on day 3 and day 10 postinjury. Data are mean ± SEM. *p < 0.05 between four groups by two-way ANOVA.

Close modal

To validate differences in gene expression between clusters observed in the scRNA-seq data, we measured mRNA expression of selected signature genes from early-stage and late-stage clusters by real-time qPCR. Because Ly6c2 was highly expressed in all three early-stage clusters (clusters 3, 4, and 6; Supplemental Fig. 2A) and we and others have shown that Ly6C+ Mos/Mϕs predominate at early stages of skin wound healing (12, 16, 62), we compared the expression of early- and late-stage cluster signature genes between live Ly6GCD11b+CD45+Ly6C+ versus live Ly6GCD11b+CD45+Ly6C Mos/Mϕs sorted from day 3 and day 10 wounds. Consistent with our scRNA-seq data, Ccl7, Ly6e, Ms4a4c, and Tmbs10, which are four signature genes of early-stage clusters 3, 4, and 6, showed significantly higher expression levels in Ly6C+ than Ly6C Mos/Mϕs regardless of the time point postwounding (Fig. 2D). In addition, Arg1, Pdpn, Bnip3, and Cd163, which are signature genes of late-stage clusters 0, 2, and 5, were significantly upregulated in Ly6C Mos/Mϕs compared with Ly6C+ cells (Fig. 2E). However, we also found genes, including Cxcl3, Isg15, Ms4a7, and C1qa, that although exhibiting a similar nonsignificant trend between time points as the scRNA-seq data, did not show clear differences between Ly6C+ and Ly6C cells (Supplemental Fig. 2B, 2C).

We next assessed expression of commonly used Mϕ markers in each cluster. As expected, the common Mϕ markers Itgam and Cd68 were highly expressed by all clusters. Adgre1 and Cd14 expression were also highly expressed by most clusters except clusters 1 and 7, where these were expressed at lower levels. Ly6c2, a marker of inflammatory Mos/Mϕs, was preferentially expressed in early-stage clusters 3, 4, and 6 (Fig. 3A), which agrees with our previous findings of increased accumulation of proinflammatory Ly6C+ Mos/Mϕs in skin wounds at the early stage of healing (16). When assessing Mo/Mϕ phenotypes, the most widely used classification scheme divides cells into classically activated M1 and alternatively activated M2 phenotypes (63), although this scheme is best suited for in vitro classification. Hence, we assessed the expression of widely used markers of M1 and M2 Mϕs in the skin wound scRNA-seq data. Interestingly, when comparing the expression of M1 versus M2 markers, there was noticeable overlap of cells that highly expressed both categories of genes, such as the M1 marker Il1β and the M2 marker Arg1, as both were coexpressed in several clusters. Similarly, coexpression of Tnf and Mrc1, as well as Ccl2 and Chil3, were also identified (Fig. 3B, 3C). Interestingly, Chil3 is commonly used as an M2 marker but was primarily expressed in early-stage clusters that demonstrated a predominantly proinflammatory phenotype. Therefore, M1/M2 markers are not good differentiators of clusters of wound Mos/Mϕs in mice, and a more multifaceted system is needed to classify Mos/Mϕs in this in vivo context.

FIGURE 3.

Skin wound Mos/Mϕs express common Mo/Mϕ markers whereas M1/M2 markers do not differentiate clusters of wound Mos/Mϕs in mice. (AC) Feature plots showing expression patterns of (A) common markers for Mos/Mϕs (Itgam, CD68, Adgre1, Cd14, and Ly6c2) and widely used markers for (B) M1 (IL1β, Tnf, Nos2, Ccl2, Cd86) and (C) M2 (Arg1, Mrc1, Cd163, Igf1, and Chil3) phenotypes. Note that common Mo/Mϕ markers were expressed as expected but that many M1 and M2 markers were distributed across several clusters and most were coexpressed with markers of the antagonistic phenotype.

FIGURE 3.

Skin wound Mos/Mϕs express common Mo/Mϕ markers whereas M1/M2 markers do not differentiate clusters of wound Mos/Mϕs in mice. (AC) Feature plots showing expression patterns of (A) common markers for Mos/Mϕs (Itgam, CD68, Adgre1, Cd14, and Ly6c2) and widely used markers for (B) M1 (IL1β, Tnf, Nos2, Ccl2, Cd86) and (C) M2 (Arg1, Mrc1, Cd163, Igf1, and Chil3) phenotypes. Note that common Mo/Mϕ markers were expressed as expected but that many M1 and M2 markers were distributed across several clusters and most were coexpressed with markers of the antagonistic phenotype.

Close modal

As APCs are a heterogeneous cell population in skin (64), we further explored the composition of APC clusters 1 and 7. We first assessed expression of common Mϕ and DC markers in these clusters. Common DC markers Cd11c/Itgax, Cd40, Cd80, and Cd86 were expressed by a large proportion of cells in clusters 1 and 7; common marker genes of Mϕs were also expressed, but at lower levels (Supplemental Fig. 2D). Interestingly, a widely used DC marker Cd207 was expressed by very few cells in either cluster. Next, we examined expression of genes used as markers for subsets of DCs, including conventional type 1 DCs (cDC1s), conventional type 2 DCs (cDC2s), Mo-derived DCs (MoDCs), and mature DCs enriched in immunoregulatory molecules (mregDCs) (6567). Although clusters 1 and 7 expressed varying levels of common Mϕ markers (Cd68, F4/80/Adgre1, Lyz2, Mac-2/Lgals3, and Cd64/Fcgr1), marker genes targeting MoDCs (Ccr2, Nrpsa, Syngr2, Itgb7, and Ahr) were highly enriched in cluster 1 with moderate enrichment of several cDC2 markers (Mgl2, Cx3cr1, Irf44, Klf4, and Cd209a), whereas cDC1 marker genes (Xcr1, Clec9a, Gpr141b/A530099J19Rik, and Cadm1) were not widely expressed (Supplemental Fig. 2E). Therefore, our data suggest that cluster 1 contains a mixture of Mϕs, MoDCs, and cDC2s. In addition, mregDC marker genes (Ccr7, Fscn1, Socs2, Relb, and Cacnb) were highly and exclusively expressed in cluster 7, indicating that mregDCs formed the bulk of cluster 7.

Taken together, our data indicate that Mos/Mϕs from early stages of wound healing tended to cluster in different proinflammatory subsets, whereas those from later stages tended to cluster in different prohealing subsets, none of which conformed well to the M1/M2 classification scheme. We also found evidence of APC subsets, including Mϕs, MoDCs, cDC2s, and mregDCs captured by our sorting scheme that tended to remain at more stable levels during the course of healing.

In addition to defining the heterogeneity of wound Mos/Mϕs throughout the course of healing, we also explored the state changes of these cells, or the dynamical changes at the transcriptional level in our scRNA-seq dataset, which could then be tested experimentally. Therefore, cell trajectory analysis was conducted using scVelo to calculate RNA velocity for each cell (35). We observed a flow of cell state transitions from the early-stage clusters (clusters 3, 4, and 6) to the late-stage ones (clusters 0, 2, and 5). Interestingly, a portion of cells in the APC clusters (clusters 1 and 7) also showed potential transition to the clusters associated with the late stage of healing (Fig. 4A). Furthermore, PAGA graph abstraction analysis confirmed such state transitions from early (cluster 3) to late (cluster 0) stages with comparatively high confidence while a weak link was also built from an early-stage cluster (cluster 6) to APCs (Supplemental Fig. 3A).

FIGURE 4.

State transitions of wound Mos/Mϕs revealed by RNA velocity and pseudotime analysis and confirmed by adoptive transfer experiments. (A) RNA velocity analysis indicated cell state transitions between early- and late-stage clusters as well as a portion of APCs to late-stage clusters. Arrows indicate the direction of transitions. (B) Latent time analysis, which is an indicator of an internal clock of the cells, revealed an initial state (low latent time values) within early-stage clusters whereas the final state of cells (high latent time values) appeared to occur in late-stage clusters. (C) Putative driver genes based on likelihood from dynamical modeling of velocity estimates from each cell. Parallel feature plots show associated RNA velocity and expression levels of driver genes. (D) Pseudotime trajectory constructed using clusters 1 and 6 as starting roots (numbers 1 and 2 in circles, respectively), thus giving those clusters low pseudotime values. Two pathways were established between early-stage clusters and APCs to late-stage clusters. Cells are colored by pseudotime. Note that primary termini were identified in late-stage clusters with high pseudotime values but that minor termini were also evident within early-stage clusters and APC cluster 1. (E) Additional pseudotime analysis with Ly6c2hi (expression > median Ly6c2) cells chosen as root cells. This analysis also showed increasing pseudotime from early-stage to late-stage clusters. (F) Adoptive transfer of naive BM live Ly6Gc-KitCD11b+CD115+Ly6Chi Mos to wounds confirming that Mos/Mϕs undergo state changes in wounds from cells expressing high levels of Ly6C to cells expressing low levels of Ly6C, especially at the later healing time point of day 6. Data are mean ± SEM. *p < 0.05 between four groups by two-way ANOVA.

FIGURE 4.

State transitions of wound Mos/Mϕs revealed by RNA velocity and pseudotime analysis and confirmed by adoptive transfer experiments. (A) RNA velocity analysis indicated cell state transitions between early- and late-stage clusters as well as a portion of APCs to late-stage clusters. Arrows indicate the direction of transitions. (B) Latent time analysis, which is an indicator of an internal clock of the cells, revealed an initial state (low latent time values) within early-stage clusters whereas the final state of cells (high latent time values) appeared to occur in late-stage clusters. (C) Putative driver genes based on likelihood from dynamical modeling of velocity estimates from each cell. Parallel feature plots show associated RNA velocity and expression levels of driver genes. (D) Pseudotime trajectory constructed using clusters 1 and 6 as starting roots (numbers 1 and 2 in circles, respectively), thus giving those clusters low pseudotime values. Two pathways were established between early-stage clusters and APCs to late-stage clusters. Cells are colored by pseudotime. Note that primary termini were identified in late-stage clusters with high pseudotime values but that minor termini were also evident within early-stage clusters and APC cluster 1. (E) Additional pseudotime analysis with Ly6c2hi (expression > median Ly6c2) cells chosen as root cells. This analysis also showed increasing pseudotime from early-stage to late-stage clusters. (F) Adoptive transfer of naive BM live Ly6Gc-KitCD11b+CD115+Ly6Chi Mos to wounds confirming that Mos/Mϕs undergo state changes in wounds from cells expressing high levels of Ly6C to cells expressing low levels of Ly6C, especially at the later healing time point of day 6. Data are mean ± SEM. *p < 0.05 between four groups by two-way ANOVA.

Close modal

Latent time, which represents the internal clock of cells, was also calculated to evaluate potential cell trajectories (35). As shown in (Fig. 4B, cluster 6 was assigned the lowest latent time value, suggesting that these cells may represent an initial state of Mos/Mϕs in the wounds. At the opposite end, cells in clusters 0 and 2 exhibited the highest latent time values, indicating that these may represent a final state of Mos/Mϕs in the wounds. In addition, among the early-stage clusters, cells in cluster 4 were assigned higher latent time values compared with the others, indicating a potential cell transition/terminal differentiation state within proinflammatory Mo/Mϕ subsets. This overall trend was also confirmed by PAGA analysis that transitions between three early-stage clusters (Supplemental Fig. 3A). Interestingly, APC clusters showed some evidence of a transition from relatively high to relatively low latent time values (e.g., from right to left in cluster 1), which may also reflect a state change of these cells. Furthermore, analysis of RNA velocity lengths demonstrated a rapid cell transition rate in clusters 4 and 0 associated with high confidence levels (Supplemental Fig. 3B, 3C). Next, we sought to identify the genes that showed the most pronounced dynamical behavior and thus may drive the cell state transitions. The expression dynamics of the top 300 driver genes resolved along latent time indicated a transcriptional cascade that may promote major transition events (Supplemental Fig. 3D). Among the top 10 likelihood-ranked genes, we noticed two major categories: 1) inflammation-related genes (Nfkbia, Nlrp3, Ppbp, Plbd1, Trim30a, and Ifrd1), and 2) Ag presentation–related genes (H2-Ab1, H2-Aa, Cd74, and Ciita), suggesting that these genes may drive transition phases along the early to late stage clusters or within the APC clusters, respectively, during the course of wound healing (Fig. 3C, Supplemental Fig. 3E) (10, 20, 6871).

Next, we ordered all cells in pseudotime to further explore Mo/Mϕ state transitions during the course of wound healing (37). RNA velocity and latent time analysis indicated that cells from clusters 6 and 1 may be in initial/early states (Fig. 4B); therefore, we chose nodes in clusters 6 and 1 as the starting “roots.” Similar to what we observed from RNA velocity and latent time analysis, there were two major trajectories connecting early- to late-stage cells and APC- to late-stage cells with the highest pseudotime values assigned to a group of cells in cluster 2. Several minor branches were also constructed within both the early- and late-stage clusters, indicating that state transitions between Mos/Mϕs may take place along multiple paths with different endpoints (Fig. 3D). To determine whether high-expressing Ly6C cells may serve as an initial state of Mos/Mϕs in wounds, we selected Ly6c2hi cells (expression level > median Ly6c2 expression) from the day 3 cell pool in cluster 3 as root cells in the pseudotime analysis. Again, we observed increased pseudotime values in a group of cells from cluster 4 and to a greater extent in the late-stage clusters, which confirmed transition of Ly6c2hi Mos/Mϕs into different proinflammatory and prohealing phenotypes (Fig. 4E).

In response to damage, Ly6Chi Mos are known to infiltrate skin wounds from peripheral blood (2). To experimentally verify that Ly6Chi cells undergo state transitions after entering the wound, we isolated live Ly6Gc-KitCD11b+CD115+Ly6Chi Mos from BM of the CD45.1 congenic mice (donor) and then injected them intradermally into the wound bed of CD45.2 recipient mice on day 2 or 5 postinjury. After 18 h, cells were dissociated from wounds and labeled with CD45.1 Ab along with Mo/Mϕ Abs F4/80 and Ly6C Abs for analysis by flow cytometry. As shown in (Fig. 4F, only 40% of transferred cells maintained surface expression of Ly6C on day 3 and 20% became Ly6CF4/80+ (mature Mϕs). This trend was even more pronounced in cells transferred to wounds at a later time point (day 6), in that only 35% of transferred cells maintained surface expression of Ly6C and more than half of these cells gained F4/80 expression. Thus, these data confirmed our findings from RNA velocity and pseudotime analysis of the scRNA-seq data (Fig. 4, Supplemental Fig. 4) that Mos/Mϕs undergo state changes in wounds from inflammation-associated Ly6ChiF4/80lo cells to healing-associated Ly6CloF4/80hi cells.

In summary, our data indicate that wound Mos/Mϕs undergo state transitions between various proinflammatory to prohealing phenotypes during the course of skin healing in mice that may influence different aspects of wound healing.

In our dataset, Ccl7 was a signature gene for two early-stage/proinflammatory clusters 3 and 4 (Fig. 2C, Supplemental Fig. 2). However, the role of Ccl7 in skin wound healing has not yet been elucidated. Therefore, we sought to determine the role of Ccl7 in regulating Mo/Mϕ accumulation during wound healing. First, Ccl7 protein levels were measured in uninjured skin and wound homogenates on days 3 and 10 after injury. As shown in (Fig. 5A, Ccl7 protein was significantly increased in wounds on day 3 postinjury, and then dropped significantly toward baseline levels on day 10 postinjury.

FIGURE 5.

Ccl7 is involved in regulation of wound Mo and Mϕ accumulation. (A) Protein levels of Ccl7 measured via custom multiplex kit via flow cytometry using homogenates of uninjured skin (n = 4) and wounds harvested on days 3 and 10 (n = 10 wounds from five mice for each time point). (B) Gating strategy for identifying Ly6C+F4/80lo/− Mos/Mϕs and Ly6CF4/80+ Mϕs in skin wounds. (C and D) Local treatments of Ccl7 recombinant protein (30 ng, black bar; 120 ng, gray bar) on a daily basis significantly increased both (C) percentages and (D) numbers of Ly6C+F4/80lo/− Mos/Mϕs and increased numbers but not percentages of Ly6CF4/80+ Mϕs in wounds on day 3 postinjury compared with those treated with PBS (open bar). n = 6 mice per group from three independent experiments. (E and F) Neutralizing Ccl7 activity in wounds using a Ccl7 blocking Ab (20 µg/wound) on days 0 and 2 lead to (E) a trend toward a lower percentage and (F) a significantly decreased number of Ly6C+F4/80lo/− Mos/Mϕs in wounds on day 3 postinjury compared with IgG control group; in contrast, accumulation of Ly6CF4/80+ Mϕs was not affected. n = 5 mice per group from three independent experiments. Data are mean ± SEM. *p < 0.05, ***p < 0.001 between treatment groups by two-way ANOVA.

FIGURE 5.

Ccl7 is involved in regulation of wound Mo and Mϕ accumulation. (A) Protein levels of Ccl7 measured via custom multiplex kit via flow cytometry using homogenates of uninjured skin (n = 4) and wounds harvested on days 3 and 10 (n = 10 wounds from five mice for each time point). (B) Gating strategy for identifying Ly6C+F4/80lo/− Mos/Mϕs and Ly6CF4/80+ Mϕs in skin wounds. (C and D) Local treatments of Ccl7 recombinant protein (30 ng, black bar; 120 ng, gray bar) on a daily basis significantly increased both (C) percentages and (D) numbers of Ly6C+F4/80lo/− Mos/Mϕs and increased numbers but not percentages of Ly6CF4/80+ Mϕs in wounds on day 3 postinjury compared with those treated with PBS (open bar). n = 6 mice per group from three independent experiments. (E and F) Neutralizing Ccl7 activity in wounds using a Ccl7 blocking Ab (20 µg/wound) on days 0 and 2 lead to (E) a trend toward a lower percentage and (F) a significantly decreased number of Ly6C+F4/80lo/− Mos/Mϕs in wounds on day 3 postinjury compared with IgG control group; in contrast, accumulation of Ly6CF4/80+ Mϕs was not affected. n = 5 mice per group from three independent experiments. Data are mean ± SEM. *p < 0.05, ***p < 0.001 between treatment groups by two-way ANOVA.

Close modal

Next, we tested the effect of exogenous Ccl7 treatment on Mo/Mϕ accumulation in wounds. Recombinant Ccl7 protein was injected intradermally into the wounds at both a low (30 ng) and high dose (120 ng) immediately after wounding and then every day until day 3 postinjury. Ccl7 treatment significantly increased accumulation of the proinflammatory Ly6C+F4/80lo/− Mos/Mϕs by both percentages and numbers. Interestingly, this effect seemed to be smaller for Ly6CF4/80+ Mϕs, as these cells did not show a significant increase in percentage but did show increased cell numbers (Fig. 5B, 5C). Additionally, Ccl7 treatment did not change the percentages of total live cells, CD11b+ myeloid cells, or Ly6G+ neutrophils (Supplemental Fig. 4C–E). The numbers of total live and CD11b+ cells were significantly increased in response to Ccl7 treatment in a dose-dependent manner, which may explain discordance between trends in percentage and number data (Supplemental Fig. 4F, 4G). Interestingly, wounds that received local treatment of Ccl7 showed a trend toward accelerated wound closure compared with the control wounds, although this difference did not reach statistical significance (Supplemental Fig. 4A).

Furthermore, we tested whether loss of Ccl7 activity at the wound site would reduce wound Mo/Mϕ accumulation. Mice received intradermal injections of NAb against Ccl7 to each wound on day 0 (after wounding) and day 2 postinjury. An equal volume of IgG control Ab was administered to the wounds in the control group. On day 3 postinjury, neutralizing Ccl7 in wounds led to a trend toward a lower percentage of Ly6C+F4/80lo/− Mos/Mϕs in wounds compared with IgG controls, whereas Ly6CF4/80+ Mϕs were not affected (Fig. 5E). In addition, the number of Ly6C+F4/80lo/− Mos/Mϕs in wounds treated with Ccl7 NAb was significantly decreased compared with those treated with IgG control. Again, the numbers of Ly6CF4/80+ cells were not affected (Fig. 5F). We did not observe significant changes of either percentage or cell numbers in total live cells, CD11b+ myeloid cells, or Ly6G+ neutrophils in response to local neutralization of Ccl7 in the wounds (Supplemental Fig. 4I–N), suggesting a specific effect on Ly6Chi cells. Moreover, wounds treated with the Ccl7 NAb tended to exhibit slower wound closure compared with IgG control-treated wounds, but again these trends were not statistically significant (Supplemental Fig. 4B).

To further assess expression of selected chemokines and their receptors in Mos/Mϕs, we compared expression of Ccl7 with that of Ccl2 and Ccl8 and assessed expression of their receptors Ccr1, Ccr2, and Ccr3 in our scRNA-seq dataset. As shown in Supplemental Fig. 4O, Ccl2 and Ccl7 were both highly expressed in early-stage Mo/Mϕ clusters whereas Ccl8 expression was limited to cells in late-stage Mo/Mϕ clusters 0 and 5. In addition, Ccr1 and Ccr2 were widely expressed in many Mo/Mϕ/DC clusters whereas Ccr3 expression was barely detected. These data indicate that Ccr1 and Ccr2 may mediate effects of Ccl7 on Ly6C+F4/80lo/− Mos/Mϕs we observed in our experiment.

Taken together, our data indicated that Ccl7 protein preferentially affects the accumulation of proinflammatory Ly6C+F4/80lo/− Mos/Mϕs in mouse skin wounds, which in turn are the primary expressors of Ccl7 in wounds. These results also highlight the power of scRNA-seq for identifying novel regulators of inflammation and wound healing that can be tested experimentally.

Mϕs are capable of performing diverse functions in many physiological processes that are associated with their ability to adopt a wide spectrum of phenotypes (72). In the context of wound healing, the remarkable plasticity of these cells enables them to orchestrate many key biological processes involved in the inflammatory, proliferative, and remodeling phases of healing (1, 2, 6). In the current study, we used scRNA-seq to comprehensively investigate the heterogeneity of skin wound Mos/Mϕs during the course of skin wound healing. Targeting Ly6GCD11b+CD45+ cells isolated from skin wounds on days 3, 6, and 10 postinjury, we identified eight subsets of Mos/Mϕs that varied in their contributions to the wound Mo/Mϕ population during the course of healing. Based on cluster composition and DE analysis, the eight clusters could be categorized into three groups that play different roles in the progression of healing: early-stage clusters (proinflammatory), late-stage clusters (prohealing), and APCs. Moreover, cell trajectory analysis by RNA velocity and pseudotime revealed dynamical state transitions not only between early and late stages of cells but also among proinflammatory Mos/Mϕs within early-stage clusters. Adoptive transfer of naive BM Ly6Chi Mos into the wound bed confirmed the transition of Ly6C+ Mos into Ly6CF4/80+ Mϕs as wound healing progresses. Finally, we demonstrated that Ccl7, a signature gene for the early stage of Mo/Mϕ clusters, preferentially induces infiltration of proinflammatory Ly6C+ Mos/Mϕs into wounds at the early stage of healing.

Our group and others have reported that Mos/Mϕs transition from proinflammatory to prohealing phenotypes and that this transition is required for efficient healing (2, 5, 12, 32). However, these studies were limited by the small numbers of markers used to characterize phenotypes. To provide a comprehensive picture of Mo/Mϕ heterogeneity during the full healing course, we performed single-cell transcriptomics of Mos/Mϕs isolated from skin wounds. Early-stage clusters showed upregulation of proinflammatory genes related to the type I IFN signaling pathway including IFN regulatory protein Irf7, IFN stimulated gene Isg15, and IFN-induced transmembrane protein Ifitm6, which likely indicates communication between inflammatory cells during the healing process (73, 74). Another study identified a similar group of Mϕs in atherosclerotic plaques, which they named IFN signaturehi Mϕs (enriched with Irf7, Isg15, Ifitm3, and Ly6e), which were associated with disease progression in mice (26). Moreover, signature genes and GO analysis of cluster 6 in our study also characterized a group of Mos/Mϕs actively involved in cytokine synthesis and endocytosis, which may contribute to the initial cascade of inflammation and cytokine storm as well as continuous destructive killing tasks in the wounds. These findings support our interpretation that early-stage Mos/Mϕs promote inflammation in wound healing. During the late stage of wound healing, Mos/Mϕs exhibited upregulation of tissue repair and remodeling genes (Pdpn, Mmp12, Arg-1, and Timp-2). Notably, genes associated with the complement system were enriched in Mos/Mϕs from late-stage clusters (clusters 0 and 5). Interestingly, a comparable population of such Mϕs was also identified in murine aortic atherosclerosis as a protective cell subset (10), and previous reports indicated beneficial effects of complement-activating components on wound healing (75), supporting the idea that the late-stage clusters have a prohealing phenotype. Furthermore, GO analysis also revealed enrichment of terms associated with angiogenesis (clusters 0 and 2), endocytosis (cluster 0), and regulation of signaling and cell communication (cluster 2), indicating prohealing and phagocytic Mos/Mϕs dominating during later stages of healing. These data are consistent with a previous report suggesting that late-stage Mϕs regulate angiogenesis by clearing endothelial cells and other debris (76).

Mϕs can function as professional APCs that assist in the propagation of inflammation (77). A distinct cluster in our data exhibited upregulation of MHC class II–associated genes (cluster 1) and had relatively equal contributions from all three healing time points. Marker gene expression indicated that there was a mixture of APC subsets including Mϕs, MoDCs, and cDC2s, which may indicate that Ag-presenting Mos/Mϕs collaborate with DCs and participate in inflammatory, proliferative, and remodeling phases of healing. In a murine atherosclerosis model, such MHC class II gene-expressing Mϕs contribute to cell communication and disease progression (10, 26). Further study on the function of each of the major Mo/Mϕ groups identified in our study is warranted and we are currently developing methods to target each group.

After infiltrating the wound bed, cells from the Mo/Mϕ lineage likely pursue different cell fates through state transitions during the healing process. RNA velocity, PAGA, and pseudotime analyses all demonstrated cell state transitions from early proinflammatory to late prohealing Mos/Mϕs. Adoptive transfer experiments confirmed state transitions of naive Ly6Chi Mos introduced into the wound bed, which showed loss of surface Ly6C expression but increased F4/80 expression, indicating a transition of Mos/Mϕs in wounds from proinflammatory to a mature phenotype. Importantly, two major gene categories appeared to drive state transitions predicted by dynamical modeling of RNA velocities: genes associated with inflammation and genes associated with Ag presentation. Among the inflammation-related genes, NLRP3 inflammasome and NF-κB signaling pathways were highlighted, and the importance of both has been well documented as key regulators of inflammation in wound healing (19, 78, 79), but the present analysis indicates that these genes may also play important roles in state transitions of early- to late-stage Mos/Mϕs during the course of healing. We also speculate that Ag-presenting genes may play important roles in state transitions between APC clusters and early- or late-stage Mos/Mϕs. These data can inform mechanistic testing of the importance of state transitions between Mo/Mϕ subsets in different aspects of wound healing as well as the drivers of these state transitions.

Finally, we tested the influence of a key signature gene of early-stage clusters, Ccl7, on wound inflammation. Ccl7 is a chemokine essential for mobilization of Ly6Chi Mos out of BM (80), but it was also reported to negatively regulate skin inflammation induced by Leishmania major infection in mice (81). In our study, excisional wounding significantly upregulated Ccl7 protein in wound homogenate peaking on day 3 postinjury, consistent with enrichment of Ccl7 gene expression in early-stage skin wound Mos/Mϕs. Daily treatment with exogenous Ccl7 recombinant protein administered locally to wounds preferentially increased proinflammatory Ly6C+F4/80lo/− Mos/Mϕs in wounds. In addition, blocking Ccl7 activity locally in wounds using a neutralizing Ab reduced the accumulation of Ly6C+F4/80lo/− Mos/Mϕs. Taken together, these data suggested a role of Ccl7 in the accumulation of proinflammatory Mos/Mϕs during the early stage of healing. These results contrast with findings from a study on L. major skin infection that Ccl7 knockout mice exhibit enhanced dermal inflammation, including increased accumulation of neutrophils, suggesting that Ccl7 plays an anti-inflammatory role in this model (81). This discrepancy may be explained by different models used, that is, infection versus excisional wounding in the current study. Interestingly, we observed a trend toward faster wound closure with added Ccl7 protein and delayed closure when blocking Ccl7 activity, suggesting that Ccl7 promotes wound closure in the early stage of healing; future studies will be designed to more fully explore the role of Ccl7 in skin wound healing.

The recent developments and implementation of scRNA-seq are beginning to expand our knowledge of cell heterogeneity during skin wound healing. Guerrero-Juarez et al. (27) described fibroblast heterogeneity in large skin wounds of mice on day 12 postinjury, a time point associated with full re-epithelialization. Their study indicated that wound fibroblasts group into 12 clusters, some of which likely represent differentiation states toward a contractile phenotype, whereas others appear to represent distinct lineages. Willenborg et al. (31) reported on the heterogeneity of wound CD11b+F4/80+ Mϕs isolated from mouse wounds on days 4 and 14 postinjury. Their study identified seven clusters of Mϕs that demonstrated heterogeneity across activation patterns and metabolic signatures. Also, recent studies have recently reported on wound cell heterogeneity in diabetic foot ulcers in humans (82, 83). Our study adds to this literature by describing the heterogeneity of Mos and Mϕs through all three stages of healing; the study by Willenborg et al. likely excluded some of the Mos and APC phenotypes from their analysis by focusing on F4/80+ cells. In addition, our cell trajectory analysis revealed state transitions of Mos/Mϕs at different stages of healing as well as the potential drivers of these transitions.

One limitation of our study is that our scRNA-seq data cannot distinguish between resident versus infiltrated Mos/Mϕs. Current understanding of wound healing indicates that most Mos/Mϕs accumulating in skin wounds are derived from circulating Mos with a minor contribution from resident cells (12). Therefore, we will plan future experiments combining lineage tracing and scRNA-seq to directly assess the contributions of resident versus infiltrating cells to Mo/Mϕ heterogeneity during wound healing. Moreover, our study identified some unique gene markers that indicated distinct characteristics and functions of Mo/Mϕ subsets during the healing process. Thus, we plan to test new protein markers based on our scRNA-seq data in future experiments to allow sorting and further characterization of the subsets identified. Finally, we note that Ccl7 does not likely attract Mos/Mϕs to wounds on its own; for example, Ccl2 is well known for its chemoattractive effects on Mos/Mϕs. In the future, we plan to conduct a comprehensive screening of chemokines, both singly and in combinations, for their impact on Mo/Mϕ infiltration and skin during healing.

In summary, our study used scRNA-seq and downstream analysis to demonstrate the complexity of skin Mo/Mϕ heterogeneity and state transitions throughout the course of wound healing. After injury, Mos/Mϕs first adopt early-stage proinflammatory phenotypes and then appear to transition to later-stage prohealing phenotypes as healing progresses. Potential transitions between APC-like cells and early- and late-stage Mos/Mϕs were also suggested by the data, and some of these state transitions were confirmed by adoptive transfer experiments. These state transitions were likely driven by genes associated with inflammation and Ag presentation, respectively. Finally, we showed that early upregulation of the chemokine Ccl7 contributes to the accumulation of proinflammatory Ly6C+F4/80lo/− Mos/Mϕs in mouse skin wounds. These data provide insight into designing future studies to mechanistically investigate factors that regulate Mo/Mϕ heterogeneity and state transitions during wound healing.

We thank Dr. Luisa A. DiPietro and Dr. Giamila Fantuzzi for input on aspects of presentation of this manuscript.

This work was supported by National Institute of General Medical Sciences Grant R35GM136228 (to T.J.K.). M.M.-C. was supported in part by National Center for Advancing Translational Sciences Grant UL1TR002003.

The online version of this article contains supplemental material.

J.P. designed the study, conducted the experiments, performed data and statistical analysis, and wrote the manuscript; M.M.-C. helped design the study, performed data and statistical analysis, and helped write the manuscript; and T.J.K. helped design the study, write the manuscript, and provided all materials.

The RNA-seq data presented in this article have been submitted to the Gene Expression Omnibus under accession number GSE203244.

Abbreviations used in this article:

     
  • ANOSIM

    analysis of similarities

  •  
  • BM

    bone marrow

  •  
  • cDC1

    conventional type 1 DC

  •  
  • cDC2

    conventional type 2 DC

  •  
  • DC

    dendritic cell

  •  
  • DE

    differentially expressed

  •  
  • DPBS

    Dulbecco’s PBS

  •  
  • GEM

    gel beads in emulsion

  •  
  • GO

    Gene Ontology

  •  
  • macrophage

  •  
  • Mo

    monocyte

  •  
  • MoDC

    Mo-derived DC

  •  
  • mregDC

    mature DC enriched in immunoregulatory molecules

  •  
  • NAb

    neutralizing Ab

  •  
  • PAGA

    partition-based graph abstraction

  •  
  • qPCR

    quantitative PCR

  •  
  • RNA-seq

    RNA sequencing

  •  
  • scRNA-seq

    single-cell RNA-seq

  •  
  • UMAP

    uniform manifold approximation and projection

1.
Eming
S. A.
,
P.
Martin
,
M.
Tomic-Canic
.
2014
.
Wound repair and regeneration: mechanisms, signaling, and translation.
Sci. Transl. Med.
6
:
265sr6
.
2.
Koh
T. J.
,
L. A.
DiPietro
.
2011
.
Inflammation and wound healing: the role of the macrophage.
Expert Rev. Mol. Med.
13
:
e23
.
3.
Leoni
G.
,
P. A.
Neumann
,
R.
Sumagin
,
T. L.
Denning
,
A.
Nusrat
.
2015
.
Wound repair: role of immune-epithelial interactions.
Mucosal Immunol.
8
:
959
968
.
4.
Gordon
S.
,
A.
Plüddemann
.
2017
.
Tissue macrophages: heterogeneity and functions.
BMC Biol.
15
:
53
.
5.
Krzyszczyk
P.
,
R.
Schloss
,
A.
Palmer
,
F.
Berthiaume
.
2018
.
The role of macrophages in acute and chronic wound healing and interventions to promote pro-wound healing phenotypes.
Front. Physiol.
9
:
419
.
6.
Novak
M. L.
,
T. J.
Koh
.
2013
.
Phenotypic transitions of macrophages orchestrate tissue repair.
Am. J. Pathol.
183
:
1352
1363
.
7.
Murray
P. J.
,
J. E.
Allen
,
S. K.
Biswas
,
E. A.
Fisher
,
D. W.
Gilroy
,
S.
Goerdt
,
S.
Gordon
,
J. A.
Hamilton
,
L. B.
Ivashkiv
,
T.
Lawrence
, et al
2014
.
Macrophage activation and polarization: nomenclature and experimental guidelines. [Published erratum appears in 2014 Immunity 41: 339–340.]
Immunity
41
:
14
20
.
8.
Martinez
F. O.
,
S.
Gordon
.
2014
.
The M1 and M2 paradigm of macrophage activation: time for reassessment.
F1000Prime Rep.
6
:
13
.
9.
Chávez-Galán
L.
,
M. L.
Olleros
,
D.
Vesin
,
I.
Garcia
.
2015
.
Much more than M1 and M2 macrophages, there are also CD169+ and TCR+ macrophages.
Front. Immunol.
6
:
263
.
10.
Cochain
C.
,
E.
Vafadarnejad
,
P.
Arampatzi
,
J.
Pelisek
,
H.
Winkels
,
K.
Ley
,
D.
Wolf
,
A. E.
Saliba
,
A.
Zernecke
.
2018
.
Single-cell RNA-seq reveals the transcriptional landscape and heterogeneity of aortic macrophages in murine atherosclerosis.
Circ. Res.
122
:
1661
1674
.
11.
Hoeksema
M. A.
,
C. K.
Glass
.
2019
.
Nature and nurture of tissue-specific macrophage phenotypes.
Atherosclerosis
281
:
159
167
.
12.
Burgess
M.
,
K.
Wicks
,
M.
Gardasevic
,
K. A.
Mace
.
2019
.
Cx3CR1 expression identifies distinct macrophage populations that contribute differentially to inflammation and repair.
Immunohorizons
3
:
262
273
.
13.
Daley
J. M.
,
S. K.
Brancato
,
A. A.
Thomay
,
J. S.
Reichner
,
J. E.
Albina
.
2010
.
The phenotype of murine wound macrophages.
J. Leukoc. Biol.
87
:
59
67
.
14.
Novak
M. L.
,
T. J.
Koh
.
2013
.
Macrophage phenotypes during tissue repair.
J. Leukoc. Biol.
93
:
875
881
.
15.
Shook
B.
,
E.
Xiao
,
Y.
Kumamoto
,
A.
Iwasaki
,
V.
Horsley
.
2016
.
CD301b+ macrophages are essential for effective skin wound healing.
J. Invest. Dermatol.
136
:
1885
1891
.
16.
Pang
J.
,
N.
Urao
,
T. J.
Koh
.
2020
.
Proliferation of Ly6C+ monocytes/macrophages contributes to their accumulation in mouse skin wounds.
J. Leukoc. Biol.
107
:
551
560
.
17.
Lurier
E. B.
,
D.
Dalton
,
W.
Dampier
,
P.
Raman
,
S.
Nassiri
,
N. M.
Ferraro
,
R.
Rajagopalan
,
M.
Sarmady
,
K. L.
Spiller
.
2017
.
Transcriptome analysis of IL-10-stimulated (M2c) macrophages by next-generation sequencing.
Immunobiology
222
:
847
856
.
18.
Das
A.
,
C. S.
Yang
,
S.
Arifuzzaman
,
S.
Kim
,
S. Y.
Kim
,
K. H.
Jung
,
Y. S.
Lee
,
Y. G.
Chai
.
2018
.
High-resolution mapping and dynamics of the transcriptome, transcription factors, and transcription co-factor networks in classically and alternatively activated macrophages.
Front. Immunol.
9
:
22
.
19.
Mirza
R. E.
,
M. M.
Fang
,
E. M.
Weinheimer-Haus
,
W. J.
Ennis
,
T. J.
Koh
.
2014
.
Sustained inflammasome activity in macrophages impairs wound healing in type 2 diabetic humans and mice.
Diabetes
63
:
1103
1114
.
20.
Beswick
E. J.
,
V. E.
Reyes
.
2009
.
CD74 in antigen presentation, inflammation, and cancers of the gastrointestinal tract.
World J. Gastroenterol.
15
:
2855
2861
.
21.
Goldman
S. L.
,
M.
MacKay
,
E.
Afshinnekoo
,
A. M.
Melnick
,
S.
Wu
,
C. E.
Mason
.
2019
.
The impact of heterogeneity on single-cell sequencing.
Front. Genet.
10
:
8
.
22.
Papalexi
E.
,
R.
Satija
.
2018
.
Single-cell RNA sequencing to explore immune cell heterogeneity.
Nat. Rev. Immunol.
18
:
35
45
.
23.
Milich
L. M.
,
J. S.
Choi
,
C.
Ryan
,
S. R.
Cerqueira
,
S.
Benavides
,
S. L.
Yahn
,
P.
Tsoulfas
,
J. K.
Lee
.
2021
.
Single-cell analysis of the cellular heterogeneity and interactions in the injured mouse spinal cord.
J. Exp. Med.
218
:
e20210040
.
24.
Mould
K. J.
,
N. D.
Jackson
,
P. M.
Henson
,
M.
Seibold
,
W. J.
Janssen
.
2019
.
Single cell RNA sequencing identifies unique inflammatory airspace macrophage subsets.
JCI Insight
4
:
e126556
.
25.
Yang
Q.
,
H.
Zhang
,
T.
Wei
,
A.
Lin
,
Y.
Sun
,
P.
Luo
,
J.
Zhang
.
2021
.
Single-cell RNA sequencing reveals the heterogeneity of tumor-associated macrophage in non-small cell lung cancer and differences between sexes.
Front. Immunol.
12
:
756722
.
26.
Lin
J. D.
,
H.
Nishi
,
J.
Poles
,
X.
Niu
,
C.
Mccauley
,
K.
Rahman
,
E. J.
Brown
,
S. T.
Yeung
,
N.
Vozhilla
,
A.
Weinstock
, et al
2019
.
Single-cell analysis of fate-mapped macrophages reveals heterogeneity, including stem-like properties, during atherosclerosis progression and regression.
JCI Insight
4
:
e124574
.
27.
Guerrero-Juarez
C. F.
,
P. H.
Dedhia
,
S.
Jin
,
R.
Ruiz-Vega
,
D.
Ma
,
Y.
Liu
,
K.
Yamaga
,
O.
Shestova
,
D. L.
Gay
,
Z.
Yang
, et al
2019
.
Single-cell analysis reveals fibroblast heterogeneity and myeloid-derived adipocyte progenitors in murine skin wounds.
Nat. Commun.
10
:
650
.
28.
Haensel
D.
,
S.
Jin
,
P.
Sun
,
R.
Cinco
,
M.
Dragan
,
Q.
Nguyen
,
Z.
Cang
,
Y.
Gong
,
R.
Vu
,
A. L.
MacLean
, et al
2020
.
Defining epidermal basal cell states during skin homeostasis and wound healing using single-cell transcriptomics.
Cell Rep.
30
:
3932
3947.e6
.
29.
Joost
S.
,
T.
Jacob
,
X.
Sun
,
K.
Annusver
,
G.
La Manno
,
I.
Sur
,
M.
Kasper
.
2018
.
Single-cell transcriptomics of traced epidermal and hair follicle stem cells reveals rapid adaptations during wound healing.
Cell Rep.
25
:
585
597.e7
.
30.
Li
D.
,
S.
Cheng
,
Y.
Pei
,
P.
Sommar
,
J.
Kärner
,
E. K.
Herter
,
M. A.
Toma
,
L.
Zhang
,
K.
Pham
,
Y. T.
Cheung
, et al
2022
.
Single-cell analysis reveals major histocompatibility complex II‒Expressing keratinocytes in pressure ulcers with worse healing outcomes.
J. Invest. Dermatol.
142
(
3 Pt A
):
705
716
.
31.
Willenborg
S.
,
D. E.
Sanin
,
A.
Jais
,
X.
Ding
,
T.
Ulas
,
J.
Nüchel
,
M.
Popović
,
T.
MacVicar
,
T.
Langer
,
J. L.
Schultze
, et al
2021
.
Mitochondrial metabolism coordinates stage-specific repair processes in macrophages during wound healing.
Cell Metab.
33
:
2398
2414.e9
.
32.
Pang
J.
,
M.
Maienschein-Cline
,
T. J.
Koh
.
2021
.
Enhanced proliferation of Ly6C+ monocytes/macrophages contributes to chronic inflammation in skin wounds of diabetic mice.
J. Immunol.
206
:
621
630
.
33.
Stuart
T.
,
A.
Butler
,
P.
Hoffman
,
C.
Hafemeister
,
E.
Papalexi
,
W. M.
Mauck
3rd
,
Y.
Hao
,
M.
Stoeckius
,
P.
Smibert
,
R.
Satija
.
2019
.
Comprehensive integration of single-cell data.
Cell
177
:
1888
1902.e21
.
34.
Hosack
D. A.
,
G.
Dennis
Jr.
,
B. T.
Sherman
,
H. C.
Lane
,
R. A.
Lempicki
.
2003
.
Identifying biological themes within lists of genes with EASE.
Genome Biol.
4
:
R70
.
35.
Bergen
V.
,
M.
Lange
,
S.
Peidli
,
F. A.
Wolf
,
F. J.
Theis
.
2020
.
Generalizing RNA velocity to transient cell states through dynamical modeling.
Nat. Biotechnol.
38
:
1408
1414
.
36.
Wolf
F. A.
,
F. K.
Hamey
,
M.
Plass
,
J.
Solana
,
J. S.
Dahlin
,
B.
Göttgens
,
N.
Rajewsky
,
L.
Simon
,
F. J.
Theis
.
2019
.
PAGA: graph abstraction reconciles clustering with trajectory inference through a topology preserving map of single cells.
Genome Biol.
20
:
59
.
37.
Cao
J.
,
M.
Spielmann
,
X.
Qiu
,
X.
Huang
,
D. M.
Ibrahim
,
A. J.
Hill
,
F.
Zhang
,
S.
Mundlos
,
L.
Christiansen
,
F. J.
Steemers
, et al
2019
.
The single-cell transcriptional landscape of mammalian organogenesis.
Nature
566
:
496
502
.
38.
Drescher
H. K.
,
E. F.
Brandt
,
P.
Fischer
,
S.
Dreschers
,
R. A.
Schwendener
,
M. A.
Kowalska
,
A.
Canbay
,
H. E.
Wasmuth
,
R.
Weiskirchen
,
C.
Trautwein
, et al
2019
.
Platelet factor 4 attenuates experimental acute liver injury in mice.
Front. Physiol.
10
:
326
.
39.
Gleissner
C. A.
,
I.
Shaked
,
K. M.
Little
,
K.
Ley
.
2010
.
CXC chemokine ligand 4 induces a unique transcriptome in monocyte-derived macrophages.
J. Immunol.
184
:
4810
4818
.
40.
Huang
Y.
,
R. W.
Mahley
.
2014
.
Apolipoprotein E: structure and function in lipid metabolism, neurobiology, and Alzheimer's diseases.
Neurobiol. Dis.
72
(
Pt A
):
3
12
.
41.
Lee
E. J.
,
H. S.
Kim
.
2014
.
The anti-inflammatory role of tissue inhibitor of metalloproteinase-2 in lipopolysaccharide-stimulated microglia.
J. Neuroinflammation
11
:
116
.
42.
Man
S. M.
,
T. D.
Kanneganti
.
2016
.
Regulation of lysosomal dynamics and autophagy by CTSB/cathepsin B.
Autophagy
12
:
2504
2505
.
43.
Sun
W.
,
Y.
Lin
,
L.
Chen
,
R.
Ma
,
J.
Cao
,
J.
Yao
,
K.
Chen
,
J.
Wan
.
2018
.
Legumain suppresses OxLDL-induced macrophage apoptosis through enhancement of the autophagy pathway.
Gene
652
:
16
24
.
44.
Jofre-Monseny
L.
,
A.
Loboda
,
A. E.
Wagner
,
P.
Huebbe
,
C.
Boesch-Saadatmandi
,
A.
Jozkowicz
,
A. M.
Minihane
,
J.
Dulak
,
G.
Rimbach
.
2007
.
Effects of apoE genotype on macrophage inflammation and heme oxygenase-1 expression.
Biochem. Biophys. Res. Commun.
357
:
319
324
.
45.
Bieniasz-Krzywiec
P.
,
R.
Martín-Pérez
,
M.
Ehling
,
M.
García-Caballero
,
S.
Pinioti
,
S.
Pretto
,
R.
Kroes
,
C.
Aldeni
,
M.
Di Matteo
,
H.
Prenen
, et al
2019
.
Podoplanin-expressing macrophages promote lymphangiogenesis and lymphoinvasion in breast cancer.
Cell Metab.
30
:
917
936.e10
.
46.
Chan
M. F.
,
J.
Li
,
A.
Bertrand
,
A. J.
Casbon
,
J. H.
Lin
,
I.
Maltseva
,
Z.
Werb
.
2013
.
Protective effects of matrix metalloproteinase-12 following corneal injury.
J. Cell Sci.
126
:
3948
3960
.
47.
Campbell
L.
,
C. R.
Saville
,
P. J.
Murray
,
S. M.
Cruickshank
,
M. J.
Hardman
.
2013
.
Local arginase 1 activity is required for cutaneous wound healing.
J. Invest. Dermatol.
133
:
2461
2470
.
48.
Hu
J. M.
,
K.
Liu
,
J. H.
Liu
,
X. L.
Jiang
,
X. L.
Wang
,
Y. Z.
Chen
,
S. G.
Li
,
H.
Zou
,
L. J.
Pang
,
C. X.
Liu
, et al
2017
.
CD163 as a marker of M2 macrophage, contribute to predicte aggressiveness and prognosis of Kazakh esophageal squamous cell carcinoma.
Oncotarget
8
:
21526
21538
.
49.
Girkin
J.
,
L.
Hatchwell
,
P.
Foster
,
S. L.
Johnston
,
N.
Bartlett
,
A.
Collison
,
J.
Mattes
.
2015
.
CCL7 and IRF-7 mediate hallmark inflammatory and IFN responses following rhinovirus 1B infection.
J. Immunol.
194
:
4924
4930
.
50.
Senatus
L.
,
R.
López-Díez
,
L.
Egaña-Gorroño
,
J.
Liu
,
J.
Hu
,
G.
Daffu
,
Q.
Li
,
K.
Rahman
,
Y.
Vengrenyuk
,
T. J.
Barrett
, et al
2020
.
RAGE impairs murine diabetic atherosclerosis regression and implicates IRF7 in macrophage inflammation and cholesterol metabolism.
JCI Insight
5
:
e137289
.
51.
Upadhyay
G.
2019
.
Emerging role of lymphocyte antigen-6 family of genes in cancer and immune cells.
Front. Immunol.
10
:
819
.
52.
Rashighi
M.
,
P.
Agarwal
,
J. M.
Richmond
,
T. H.
Harris
,
K.
Dresser
,
M. W.
Su
,
Y.
Zhou
,
A.
Deng
,
C. A.
Hunter
,
A. D.
Luster
,
J. E.
Harris
.
2014
.
CXCL10 is critical for the progression and maintenance of depigmentation in a mouse model of vitiligo.
Sci. Transl. Med.
6
:
223ra23
.
53.
Carrero
R.
,
I.
Cerrada
,
E.
Lledó
,
J.
Dopazo
,
F.
García-García
,
M. P.
Rubio
,
C.
Trigueros
,
A.
Dorronsoro
,
A.
Ruiz-Sauri
,
J. A.
Montero
,
P.
Sepúlveda
.
2012
.
IL1β induces mesenchymal stem cells migration and leucocyte chemotaxis through NF-κB.
Stem Cell Rev. Rep.
8
:
905
916
.
54.
Kyriakides
T. R.
,
S.
Maclauchlan
.
2009
.
The role of thrombospondins in wound healing, ischemia, and the foreign body reaction.
J. Cell Commun. Signal.
3
:
215
225
.
55.
Wight
T. N.
,
I.
Kang
,
S. P.
Evanko
,
I. A.
Harten
,
M. Y.
Chang
,
O. M. T.
Pearce
,
C. E.
Allen
,
C. W.
Frevert
.
2020
.
Versican—a critical extracellular matrix regulator of immunity and inflammation.
Front. Immunol.
11
:
512
.
56.
Han
J. H.
,
S.
Lee
,
Y. S.
Park
,
J. S.
Park
,
K. Y.
Kim
,
J. S.
Lim
,
K. S.
Oh
,
Y.
Yang
.
2011
.
IFITM6 expression is increased in macrophages of tumor-bearing mice.
Oncol. Rep.
25
:
531
536
.
57.
Asleh
R.
,
J.
Ward
,
N. S.
Levy
,
S.
Safuri
,
D.
Aronson
,
A. P.
Levy
.
2014
.
Haptoglobin genotype-dependent differences in macrophage lysosomal oxidative injury.
J. Biol. Chem.
289
:
16313
16325
.
58.
Shi
J.
,
H. L.
Karlsson
,
K.
Johansson
,
V.
Gogvadze
,
L.
Xiao
,
J.
Li
,
T.
Burks
,
A.
Garcia-Bennett
,
A.
Uheida
,
M.
Muhammed
, et al
2012
.
Microsomal glutathione transferase 1 protects against toxicity induced by silica nanoparticles but not by zinc oxide nanoparticles.
ACS Nano
6
:
1925
1938
.
59.
Lee
B. C.
,
S. G.
Lee
,
M. K.
Choo
,
J. H.
Kim
,
H. M.
Lee
,
S.
Kim
,
D. E.
Fomenko
,
H. Y.
Kim
,
J. M.
Park
,
V. N.
Gladyshev
.
2017
.
Selenoprotein MsrB1 promotes anti-inflammatory cytokine gene expression in macrophages and controls immune response in vivo.
Sci. Rep.
7
:
5119
.
60.
Stables
M. J.
,
S.
Shah
,
E. B.
Camon
,
R. C.
Lovering
,
J.
Newson
,
J.
Bystrom
,
S.
Farrow
,
D. W.
Gilroy
.
2011
.
Transcriptomic analyses of murine resolution-phase macrophages.
Blood
118
:
e192
e208
.
61.
Yamakita
Y.
,
F.
Matsumura
,
M. W.
Lipscomb
,
P. C.
Chou
,
G.
Werlen
,
J. K.
Burkhardt
,
S.
Yamashiro
.
2011
.
Fascin1 promotes cell migration of mature dendritic cells.
J. Immunol.
186
:
2850
2859
.
62.
Kimball
A.
,
M.
Schaller
,
A.
Joshi
,
F. M.
Davis
,
A.
denDekker
,
A.
Boniakowski
,
J.
Bermick
,
A.
Obi
,
B.
Moore
,
P. K.
Henke
, et al
2018
.
Ly6CHi blood monocyte/macrophage drive chronic inflammation and impair wound healing in diabetes mellitus.
Arterioscler. Thromb. Vasc. Biol.
38
:
1102
1114
.
63.
Mantovani
A.
,
A.
Sica
,
M.
Locati
.
2005
.
Macrophage polarization comes of age.
Immunity
23
:
344
346
.
64.
Kashem
S. W.
,
M.
Haniffa
,
D. H.
Kaplan
.
2017
.
Antigen-presenting cells in the skin.
Annu. Rev. Immunol.
35
:
469
499
.
65.
Maier
B.
,
A. M.
Leader
,
S. T.
Chen
,
N.
Tung
,
C.
Chang
,
J.
LeBerichel
,
A.
Chudnovskiy
,
S.
Maskey
,
L.
Walker
,
J. P.
Finnigan
, et al
2020
.
A conserved dendritic-cell regulatory program limits antitumour immunity. [Published erratum appears in 2020 Nature 582: E17.]
Nature
580
:
257
262
.
66.
Malissen
B.
,
S.
Tamoutounour
,
S.
Henri
.
2014
.
The origins and functions of dendritic cells and macrophages in the skin.
Nat. Rev. Immunol.
14
:
417
428
.
67.
Anderson
III
D. A.
,
C. A.
Dutertre
,
F.
Ginhoux
,
K. M.
Murphy
.
2021
.
Genetic models of human and mouse dendritic cell development and function.
Nat. Rev. Immunol.
21
:
101
115
.
68.
Yeo
L.
,
N.
Adlard
,
M.
Biehl
,
M.
Juarez
,
T.
Smallie
,
M.
Snow
,
C. D.
Buckley
,
K.
Raza
,
A.
Filer
,
D.
Scheel-Toellner
.
2016
.
Expression of chemokines CXCL4 and CXCL7 by synovial macrophages defines an early stage of rheumatoid arthritis.
Ann. Rheum. Dis.
75
:
763
771
.
69.
Durocher
M.
,
B. P.
Ander
,
G.
Jickling
,
F.
Hamade
,
H.
Hull
,
B.
Knepp
,
D. Z.
Liu
,
X.
Zhan
,
A.
Tran
,
X.
Cheng
, et al
2019
.
Inflammatory, regulatory, and autophagy co-expression modules and hub genes underlie the peripheral immune response to human intracerebral hemorrhage.
J. Neuroinflammation
16
:
56
.
70.
Shi
M.
,
W.
Deng
,
E.
Bi
,
K.
Mao
,
Y.
Ji
,
G.
Lin
,
X.
Wu
,
Z.
Tao
,
Z.
Li
,
X.
Cai
, et al
2008
.
TRIM30α negatively regulates TLR-mediated NF-κB activation by targeting TAB2 and TAB3 for degradation.
Nat. Immunol.
9
:
369
377
.
71.
Buxadé
M.
,
H.
Huerga Encabo
,
M.
Riera-Borrull
,
L.
Quintana-Gallardo
,
P.
López-Cotarelo
,
M.
Tellechea
,
S.
Martínez-Martínez
,
J. M.
Redondo
,
J.
Martín-Caballero
,
J. M.
Flores
, et al
2018
.
Macrophage-specific MHCII expression is regulated by a remote Ciita enhancer controlled by NFAT5.
J. Exp. Med.
215
:
2901
2918
.
72.
Tabas
I.
,
K. E.
Bornfeldt
.
2016
.
Macrophage phenotype and function in different stages of atherosclerosis.
Circ. Res.
118
:
653
667
.
73.
Kanno
E.
,
H.
Tanno
,
A.
Masaki
,
A.
Sasaki
,
N.
Sato
,
M.
Goto
,
M.
Shisai
,
K.
Yamaguchi
,
N.
Takagi
,
M.
Shoji
, et al
2019
.
Defect of interferon γ leads to impaired wound healing through prolonged neutrophilic inflammatory response and enhanced MMP-2 activation.
Int. J. Mol. Sci.
20
:
5657
.
74.
Ishida
Y.
,
T.
Kondo
,
T.
Takayasu
,
Y.
Iwakura
,
N.
Mukaida
.
2004
.
The essential involvement of cross-talk between IFN-γ and TGF-β in the skin wound-healing process.
J. Immunol.
172
:
1848
1855
.
75.
Sinno
H.
,
S.
Prakash
.
2013
.
Complements and the wound healing cascade: an updated review.
Plast. Surg. Int.
2013
:
146764
.
76.
Gurevich
D. B.
,
C. E.
Severn
,
C.
Twomey
,
A.
Greenhough
,
J.
Cash
,
A. M.
Toye
,
H.
Mellor
,
P.
Martin
.
2018
.
Live imaging of wound angiogenesis reveals macrophage orchestrated vessel sprouting and regression.
EMBO J.
37
:
e97786
.
77.
Barker
R. N.
,
L. P.
Erwig
,
K. S.
Hill
,
A.
Devine
,
W. P.
Pearce
,
A. J.
Rees
.
2002
.
Antigen presentation by macrophages is enhanced by the uptake of necrotic, but not apoptotic, cells.
Clin. Exp. Immunol.
127
:
220
225
.
78.
Ito
H.
,
A.
Kanbe
,
H.
Sakai
,
M.
Seishima
.
2018
.
Activation of NLRP3 signalling accelerates skin wound healing.
Exp. Dermatol.
27
:
80
86
.
79.
Wullaert
A.
,
M. C.
Bonnet
,
M.
Pasparakis
.
2011
.
NF-κB in the regulation of epithelial homeostasis and inflammation.
Cell Res.
21
:
146
158
.
80.
Kratofil
R. M.
,
P.
Kubes
,
J. F.
Deniset
.
2017
.
Monocyte conversion during inflammation and injury.
Arterioscler. Thromb. Vasc. Biol.
37
:
35
42
.
81.
Ford
J.
,
A.
Hughson
,
K.
Lim
,
S. V.
Bardina
,
W.
Lu
,
I. F.
Charo
,
J. K.
Lim
,
D. J.
Fowell
.
2019
.
CCL7 is a negative regulator of cutaneous inflammation following Leishmania major infection.
Front. Immunol.
9
:
3063
.
82.
Theocharidis
G.
,
B. E.
Thomas
,
D.
Sarkar
,
H. L.
Mumme
,
W. J. R.
Pilcher
,
B.
Dwivedi
,
T.
Sandoval-Schaefer
,
R. F.
Sîrbulescu
,
A.
Kafanas
,
I.
Mezghani
, et al
2022
.
Single cell transcriptomic landscape of diabetic foot ulcers.
Nat. Commun.
13
:
181
.
83.
Januszyk
M.
,
K.
Chen
,
D.
Henn
,
D. S.
Foster
,
M. R.
Borrelli
,
C. A.
Bonham
,
D.
Sivaraj
,
D.
Wagh
,
M. T.
Longaker
,
D. C.
Wan
,
G. C.
Gurtner
.
2020
.
Characterization of diabetic and non-diabetic foot ulcers using single-cell RNA-sequencing.
Micromachines (Basel)
11
:
815
.

The authors have no financial conflicts of interest.

Supplementary data