CD201+ fascia progenitors choreograph injury repair

Optimal tissue recovery and organismal survival are achieved by spatiotemporal tuning of tissue inflammation, contraction and scar formation1. Here we identify a multipotent fibroblast progenitor marked by CD201 expression in the fascia, the deepest connective tissue layer of the skin. Using skin injury models in mice, single-cell transcriptomics and genetic lineage tracing, ablation and gene deletion models, we demonstrate that CD201+ progenitors control the pace of wound healing by generating multiple specialized cell types, from proinflammatory fibroblasts to myofibroblasts, in a spatiotemporally tuned sequence. We identified retinoic acid and hypoxia signalling as the entry checkpoints into proinflammatory and myofibroblast states. Modulating CD201+ progenitor differentiation impaired the spatiotemporal appearances of fibroblasts and chronically delayed wound healing. The discovery of proinflammatory and myofibroblast progenitors and their differentiation pathways provide a new roadmap to understand and clinically treat impaired wound healing.

Optimal tissue recovery and organismal survival are achieved by spatiotemporal tuning of tissue inflammation, contraction and scar formation 1 .Here we identify a multipotent fibroblast progenitor marked by CD201 expression in the fascia, the deepest connective tissue layer of the skin.Using skin injury models in mice, single-cell transcriptomics and genetic lineage tracing, ablation and gene deletion models, we demonstrate that CD201 + progenitors control the pace of wound healing by generating multiple specialized cell types, from proinflammatory fibroblasts to myofibroblasts, in a spatiotemporally tuned sequence.We identified retinoic acid and hypoxia signalling as the entry checkpoints into proinflammatory and myofibroblast states.Modulating CD201 + progenitor differentiation impaired the spatiotemporal appearances of fibroblasts and chronically delayed wound healing.The discovery of proinflammatory and myofibroblast progenitors and their differentiation pathways provide a new roadmap to understand and clinically treat impaired wound healing.
Wound healing is a conserved process that requires tight spatiotemporal tuning of consecutive phases 1 .In the inflammatory phase, local signals trigger the innate immune system around the wound 2,3 .Timely wound healing requires the transition of the inflammatory phase into a proliferation phase in which specialized fibroblasts called myofibroblasts 4 produce the wound extracellular matrix (ECM).The final remodelling phase of wound healing is characterized by the lasting action of myofibroblasts in the restored tissue, in which further ECM production forms scar tissue.
The coordinated and timely initiation, persistence and resolution of each of the phases is key to effective healing.Protracted inflammatory phase or delayed proliferative phase can lead to chronic wounds observed in conditions such as ageing, diabetes, surgical or infectious wounds, vascular and heart diseases, or even cancer 1,3,5 .By contrast, a protracted remodelling phase leads to excessive tissue contraction and fibrosis, as observed in Dupuytren's contractures, hypertrophic scars and keloid lesions 1 .In sum, impaired wound healing imposes an enormous burden.Therefore, better understanding the choreography of tissue recovery is a matter of tremendous clinical and research interest.
We previously showed that wound fibroblasts in the dorsal skin originate from a single embryonic cell lineage that expresses the engrailed 1 (En1) gene during development, termed En1-past fibroblasts 6,7 (EPFs), and more recently we showed that EPFs from the fascia possess a higher scar-forming potential 8 .
Here we perform longitudinal stromal-enriched scRNA-seq combined with in vivo genetic lineage tracing, cell ablation and gene deletion approaches to generate a high-resolution pedigree study of all fibroblasts during skin wound healing.We identify a fibroblast progenitor in the subcutaneous fascia that controls the spatiotemporal features of wound healing by undergoing sequential differentiation into specialized fibroblastic cell types.

Fibroblast heterogeneity in wounds
To reveal the biological transition of fibroblasts during wound healing, we combined scRNA-seq with genetic lineage tracing of EPFs by crossing mice expressing Cre recombinase under the En1 promoter (En1 cre ) with the fluorescent Rosa26-mTmG (R26 mTmG ) reporter mouse line (Extended Data Fig. 1a).The resulting En1 cre R26 mTmG line tags EPFs permanently with the GFP.We then performed bilateral, full-thickness excisional wounds on the back skin and, at relevant days post-injury (dpi), we sequenced the mesenchymal-enriched fraction from the inflammatory (1 and 3 dpi), proliferation (5 and 7 dpi) and remodelling phases of wound healing (14 and 27 dpi) (Extended Data Fig. 1b).The initial sequencing approach was unbiased towards EPFs, and therefore included both GFP + and GFP − fibroblasts.
We detected seven distinct fibroblast clusters (Extended Data Fig. 1c).Clusters 0 to 2 were present in both uninjured and injured skin at all timepoints, whereas clusters 3 to 6 were only found in injured tissue (Extended Data Fig. 1d).To annotate previously known populations, we scored the expression of a curated set of published markers for resident papillary, reticular and fascia fibroblasts, and for classical myofibroblasts.Cluster 0 scored highest for the papillary profile, marked by Sparc 17 , Dcn and Lum.Cluster 2 scored highest for the reticular profile, marked by Cxcl12 17 , Cygb 17 , and Mgp.Cluster 1 scored the highest for the fascia profile, marked by Sca1 8, 24 (also known as Ly6a) and Plac8 17 , as well as Pi16 and Dpt, both previously reported as universal fibroblast markers 25 .We also identified Cd201 (also known as Procr), a marker of progenitor mesenchymal cells in spleen 26 and muscle 27 , as a specific fibroblast marker of the fascia (Extended Data Fig. 1e,f).
Next, we sought to trace the most likely origin of myofibroblasts by using partition-based graph abstraction (PAGA) for trajectory inference, which calculates the connections (for example, differentiation steps) between clusters.This initial unbiased analysis revealed a complex interconnection between all fibroblasts and corroborated previous observations that papillary fibroblasts do not contribute to the wound fibroblast pool 24 (Fig. 1a), and that the most likely fate of reticular fibroblasts is the SFRP2 + fibroblast (Fig. 1b).Furthermore, we identified the main myofibroblast source to be the proto-myofibroblast, which in turn branched off from the proinflammatory fibroblast.Notably, the proinflammatory fibroblast cluster derived exclusively from the fascia cluster, suggesting that proinflammatory, proto-myofibroblasts and myofibroblasts follow a sequential lineage trajectory that emerges from the fascia fibroblasts (Fig. 1b).As EPFs are the predominant lineage in wounds 6 , we aimed to place this lineage in the context of this new differentiation trajectory.We sorted the GFP + EPFs in our scRNA-seq dataset and reanalysed them with a higher clustering resolution (Extended Data Fig. 2a).EPFs were present across all homeostatic populations as well as in the proinflammatory, proto-myofibroblast and myofibroblast clusters (Extended Data Fig. 2b-d).
Our analysis also revealed an additional intermediate fascia EPF cluster that is present in injured samples but is absent from uninjured skin (Extended Data Fig. 2c,d).These injured fascia EPFs express the naive fascia markers Pi16 and Plac8, while simultaneously upregulating the proinflammatory fibroblast markers Pdpn, Ccl2 and Cxcl1 (Extended Data Fig. 2e).This reinforces the notion that the three injury-related cell clusters (proinflammatory, proto-myofibroblast and myofibroblast) are most probably derived from fascia-resident EPFs.
Based on our initial trajectory analysis and previous evidence of the contribution of fascia fibroblast into the myofibroblast pool 8 , we recalculated the sequential trajectory sprouting from the fascia cluster, passing through the proinflammatory and proto-myofibroblast states and ending in the myofibroblast state (Fig. 1c) for further study.

Genetic programmes along the trajectory
We next explored the genetic programmes of each fibroblast cluster from our trajectory.To do so, we mapped the transcription factor 'programme', their respective target genes ('regulon') and GO terms analysis of the same revealed discrete functions for each fibroblast cluster (Extended Data Fig. 3a).
The source of the trajectory, fascia fibroblasts, maintained a Fos, Jun and Egr1 programme, which controlled a 263-gene regulon that promotes cell survival and growth (Extended Data Fig. 3b), fitting its role as a surveillance system for injury repair.
The proinflammatory fibroblasts expressed a 116-gene regulon controlled by the Nfe2l2 and Bach1 programme, which directs oxygen sensing (Extended Data Fig. 3c), indicating that, besides immunomodulation, oxygen sensing is an important function acquired by fascia fibroblasts when becoming proinflammatory fibroblasts.
For the second step in the trajectory, the proto-myofibroblast programme comprised Stat3, hypoxia-inducible factor-1-alpha (Hif1a), En1, Ets2, Cebpb, Rora and Fosl1, which controls a 251-gene regulon involved in migration, oxygen sensing, adhesion and collagen production (Extended Data Fig. 3d).These programme changes reveal a progressive phenotype switch from immunomodulatory activities from the proinflammatory fibroblasts towards matrix-oriented functions, mimicking the transition from the inflammatory to proliferation phases of wound healing.
Meanwhile in the last step of the trajectory, the myofibroblast maintained a Tcf4, Runx1 and Runx2 programme that controlled a 21 genes regulon that directs terminal differentiation processes (Extended Data Fig. 3e).Together, our inferred trajectory and discrete phenotypes indicate that fascia fibroblast differentiation into proinflammatory, proto-myofibroblasts and myofibroblasts mirror the transition between the different wound healing phases (Extended Data Fig. 3f).

Spatiotemporal differentiation of CD201 + cells
We then sought to frame the cell differentiation steps of our trajectory within the wound healing phases in vivo.First, we confirmed the temporal expression patterns of our cluster-specific markers (Extended Data Fig. 4a) in mouse skin wounds at 3 and 7 dpi: two critical timepoints marking the transition from the inflammatory to proliferation phases and from the proliferation to remodelling phases of wound healing, respectively.PDPN + proinflammatory and phosphorylated activated STAT3 + (pSTAT3 + ) proto-myofibroblasts peaked during the first transition at 3 dpi and declined at later timepoints (Extended Data Fig. 4b,c), whereas RUNX2 + myofibroblasts were prominent at the second transition into the remodelling phases of wound healing at 7 dpi (Extended Data Fig. 4d).
We next analysed the spatial distribution of proinflammatory, proto-myofibroblast and myofibroblast clusters in three distinct wound compartments (Extended Data Fig. 4e).At 3 dpi, PDPN + proinflammatory fibroblasts were confined to the wound bed and upper wounds (Extended Data Fig. 4f), whereas pSTAT3 + proto-myofibroblasts showed no preference among these areas (Extended Data Fig. 4g).By contrast, at 7 dpi, RUNX2 + myofibroblasts clearly localized in upper wound compartments, with minimal presence in deep wound areas (Extended Data Fig. 4h).This gradual spatial bottom-up differentiation of fibroblasts Article towards upper wound compartments indicates a spatial distribution of guiding signals promoting the sequential differentiation of fascia fibroblasts (Extended Data Fig. 4i).
To validate that fascia progenitors differentiate into all three distinct fibroblast states in vivo, we performed genetic lineage tracing using an inducible CreER recombinase mouse line driven by the expression of the fascia marker CD201 (Cd201 creER ).In our scRNA-seq dataset, Cd201 expression enriched the fascia fraction 1.7-fold versus the dermal fraction in relation to the pan-fibroblast marker Pdgfra (60% fascia and 39% dermal), whereas Sca1 showed a lower enrichment of 1.24-fold (44.6% fascia and 55.4% dermal) (Extended Data Fig. 5a), indicating that Cd201 is a better candidate to specifically target fascia fibroblasts.We then crossed the Cd201 creER line with the fluorescent reporter line R26 Ai14 , in which tamoxifen administration induces the enduring expression of tdTomato in cells that express CD201 at the time of exposure in the Cd201 creER R26 Ai14 offspring (Fig. 1d).
In uninjured skin, our system almost exclusively labelled fibroblasts from the fascia compartment and not endothelial, epithelial or immune cells (Fig. 1e and Extended Data Fig. 5b-d).Similarly, tdTomato + cells in 7 dpi wounds were predominantly fibroblasts and not endothelial, epithelial, or immune cells (Extended Data Fig. 5b,c).Furthermore, tdTomato + cells contributed to 83% of all αSMA + myofibroblasts (Fig. 1e).Although CD201 has been reported as a marker of haematopoietic stem cells 30 , we observed no tdTomato expression on immune cells in the bone marrow or skin wounds (Extended Data Fig. 5b,c), discarding their potential contribution to the tdTomato + pool in skin wounds.We also observed that tdTomato + fibroblasts in the uninjured fascia spread along the connective tissue and were not restricted to the perivascular space (Fig. 1f).
Next, we tested whether CD201 + fascia fibroblasts differentiate into proinflammatory, proto-myofibroblasts and myofibroblasts upon injury.During the progression from inflammatory-proliferation to proliferation-remodelling transitions, CD201 + progenitors consistently downregulated proinflammatory and proto-myofibroblast markers: 71.5 down to 46.4% for PDPN and 88.2 down to 38.4% for pSTAT3 from 3 to 7 dpi.By contrast, the expression of the mature myofibroblast marker RUNX2 increased up to 73.5% in the CD201 fibroblast lineage at 7 dpi (Fig. 1g and Extended Data Fig. 6a).The temporal shift in marker gene expression of CD201 + fascia fibroblasts indicates that their transition into proinflammatory and to myofibroblasts occurs in synchrony to the progression of the wound healing phases.
To validate that the CD201 + fibroblast differentiation into myofibroblasts is preceded by the intermediate proinflammatory state, we generated a complementary system in which the expression of CreER recombinase is controlled by the proinflammatory marker Pdpn (Pdpn creER ).These mice were crossed with the R26 mTmG reporter mouse line, enabling lineage tracing from the proinflammatory state onwards (Fig. 1h).In contrast to our Cd201 creER fascia tracing system, Pdpn creER R26 mTmG mice show no GFP expression in any skin fibroblast under homeostatic conditions (Fig. 1i).At 7 dpi, αSMA + GFP + (proinflammatory-derived) myofibroblasts were detected abundantly in upper wounds (Fig. 1i).Furthermore, co-expression of PDPN, pSTAT3 and RUNX2 markers in GFP + cells at relevant times after injury further indicate that proinflammatory fibroblasts follow the same sequential trajectory into mature myofibroblasts (Extended Data Fig. 6b).Given that PDPN marks lymphatic vessels 31 , we determined the relative contribution of LYVE1 + lymphatic cells to the GFP + pool.We observed that only one-quarter of all GFP + cells were lymphatic cells and only one-third of all lymphatic cells were tagged with GFP + (Extended Data Fig. 6c), making most of the labelled cells descendants of the proinflammatory fibroblasts.These results confirm that proinflammatory fibroblasts transition into myofibroblasts during skin wound healing and that our Cd201 creER and Pdpn CreER transgenic systems faithfully trace their conversion steps.
In sum, in silico trajectory inference, spatiotemporal expression patterns and complementary genetic lineage tracing methods, detail a spatiotemporal coordination of fibroblast differentiation taking place during wound healing, which arises from CD201 + progenitors differentiating into proinflammatory fibroblasts within wound beds that finally mature into myofibroblasts in the upper wound region.

Targeting the trajectory
We next tested the relevance of the multi-step differentiation of fascia into myofibroblast for wound healing.As myofibroblasts are contractile cells, we first correlated tissue changes to fascia fibroblast differentiation in a mouse back-skin fascia ex vivo model (Fig. 2a).During culture and similar to mammalian wounds, the translucent fascia explants gradually contracted down to 50% of their original area until they formed an opaque sphere of scar-like tissue (Fig. 2b).
We then performed genetic cell ablation of PDPN + proinflammatory fibroblasts in vivo by crossing our Pdpn creER mice with the R26 DTA line, in which Cre-mediated recombination enables expression of the diphtheria toxin protein, leading to cell death.Thus, upon exposure to tamoxifen in Pdpn creER R26 DTA mice, PDPN + proinflammatory fibroblasts die before transiting into contractile proto-and myofibroblasts.We used the splinted wound model, which better mimics the human healing process by mitigating the action of other contractile tissues such as muscle.Splinted wounds in control mice (Pdpn creER R26 WT and Pdpn creER R26 mTmG ) showed a clear wound closure between 7 and 14 dpi, whereas the ablation of proinflammatory fibroblasts caused a significant delay in wound closure and contraction in Pdpn creER R26 DTA mice, whose wounds did not fully close by 14 dpi (Fig. 2e).In the experimental group, 7 dpi wounds showed a minimal wound bed with poor matrix content compared with control mice (Fig. 2f), indicative of a null wound healing progression.Notably, ablation of proinflammatory fibroblasts caused a significant decrease in the number of contractile αSMA + myofibroblasts in 7 dpi wounds (Fig. 2f).Together, our results indicate that fascia-derived proinflammatory fibroblasts have a pivotal role in the proper progression of the wound healing phases.

Retinoic acid gates the first transition step
The dynamic action of biochemical and biomechanical signals ensues the progression of the wound healing phases.To explore the connections between CD201 + progenitor differentiation and relevant signals, we scored the expression of genes from seven signalling pathways in our scRNA-seq dataset throughout the fascia-to-myofibroblast trajectory.Classical pathways associated with myofibroblast differentiation 4 -and thus, the proliferation and remodelling phases-such as TGFβ, Wnt and ECM mechanotransduction peaked at the myofibroblast state.Inflammatory-related signals, such as chemokines peaked at the proinflammatory fibroblast state (Fig. 3a).
The spatial link between RA synthesis and degradation, and fibroblast differentiation suggests that a RA gradient, originating within the wound bed and decreasing in the upper wound regions, enables the differentiation of progenitor fibroblasts to myofibroblasts.To test this, we either added exogenous RA or increased the endogenous RA levels using a CYP26B1 inhibitor in fascia explants.Both treatments prevented tissue contraction (Fig. 3g), indicating that RA is actively metabolized in culture and that high RA concentrations prevent fascia Article tissue contraction.Similarly, agonists for each of the three RA receptors phenocopied the contraction obstruction, and a pan-RA receptor antagonist significantly increased the overall tissue contraction (Extended Data Fig. 7a), showing that the RA signalling acts through its canonical transcriptional activity.
In fascia explants, treatment with an agonist of RARγ, the most highly expressed receptor in fibroblasts (Fig. 3b), increased the number of PDPN + proinflammatory fibroblasts and CCL2 expression, and significantly decreased the number of RUNX2 + and αSMA + myofibroblasts after 6 days of culture (Fig. 3h), indicating that RA overactivation promotes the proinflammatory state and prevents their conversion to myofibroblasts, leading to impaired tissue contraction.Notably, the pan-RAR antagonist did not increase myofibroblast numbers (Fig. 3h), indicating that other signals cause the transition from proinflammatory fibroblast to myofibroblast.
To explore the potential role of RA overactivation in sustaining the proinflammatory state, and thus the inflammatory phase, we treated skin wounds with RARγ agonist.Although subcutaneous exposure to the agonist had no significant effect on the wound closure rate (Extended Data Fig. 7b), treated wounds had a stronger and protracted inflammatory phase, revealed by a significant increase of monocyte and macrophage numbers in wounds at 3 and 7 dpi (Fig. 3i and Extended

Article
Data Fig. 7c,d).Conversely, treated wounds had significantly fewer αSMA + myofibroblasts associated with smaller scars at 14 dpi (Fig. 3j).These observations confirm that in vivo RA activity has a role in the proinflammatory fibroblast state and its overactivation limits the subsequent transition to myofibroblasts, reducing scar formation as result.
To expand on the effects of overactivation of RA signalling observed in vivo, namely monocyte recruitment and myofibroblast differentiation prevention, we cultured CD201 + fascia fibroblasts in inflammatory or proliferation phase-simulating media in combination with exogenous RA (Extended Data Fig. 7e).Inflammatory-simulating medium increased the expression of Ccl2 and Cxcl1 chemokine genes, whereas the proliferation-simulating medium only increased the expression of Acta2 (Extended Data Fig. 7e).Addition of RA to the inflammatory medium further increased Pdpn and Ccl2 expression and decreased Acta2 and Cxcl1 expression (Extended Data Fig. 7f).Addition of RA to proliferation medium caused a concentration-dependent decrease in Acta2 expression (Extended Data Fig. 7g), in line with the effects limiting myofibroblast differentiation observed in vivo and ex vivo.
Together, our results indicate that RA has a supporting role during the inflammatory phase of wound healing by promoting a monocyte and macrophage recruiting phenotype in the fascia-derived proinflammatory fibroblasts via expression of the monocyte-chemoattractant CCL2.Even though overactivation of the pathway effectively limits the amount of myofibroblasts, downregulation of the RA signal is insufficient to trigger the transition to the contractile states, and thus other signals are required.

Hypoxia gates the second transition step
Previous analyses highlighted the relevance of hypoxia during the transition from proinflammatory fibroblast to proto-myofibroblast, preceding classical myofibroblast-inductive signals such as TGFβ and ECM mechanotransduction mediated by YAP-TAZ 4 (Fig. 3a and Extended Data Fig. 1g).Notably, HIF1α-the master regulator of hypoxia-is part of the transcriptional programme of proto-myofibroblasts (Extended Data Fig. 3d).Furthermore, Hif1a transcriptional activity directly correlates with its expression in proto-and myofibroblasts (Fig. 4a) and HIF1α + cells were indeed detected within CD201 + -and PDPN + -derived lineages as well as in fascia explants at times when proinflammatory fibroblasts transition to proto-myofibroblasts (Fig. 4b).
Inhibition of HIF1α activity in explants impaired tissue contraction (Extended Data Fig. 8a), increased numbers of PDPN + proinflammatory fibroblasts and decreased numbers of αSMA + myofibroblasts (Extended Data Fig. 8b).Notably, HIF1α inhibition prevented tissue contraction in the presence of either TGFβ or YAP-TAZ activators; comparatively, exogenous RA was bypassed by the action of the stronger TGFβ signal (Extended Data Fig. 8c,d).Furthermore, HIF1α inhibition decreased the activation (via phosphorylation) of SMAD2, indicating a minimal TGFβ signal activity upon treatment (Extended Data Fig. 8e).To expand on the molecular connections between hypoxia and TGFβ pathways, we examined the correlation between gene expression of the HIF1α regulon and the inferred hypoxia signal activity in our scRNA-seq dataset (Extended Data Fig. 8f).Directly regulated genes included external modulators of TGFβ (Lgals1 32 and Serpine1 33 ), as well as internal modulators such as Eno1 34 , Ldha 35 , Pkm 36 , Ddit4 37 and En1 38 , which were all enriched in the proto-and myofibroblast clusters (Extended Data Fig. 8g,h).These data indicate that HIF1α acts as an upstream regulator of both YAP-TAZ mechanotransduction and TGFβ pathways.
In vivo chemical inhibition of HIF1α significantly delayed wound closure (Extended Data Fig. 8i) and these wounds had a marginal wound bed, similar to inflammatory phase wounds, whereas control wounds were completely closed and exhibited rich wound beds (Extended Data Fig. 8j).Notably, 3 dpi treated wounds showed similar leukocyte infiltration to control wounds (Extended Data Fig. 8k), confirming that HIF1α has a role only during the later contraction-related phases of wound healing.Furthermore, HIF1α inhibition in Pdpn creER R26 mTmG wounds significantly decreased proinflammatory-derived proto-myofibroblast (GFP + pSTAT3 + ) and myofibroblast (GFP + RUNX2 + ) numbers at 7 dpi (Fig. 4c), confirming that HIF1α activity is necessary for the transition from proinflammatory fibroblast to proto-and myofibroblast states, which in turn induces wound closure, tissue contraction and scar formation.
To unambiguously demonstrate the role of HIF1α during fasciato-myofibroblast differentiation and its effect in the progression to contraction phases of wound healing, we made mice deficient in HIF1α in fascia or proinflammatory fibroblasts by crossing Cd201 creER and Pdpn creER mice to a line harbouring loxP sites flanking exon 2 of the HIF1α gene (Hif1a flox ), which results in a null allele in response to tamoxifen-induced recombination 39 .Wounds on mice with single allele deletion in CD201 + fascia progenitors (Cd201 creER Hif1a fl/+ ) showed equivalent wound closure dynamics to wild-type mice, in which 91.7% closure of the initial wound size takes place around 9 dpi (Fig. 4d).However, bi-allele deletion of Hif1a in CD201 + fascia fibroblasts (Cd201 creER Hif1a fl/fl ) resulted in a significant delay in wound healing (Fig. 4d), with a closure rate of 72.7% at 9 dpi.We also observed impaired wound closure in mice that were deficient in HIF1α in PDPN + proinflammatory fibroblasts (Pdpn creER Hif1a fl/fl ), which had a wound closure rate of 62.5% at 9 dpi (Fig. 4d).At the same timepoint, wounds of mice deficient in HIF1α in CD201 + fascia and PDPN + proinflammatory fibroblasts exhibited an immature wound bed (Fig. 4e), similar to wounds in the inflammatory phase.In both systems, wound beds did not contract and develop (Fig. 4e), revealing a limited action of contractile proto-myofibroblast and myofibroblast states, which severely impairs the progression to the final phases of wound healing.
We previously showed that fascia fibroblasts transfer pre-made ECM into wounds, which contracts and contributes to the scar matrix 8 .To understand how the transition from immunomodulatory (proinflammatory fibroblast) to contractile and matrix-modifying states (proto-myofibroblast and myofibroblast) affects ECM dynamics in vivo, we labelled the fascia ECM prior to injury (Extended Data Fig. 9a).At 3 dpi, transferred pre-made fascia ECM covered the wound bed, where PDPN + proinflammatory fibroblasts were preferably allocated on thinner, looser fibres, whereas pSTAT3 + proto-myofibroblasts were associated with thicker, woven fibre bundles (Extended Data Fig. 9b).These associations support the notion that differentiation of CD201 + progenitors into ECM-modifying states (proto-myofibroblast and myofibroblast) mediates matrix malleability and tissue contraction.To test this, we blocked the differentiation into proto-myofibroblast and myofibroblasts with HIF1α inhibitor in mice with labelled fascia ECM.Again, control 3 dpi wounds showed fibre bundles with a condensed organization, whereas treated mice retained slackened and weakly organized fibres (Fig. 4f).We assessed the contraction of the traced matrix by measuring its porosity and complexity using fractal analysis.This analysis revealed a significant reduction in the complexity and an increase in porosity in the HIF1α-inhibition group compared with controls (Fig. 4g), indicative of a relaxed fibre organization compared wuth the more condensed conformation in controls.
Together, our comprehensive scRNA-seq, ex vivo and in vivo chemical inhibition experiments, genetic lineage tracing, gene deletion and matrix fate map studies showed that in mouse, HIF1α gates the second differentiation step from proinflammatory fibroblast to proto-myofibroblast and myofibroblast that controls the timing of the progression of the inflammatory and contraction (proliferation and remodelling) phases of wound healing.

Trajectory in human skin pathologies
To explore whether fascia progenitor differentiation into myofibroblasts is conserved in humans, we merged scRNA-seq datasets of PDGFRA + fibroblasts from fibrotic conditions hallmarked by the action of myofibroblasts, such as Keloid lesions 18 , hypertrophic scars 40 and scleroderma 41 .We also included psoriatic skin, an autoimmune disease in which fibroblasts acquire a proinflammatory phenotype 42 (Extended Data Fig. 10a,b).We then used transfer learning to annotate the human fibroblast supercluster using our mouse dataset as reference.As well as the three homeostatic populations, we detected the proinflammatory fibroblast-like, proto-myofibroblast-like and myofibroblast-like populations in the human skin pathologies (Fig. 5a).Comparing fibroblast compositions between skin pathologies, we noted that fascia fibroblasts were four times more abundant in hypertrophic scars, and papillary fibroblasts were twice as common in scleroderma samples than other conditions.Proinflammatory fibroblasts

Article
and myofibroblasts were present in all conditions.Proinflammatory fibroblasts were more prominent in keloids and scleroderma, whereas myofibroblasts were more abundant in psoriatic skin.Notably, proto-myofibroblasts were only present in hypertrophic scars and psoriatic skin samples (Extended Data Fig. 10c).We next explored the conservation of the population-specific markers by comparing the expression of the top markers of the human clusters to their mouse counterpart.The most conserved marker for papillary fibroblasts across mouse and human was CD9, and APOE, CFH and IGFBP7 marked the reticular population in both species (Extended Data Fig. 10d).Conserved markers for proinflammatory fibroblasts included CCL2, CXCL2 and NFKBIA, and COL12A1 marked the proto-myofibroblast cluster (Extended Data Fig. 10d,e).The fascia clusters had PI16 as the most specific marker in both human and mouse datasets (Extended Data Fig. 10d,e).Immunolabelling of PI16 in human keloid and healthy dermis and hypodermis corroborated that PI16 protein expression was restricted to fascia connective tissue in the hypodermis and was absent in dermis or keloid lesions (Fig. 5b).By contrast, myofibroblasts showed the highest conservation of marker genes, including POSTN, several collagen genes and RUNX2 (Extended Data Fig. 10d,e).As RUNX2 is not a classical myofibroblast marker, we confirmed its presence within keloid lesions in αSMA + cells and its absence from healthy tissue (Fig. 5b).These results strongly indicate that homeostatic and wound fibroblasts from mouse wound healing strongly resemble human fibroblasts.
We next analysed the trajectories between the human fibroblast clusters.Given that proto-myofibroblasts were present in only half of the conditions, we omitted them to obtain a general trajectory roadmap applicable for all pathologies.In this initial trajectory analysis, we observed stronger connectivities within all homeostatic populations and into the proinflammatory and myofibroblast populations as terminal states (Fig. 5c), indicating variable contributions from the homeostatic populations to the wound fibroblast pools.
We then assessed the conservation of the fascia fibroblast differentiation trajectory to myofibroblasts in each disease.All pathologies showed an initial transition to proinflammatory fibroblasts before transiting to myofibroblasts; scleroderma was the sole exception that showed the proinflammatory fibroblast as the terminal state (Fig. 5d).These analyses indicate that human fascia fibroblasts transition to myofibroblasts through a proinflammatory state, as in mouse wound healing, and that this differentiation process can be disturbed in skin pathologies.
We then investigated the conservation of the two signalling gates, RA and hypoxia.Comparative analysis of the chronology of these signalling pathways with the human datasets showed that RA activity peaked in the proinflammatory state across all diseases with the exception of psoriatic skin, in which RA activity peaked in the fascia fibroblasts (Fig. 5e).Conversely, hypoxia signal activity closely correlated with the TGFβ signalling pathway along the transition from the proinflammatory to myofibroblast state in all human pathologies (Fig. 5e).Furthermore, human fascia explant contraction halts when exposed to exogenous RA or HIF1α inhibitor (Fig. 5f), indicating a conserved role of these two signalling gates during fascia fibroblast differentiation to myofibroblasts.
Together, our analyses indicate that similar to mouse, human fascia fibroblasts undergo a differentiation trajectory to proinflammatory fibroblasts, partially supported by RA signalling, and to later contractile myofibroblasts, mediated by HIF1α-mediated hypoxia signalling.This differentiation trajectory has a major impact in the progression from the inflammatory phase to the tissue contraction phases (proliferation and remodelling) of wound healing.

Discussion
This study revealed a conserved differentiation process with clinical implications for injury repair.At the root of this process are fascia progenitors that give rise to specialized fibroblast types and two regulatory gates needed for the timely progression of the wound healing phases.
The first transition of fascia progenitors to proinflammatory fibroblasts has an important effect on the initial inflammatory phase.Our data revealed that RA supports this phase by promoting a monocyte recruitment phenotype and by limiting their subsequent conversion into myofibroblasts, accounting for the observed anti-fibrotic effects of exogenous RA 11 .The second differentiation step from proinflammatory fibroblast to proto-myofibroblast, marked by the acquisition of a contractile phenotype, enables the transition from the inflammatory to the proliferation and remodelling phases of wound healing.Previous studies highlight the importance of HIF1α in fibroblasts for timely wound healing 43 .Besides its classical roles of reactive oxygen species clearance and angiogenesis regulation 44 , HIF1α has previously been only indirectly linked to myofibroblast transition via its synergy with the TGFβ signalling in vitro 45 ; our data provide direct evidence for upstream regulation of the TGFβ pathway by HIF1α.
An important caveat to consider in this present study, is the limited number of biological replicates in our in vivo experiments.Nonetheless, multiple lines of evidence, including cross-species in silico analyses, novel ex vivo systems and the use of several complementary chemical and genetic mouse models, point towards the conclusion that the conserved differentiation process of fascia fibroblasts to myofibroblasts orchestrates the progression of the wound healing phases during skin repair.
Despite early reports showing a prevalence of myofibroblasts in upper wound regions 46 , little attention has been paid to the spatial distribution of specialized fibroblast populations within wounds.Our results draw an important distinction between wound compartments in several ways.CD201 + progenitors are found within the very deep layers of the skin, where they follow spatial bottom-up cues during their differentiation road upon injury.The first spatial cue occurs in the wound bed where progenitors differentiate into proinflammatory fibroblasts and RA signalling is at its highest.In parallel, the second spatial cue appears within the oxygen deprived wound bed in the form of an intrinsic hypoxia response in proinflammatory fibroblasts.The combination of these two spatial cues promotes the controlled transition into proto-myofibroblasts while limiting their terminal differentiation into myofibroblasts due to the activity of RA.The final spatial cue occurs in upper wound regions during the transition to the proliferation phase, where RA is degraded and terminal differentiation into myofibroblasts is granted.Thus, our data indicate that the interplay between RA gradients in hypoxic niches regulate fibroblast differentiation in response to injury.
These results have important clinical implications, as alterations in this differentiation process could result in pathologies such as impaired wound healing and fibrotic scarring.In particular, the manipulation of individual steps along the trajectory opens new roads to modulate different wound fibroblasts to obtain disease-tailored therapies.Treatments of wounds (such as surgical and keloid wounds, and burns) could be modulated to restrain myofibroblasts to prevent excessive scarring and contraction without compromising the proinflammatory role needed during the inflammatory phase, and treatments tailored for modulating the early progenitors would be beneficial for inflammatory skin diseases.

Online content
Any methods, additional references, Nature Portfolio reporting summaries, source data, extended data, supplementary information, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at https://doi.org/10.1038/s41586-023-06725-x.
Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made.The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material.If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

In vivo wounds, tamoxifen induction and treatments
Male and female mice, equally distributed between time groups, were used for descriptive experiments.For functional studies, littermates or age-matched mice were randomly assigned into the different experimental groups.To decrease variability, only males were used for the RARγ agonist and HIF1α inhibition treatments (7 dpi), and for the Hif1α genetic deletion experiment.Females and males, equally distributed between treatment groups, were used for the PDPN + cell ablation studies and for the 3 dpi HIF1α inhibition experiment.Back-skin full-thickness wounds on adult (8-to 12-week-old) mice were performed as previously 8 .In brief, bilateral 5-mm circular wounds were performed by piercing through a skin fold along the midline of the back using biopsy punches (Stiefel 10005).For the splinted wound model 47 , immediately after making the excisional wounds, circular silicon rings were glued (Gracebio 664581) and sutured around the open wounds.For the genetical lineage tracing of CD201 + fascia fibroblasts and Hif1αa gene deletion, Cd201 creER R26 Ai14 , Cd201 creER Hif1a fl/+ , Cd201 creER Hif1a fl/fl , and Pdpn creER Hif1a fl/fl mice received three intraperitoneal injections of 2 mg (Z)−4-Hydroxytamoxifen (Sigma Aldrich H7904) diluted in 0.2 ml miglyol-812 (Caesar & Lorentz GMBH CSLO3274) two days and one day before wounding, and at the day of wounding.For the genetic lineage tracing and ablation of PDPN + proinflammatory fibroblasts, Pdpn creER R26 mTmG or Pdpn CreER R26 DTA mice received daily intraperitoneal injections of 1 mg (Z)−4-hydroxytamoxifen diluted in 0.1 ml miglyol-812 from day four before wounding to the day of injury.For fascia matrix fate mapping 8,48 , 40 μl of 2 mg ml −1 Pacific Blue succinimidyl Ester (Invitrogen P10163) with 0.1 M sodium bicarbonate in saline were subcutaneously injected on the mouse back three-and one day prior to wounding to label the fascia matrix.HIF1α inhibitor (echinomycin, Sigma Aldrich sml0477), at a concentration of 100 μg kg −1 body weight, RARγ agonist (BMS 961, Tocris 3410), at a concentration of 100 μg kg −1 body weight, or vehicle control dimethyl sulfoxide (DMSO, Sigma Aldrich D5879) were subcutaneously injected near and below the wounds immediately after injury and every other day to have a more direct effect on fibroblasts.

Image Analysis
All image analysis was performed using Fiji (v.1.53c).For assessment of tissue contraction in vitro and wound closure in vivo, the explant or wound area from different timepoints was selected using the 'selection brush' tool and measured using the 'measure' function.Area values were normalized in percentage to the initial (day 0, A 0 ) as 100%.In vitro final contraction (C) was defined as the difference between the highest area value of the two initial days after culture (day 0-2, A max ) minus the lowest area value from the last two days of culture (day 4-6 A min ); thus, C = A max − A min .Wound length was determined as the mean length from three individual measurements in each section by tracing lines covering along (1) the new epidermis, (2) distance between flanking hair follicles, and (3) between the branched Panniculus carnosus muscle using the line measurement tool.For cell marker-positive quantifications, individual channels were pre-processed using the functions 'substract background' (rolling=20px), 'median' (radius=1), 'unsharp mask' (radius=1, mask=0.60),'despeckle', 'enhance contrast' (saturated=0.1 normalize), and 'auto threshold' methods were empirically selected for individual markers: Moments (for DAPI, pSTAT3, YAP1, pSMAD2, CD45, and GFP signal), Intermodes (RUNX2, tdTomato, PDGFRα, and LYVE1), Yen (CCL2), IsoData (PDPN), or RenyiEntropy (αSMA).The wound region of interest (ROI) was manually selected using the 'selection brush' tool taking as reference the breached P. carnosus below the dermis and the flanking hair follicles in the epidermis above.Total wound cell numbers (within the ROI) from binary images were measured using the DAPI channel and using the functions 'fill holes', 'watershed', and 'analyze particles' (size=70-1400px) in that order.Marker-positive cell numbers were measured by using the mask of the marker with the DAPI channel using the 'image calculator_AND', 'dilate' twice, 'close', 'fill holes', 'watershed', and 'analyze particles' (size=70-1400px) functions.For double marker-positive cell calculations, an initial mask from two marker channels was created and processed as before.For wound region quantifications, 3 ROIs were manually selected with the 'selection brush' tool.The wound lateral limits were delimited as stated before, the upper wound region extended vertically 200 microns from the wound surface (at 3 dpi) or the wound epidermis (at 7 dpi).The fascia was measured as the area extending below the breached P. carnosus muscle down to the dorsal muscles, while the wound core covered the area between the wound fascia and the upper wound.Non-mesenchymal areas (such as hair follicles and wound epidermis) were deleted from the ROIs.Marker channels were pre-processed as before and the mean grey value (MGV) for each ROI was then quantified with the 'measure' function.For visual aid, marker signal intensity was represented with the 'fire' or '16-color' pseudo-colours.ECM area density was determined by converting Masson's trichrome-stained RGB images into CMYK using the RGB-to-CMYK plugin.Then, the cyan channel (ECM) was extracted and the MGV measured before being pre-processed as above.ECM density was defined as the ECM (cyan) area in the wound covering the total wound area.ECM density and ECM MGV in wounds were normalized to the healthy adjacent dermis to avoid batch effects between stained sections.Fractal analysis was performed as before 7 .In brief, the cyan channel was pre-processed as before, and binary images were generated using 'auto threshold' ('Default' method) function.Binary images were then analysed with the FracLac plugin using the same settings as before.All quantification plots were generated using the Python-based packages Matplotlib and Seaborn.

Dropseq libraries preparation and sequencing
Stromal-enriched cell suspensions were separated using the Dropseq 49 microfluidic-based method following the McCaroll lab protocol with a few adaptations.In brief, using a microfluidic PDMS device (Nanoshift), cells were encapsulated in droplets with barcoded beads (120 per μl, ChemGenes Corporation, Wilmington, MA) at a rate of 4,000 μl h −1 .Droplet emulsions were collected for 15 min each prior to droplet breakage with perfluoro-octanol (Sigma-Aldrich).After breakage, beads were collected, and the hybridized mRNA transcripts were reverse transcribed (Maxima RT, Thermo Fisher).Unused primers were removed by the addition of exonuclease I (New England Biolabs).Beads were then washed, counted, and aliquoted for pre-amplification (2,000 beads per reaction, equal to ~100 cells per reaction) with 12 PCR cycles (with the recommended settings).PCR products were pooled and purified twice by 0.6× clean-up beads (CleanNA).Prior to tagmentation, cDNA samples were loaded on a DNA High Sensitivity Chip on the 2100 Bioanalyzer (Agilent) to ensure transcript integrity, purity, and amount.For each sample, 1 ng of pre-amplified cDNA, from an estimated 1,000 cells, was tagmented by Nextera XT (Illumina) with a custom P5 primer (Integrated DNA Technologies).Single-cell libraries were sequenced in a 100 bp paired-end run on an Illumina HiSeq4000 using 0.2 nM denatured sample and 5% PhiX spike-in.For priming of the first read, 0.5 μM Read1CustSeqB (primer sequence: 5′-GCCTGTCCGCGGAAGCAGTGGTATCAACGCAGAGTAC-3′) was used.STAR (version 2.5.2a) was used for mapping the reads and to align them to the mm10 genome reference (provided by Drop-seq group, GSE63269) that was tailored to include the eGFP cDNA transcript.

scRNA-seq data analysis
All the analyses were performed using the phyton toolkit Scanpy 50 and complementary tools under its ecosystem.Matrices of individual samples were concatenated, quality checked by filtering cells with fewer than 100 genes.Genes that appeared in less than 50 cells and cells containing less than 0.5% of the total genes in the dataset were also filtered out.To decrease batch-and cell cycle effects, combat and cell cycle regression algorithms were used as recommended by the developers.The UMAP algorithm was used as the preferred dimensional reduction method.
Human skin diseases datasets were obtained from the GEO repository (GSE163973, GSE162183, GSE156326 and GSE160536), merged, and processed as before.Annotation of human fibroblast clusters was performed using scArches 54 .In brief, after training on the reference dataset (mouse) the node weights were transferred to a new model.Human dataset was then tested on the trained model and cells were annotated accordingly.Parameters were initialized in accordance with the scArches optimization.
Mouse trajectories from all fibroblasts and from fascia-tomyofibroblast trajectory clusters were inferred by PAGA with RNA velocity-directed edges using the scvelo toolkit 55,56 .'Dynamic modeling' was used under standard settings to calculate velocities.Cells were arranged under the inferred trajectory using the 'velocity pseudotime' function rooting the trajectory origin in the most distant cells within the fascia cluster from the rest of the clusters.Human trajectories (unbiased and from fascia-to-myofibroblast trajectory clusters) were inferred using PAGA without RNA velocity due to the non-matched, monotemporal nature of the human samples.Arragment of the trajectories were performed then using diffusion pseudotime instead of velocity pseudotime.
Regulon analysis was performed by mapping the expression of 1623 mammalian transcription factors from the RIKEN database and statistically correlating (Pearson R) their expression with the fibroblastic state signatures.Programmes were selected from the list of transcription factors whose expression positively correlated with a signature that had validated downstream targets in the TRANSFAC and ENCODE databases.Individual programme scores were calculated according to the expression of the selected transcription factors using the score gene function.To identify the complete regulon, the combined list of downstream targets for each programme was mapped and correlated with the programme scores.Positively correlated genes, whose expression matched the upstream regulating programmes, were used for the regulon score using the gene score function and for the GO term overrepresentation test using the PANTHERdb 57 web tool.
Transcription factor activity was inferred using the DoRothEA toolkit 58 .Signalling pathway activity was predicted by mapping the expression of genes from a manually curated list for relevant factors of the TGFβ, Wnt, mechanotransduction, cell-cell adhesion, chemokine, and RA signalling pathways.Complementarily, hypoxia, Wnt, and TGFβ pathway activity was predicted using the PROGENY tool 58 , showing similar results, for TGFβ and Wnt signal activity, to our manual method.The HIF1α regulon positively regulated during hypoxia signal was identified by correlating the validated downstream targets of HIF1α from the ENCODE and DoRothEA databases to the hypoxia signal activity from PROGENY.

Statistics
Statistical analyses were performed using the Python toolkit Scipy.No methods were used to predetermine sample size.No blinding or randomization was used for data analysis.Two-tailed t-tests, one-way or two-way ANOVAs, and Pearson's R were implemented with standard settings.Significant differences were considered with P values below 0.05 and are indicated in the figures.In all plots, error bars depict the standard deviation.

Reporting summary
Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Fig. 1 |
Fig. 1 | CD201 + fascia fibroblasts differentiate into specialized states during wound healing.a, Uniform manifold approximation and projection (UMAP) analysis of all fibroblast clusters embedded with PAGA connectivities for trajectory inference.FB, fibroblast; MyoFB, myofibroblast; proinf., proinflammatory.b, Left, PAGA connectivity values with potential sources (rows) and fates (column) colour-coded as most (green), intermediate (yellow) and least probable (orange).Right, schematic of the potential trajectories highlighting the fascia-to-myofibroblast trajectory.c, UMAP analyses of fibroblasts from the fascia-to-myofibroblast trajectory, colour-coded for individual clusters (top) and velocity pseudotime score (bottom).d, Genetic lineage system to trace the fate of CD201 + fascia fibroblasts during wound healing.e, Representative histology of uninjured (left) Cd201 creER R26 Ai14 and 7 dpi wound (middle) immunolabelled for αSMA.Right, the contribution of CD201 + fibroblasts to the myofibroblast pool was determined for 6 wounds

Fig. 2 |Fig. 3 |
Fig. 2 | Targeting fascia differentiation impairs tissue contraction and wound closure.a, The fascia explant culture method.b, Representative bright-field images of fascia explants at indicated timepoints of culture (left) and area versus time contraction measurements (right).n = 5 from 3 biological replicates.c, Representative micrographs of fascia explants immunolabelled for PDPN, pSTAT3, RUNX2 or αSMA (left) and quantification of expressing cell fraction at 1, 3 and 6 days of culture normalized to day 1 values (right).n = 3 biological replicates.d, Representative images of fascia explants at 6 days of culture (left), contraction dynamics (top right) and total contraction values (bottom right) showing the effects of cell ablation (freeze-thawing) and different serum conditions.n = 3 (full serum and freeze-thawing) and 4 (heat-inactivated high and low serum) biological replicates.HI, heat-inactivated.Colours of the image outlines indicate treatments in the graphs.e, Representative photographs of splinted wounds on PDPN-related transgenic lines at indicated timepoints after injury (top), wound area quantification (bottom left) and contraction percentage at the indicated dpi (bottom right), showing that ablation of proinflammatory fibroblasts delays wound closure.n = 6 wounds from 3 biological replicates.WT, wild type.f, Masson's trichrome staining of control (top left) and ablated (bottom left) 7 dpi wounds.The dotted line delimits the wound region.Arrowheads indicate the original injury site.High-magnification images of control (top middle) and ablated (bottom middle) 7 dpi wounds immunolabelled for αSMA.Right, myofibroblast ratio in wounds.n = 4 (control) and 5 wounds from 3 biological replicates.Colours of the image outlines indicate treatments in the graph.Two-tailed t-tests.Scale bars: 2 mm (b,d), 500 μm (f, left) and 50 μm (c and f, middle).

6 Fig. 4 |
Fig. 4 | HIF1α mediates the transition to proto-myofibroblasts and tissue contraction.a, HIF1α activity versus expression correlations in the different fibroblast clusters from the scRNA-seq data.Pearson's R coefficient is shown.b, Representative micrographs of wounds in Cd201 creER R26 Ai14 and Pdpn creER R26 mTmG mice, and fascia explant at indicated timepoints after injury or culture showing the expression of HIF1α.c, Representative high-magnification images of control and HIF1α inhibitor-treated Pdpn creER R26 mTmG wounds showing the expression of pSTAT3 or RUNX2 in GFP-traced cells (left) and marker-expressing cell ratios (right).n = 6 wounds from 3 biological replicates.d, Representative photographs of control, fascia or proinflammatory HIF1α null wounds at 9 dpi (left).Wound area quantification over time (middle), and contraction at 7 dpi (top right) and at 9 dpi (bottom right).Colours of the image outlines indicate treatments in the graph.n = 6 wounds from 3 biological replicates for each genotype.e, Masson's trichrome staining (top) and wound maturation-related measurements (bottom) of 9 dpi wounds from indicated genotypes.n = 12 images from 3 biological replicates per each genotype.f,g, Fascia ECM fate mapping in 3 dpi wounds treated with HIF1α inhibitor.f, Representative low-magnification (left) and high-magnification (right) images of control or HIF1α inhibitor-treated wounds.g, Lacunarity and fractal dimension density plots to assess porosity (lacunarity) and shape complexity (fractal dimension) differences in labelled ECM from control or HIF1α inhibitor-treated wounds (left) and individual comparisons (right).n = 11 images from 3 biological replicates for each condition.The dotted line delimits the wound region.Arrowheads indicate the original injury site.All P values (except in a) were obtained from two-tailed t-tests.Scale bars: 500 μm (e and f, left) and 50 μm (b,c, and f, right).

Fig. 5 |
Fig. 5 | Fascia-to-myofibroblast trajectory in human skin pathologies.a, 'Transfer learning' for the mapping of human fibroblasts using the mouse atlas as a reference.UMAP representation of mouse (left) and human fibroblasts (right).b, Photograph (top left) and low-magnification image of trichrome-stained section (bottom left) of a human keloid lesion.Right, high-magnification micrographs of indicated regions of healthy and keloid tissue that were trichrome-stained or immunolabelled for the fascia marker PI16 or the myofibroblast markers αSMA and RUNX2.c, PAGA connectivities embedded onto the UMAP graph of the human fibroblast supercluster (left) and the diffusion pseudotime arrangement (right).d, UMAP analyses of subclusters involved in the fascia-to-myofibroblast trajectories in mouse wound healing and indicated human skin pathologies (top) and their

Extended Data Fig. 1 |
Longitudinal single-cell RNA sequencing of stromal cells during murine wound healing.a. Experimental approach for the Engrailed-1 lineage tracing coupled with stromal cell enrichment and Dropseq method for single-cell sequencing.b.UMAP representations of complete dataset color-coded for samples (time points and cell number, left) and Pdgfra expression (right) to identify fibroblasts.c.UMAP representation of fibroblast clusters (left) and the cluster composition of the different samples (right) revealing homeostatic (always present) and injury-specific populations (absent in uninjured sample).d.Fraction in % of each injury-specific population across samples (time-points).e. Ridge plots depicting profile scores from indicated fibroblast subtypes in each cluster.Gray lines indicate the overall mean for each profile score and the white lines indicate the cluster mean.Colour arrows indicate clusters with the highest profile for each fibroblast subtype.f.Gene expression dot-plots of top markers for homeostatic (left) and injury-related (right) clusters.g.Overrepresented GO terms in indicated injury cluster.Scatter plots show the fold enrichment (x-axis) of the term (y-axis), genes expressed/genes in term coverage in % (dot size), and the false discovery rate value (colour gradient).Extended Data Fig. 2 | Fascia-to-myofibroblast differentiation occurs within the En1-lineage.a. UMAP representation of fibroblast showing GFP expression (EPFs) that were sorted and subclustered.b.UMAPs showing the skin fibroblast profile scores (top) and their associated EPF subclusters (bottom).c.UMAP of EPFs subclusters.d-e.An intermediate state with both naïve fascia and proinflammatory profiles supports the transition sprouting from fascia fibroblasts.d.Sample composition of EPF subclusters from the fascia-to-myofibroblast trajectory showing that the intermediate state appears in early wounds.e.Expression dot-plot of naïve fascia and proinflammatory markers being present in the intermediate state.Extended Data Fig. 3 | State-specific regulons reveal a stepwise differentiation linked to the wound healing phases.a. Analysis pipeline to determine the collection of transcription factors (programme) and their associated regulons for each fibroblastic state.b-e.UMAPs of programme (left) and regulon scores (centre), and their overrepresented GO terms (right) for the fascia (b), proinflammatory (c), proto-myofibroblast (d), and myofibroblast states (e).f.Scheme of the fibroblastic states along the trajectory and their phenotype revealed from scRNAseq analysis.Extended Data Fig. 4 | Spatiotemporal patterns during the transition.a. UMAPs of fibroblast from trajectory, colour-coded for individual clusters (left) and marker genes used in this study (right).b-d.Representative immunostainings (left) for PDPN (a), pSTAT3 (b), or RUNX2 (c) from 3-(a-b) or 7-dpi wounds (c) and marker-expressing cell quantifications at different timepoints after injury (right).N = 16 (3 dpi), 12 (7 dpi), and 3 (14 dpi) images from PDPN-stained (a) sections.6 (3 dpi) and 4 (7 and 14 dpi) images from pSTAT3-stained (b) sections.3 (3 dpi) and 4 (7 dpi) images from RUNX2-stained (c) sections.All images were obtained from 3 biological replicates for each timepoint.High magnification images (bottom) from rectangle insert in lowpower images (top).Dotted line delimits the wound region.Arrowheads indicate the original injury site.e. Schematics of the wound regions analysed and the spatial cues that defined them.f-h."Fire" pseudo-colouring representations (top) to depict expression intensity of PDPN (e), pSTAT3 (f), or RUNX2 (g) within the different wound regions at indicated dpi.Mean Gray Value (MGV) expression intensity quantifications (bottom) at indicated dpi in the different wound regions.N = 12 (upper wound and wound core) and 13 (wound fascia) values from PDPN-stained (e) sections.7 (upper wound and wound core) and 8 (wound fascia) values from pSTAT3-stained (f) sections at 3 dpi.4 values from pSTAT3-stained (f) sections at 7 dpi.3 (upper wound) and 4 (wound core and wound fascia) values from RUNX2-stained (g) sections.All images were obtained from 3 biological replicates for each timepoint.Dotted line delimits each wound region.All p-values (p) indicated were obtained from two-tailed T-tests.Scale bars: 500 microns in low magnification (a-c and e-g), and 50 microns in high magnification micrographs (a-c).i. Spatial distribution of the fibroblast states in discrete wound regions.Extended Data Fig. 5 | CD201 CreER R26 Ai14 genetic lineage tracing system to fate map fascia fibroblasts.a. Homeostatic fibroblast fractions in fascia marker-expressing cells within the scRNAseq dataset.b.Flow cytometry strategy to identify proportions of different cell types labelled with TdTomato in the uninjured dermis and fascia, as well from wound and bone marrow from 7 dpi CD201 CreER R26 Ai14 mice.c.TdTomato-labelled fraction (left) from total fibroblasts (PDGFRα + ), epithelial (EPCAM + ), leukocytes (PTPRC + ), and endothelial cells (PECAM1 + ) present in dermis, fascia, wounds, or bone marrow.Insert: labelled fibroblast fraction in wildtype and CreER positive animals shows minimal leaking of the genetic system.All p-values (p) indicated were obtained from two-tailed T-tests.2-way ANOVA results from TdTomato-positive fractions comparison (right).N = 5 biological replicates.Representative micrographs immunolabeled for fibroblast (PDGFRα), keratinocytes (KRT14), leukocytes (PTPRC), and endothelial cells (PECAM1) markers on CD201 CreER R26 Ai14 skin sections validating the flow cytometry observations.Scale bars: 50 microns.Extended Data Fig. 6 | Expression of wound fibroblast markers on the CD201 CreER R26 Ai14 and Pdpn CreER R26 mTmG genetic lineage tracing systems.a. Representative low and high magnification (inserts) images showing the expression of PDPN (top), pSTAT3 (middle), and RUNX2 (bottom) in TdTomatolabelled fascia-derived cells on 3 (left) and 7 dpi (right) wounds.b.Representative low (left) and high magnification (right) images showing the expression of PDPN (top), pSTAT3 (middle), and RUNX2 (bottom) in GFP-labelled proinflammatoryderived cells on wounds at indicated timepoints.c.Representative low (left) and high magnification (middle) images showing the residual labelling with GFP in lymphatic vessels (LYVE1 + ) in wounds at indicated timepoints.Cell fractions of labelled-lymphatic endothelial cells from total traced and from total lymphatic vessels in wounds at indicated timepoints.N = 3 biological replicates.Dotted line delimits the wound region.Arrowheads indicate the original injury site.Scale bars: 500 microns in low magnification, and 50 microns in high magnification micrographs.Extended Data Fig. 7 | RA supports the proinflammatory state by inducing a monocyte recruitment phenotype.a. Representative brightfield images at day 6 of culture (left), explant area quantifications over time (middle), and overall contraction quantifications (right) from fascia explants treated with RAR-selective agonists and pan-antagonist.N = 4 (area vs time) and 3 (total contraction) biological replicates.Two-tailed T-tests.Scale bars: 2 mm.b.Representative photos of control or RARγ agonist-treated wounds at 6 dpi (left), wound closure measurements (middle), and total wound contraction by 6 dpi (right).N = 6 wounds from 3 biological replicates.c.Flow cytometry strategy to identify treatment-induced changes in the recruitment of different immune cell types in wounds.d.Total cell type fractions at 3 (top) and 7 dpi (bottom) of general leukocytes (PTPRC + ), monocytes/macrophages (ADGRE1 + ), neutrophiles (LY6G + ), T (CD3 + ), and B lymphocytes (CD19 + ).N = 3 biological replicates.Two-tailed T-tests.e. Strategy for fascia fibroblast purification and culture (left).Expression changes (right) of proinflammatory (Ccl2 and Cxcl1) and myofibroblast markers (Acta2) in IL1β-(inflammation-inducing) or TGFβ1containing media (myofibroblast-inducing).N = 3 technical replicates.Expression changes normalized to control medium.p values on bars from two-tailed T-tests (top) and 1-way ANOVA comparisons between treatments (bottom).f-g.Expression changes of indicated markers in inflammation-(f) and myofibroblast-inducing media (g) exposed to exogenous RA at indicated concentrations.N = 3 technical replicates.Expression changes normalized to control IL1β-(f) or TGFβ1-containing medium (g).p values on bars from twotailed T-tests (top) and 1-way ANOVA comparisons between treatments (bottom).Extended Data Fig. 8 | HIF1α instructs the transition into (proto) myofibroblasts.a.Control (top) or HIF1α inhibitor-treated (bottom) explants after 6 days of culture (left), area over time (middle), and contraction (right) measurements.N = 6 (control) and 3 (inhibitor) replicates.b.Immunolabeling (left) and marker-positive cell ratios (right) of proinflammatory fibroblast (PDPN + ) or myofibroblasts (αSMA + ) in control (top) or inhibitor-treated (bottom) explants.N = 3 biological replicates.c.Immunolabeling (top left) and YAP1 + cell ratios (bottom left) upon chemical TAZ agonism.N = 3 (controls) and 4 (agonist) biological replicates.TAZ agonist-treated (centre top) and combined treatments with the inhibitor (centre middle) or RA (centre bottom) explants after 6 days of culture.Area over time (top right), and contraction (bottom right) measurements.N = 3 biological replicates.d.TGFβ1-treated (top left) and combined treatments with the inhibitor (centre left) or RA (bottom left) explants after 6 days of culture.Area over time (top right), and contraction (bottom right) measurements.N = 3 biological replicates.e. Immunolabeling (top) and pSMAD2 + cell ratios (bottom).N = 3 biological replicates.f.Correlation of HIF1α regulon vs hypoxia activity, positively regulated genes involved in TGFβ pathway (bottom).g-h.UMAPs of gene expression and schematic representation of HIF1α-regulated TGFβ modulators.i. Representative photographs (left) of control (top) or inhibitor-treated wounds (bottom) at 7 dpi.Wound area measurements over time (top right) and overall contraction at 7 dpi (bottom right).N = 3 biological replicates.j.Trichromestained micrographs of control (top) and treated wounds (bottom) at 7 dpi.k.Micrographs of control (top left) and treated 3 dpi wounds (bottom left) immunolabeled for CD45/PTPRC and positive cell ratios (right).N = 6 (control) and 5 (treated) wounds from 3 biological replicates.All p-values from twotailed T-tests.Scale bars: 2 mm in explants photographs (a, c-d), 500 microns (g), and 50 microns (b-c, e, and h).Extended Data Fig. 9 | Fate mapping of fascia extracellular matrix.a. Schematics of fascia-matrix labelling using NHS-Pacific Blue.b.Representative low magnification (top) and high magnification (bottom) images of control wounds immunolabeled for PDPN or pSTAT3, showing association with traced extracellular matrix.Scale bars: 500 microns in low magnification and 50 microns in high magnification micrographs.Extended Data Fig. 10 | Comparative analysis of human and mouse fibroblasts subtypes.a-b.UMAP representation of scRNA-seq datasets of denoted human skin pathologies color-coded for sample of origin (a) and the pan-fibroblast marker PDGFRA expression levels (b).c.Composition charts depicting the proportion of the different fibroblast subtypes for each dataset.

d.
Expression comparison between merged human (red) and mouse datasets (blue) of top marker genes for the human fibroblast clusters.Expression of analogous fibroblast clusters is denoted in red rectangles.e. UMAP representations the expression levels of highly conserved markers in the merged human (top) and mouse (bottom) datasets.