Almost all tumours carry one or more cancer driver mutations, which are essential for cell transformation. However, recent advances in cancer genomics have demonstrated that normal human tissues contain millions of cells carrying known driver mutations, while preserving homeostasis. Most of these mutated cells will never transform into tumours. Moreover, studies of known or suspected human carcinogens have shown that the majority are not mutagens. These observations suggest that exogenous carcinogenic exposures might increase cancer risk by modifying selective constraints, promoting the expansion of pre-existing clones carrying specific oncogenic mutations. In this Review, we propose a synthesis between ideas put forward almost a century ago based on seminal experiments on carcinogen-induced tumours in mice, observations made by cancer epidemiologists over several decades, and the recent revelation that normal human tissues are a patchwork of mutant clones. The repeated interplay between variation and selection-the first principles of Darwinian evolution-underlies the clonal selection leading to tumorigenesis. A deeper understanding of these processes can enhance prospects for cancer prevention by eliminating or mitigating the effects of environmental or endogenous tumour promoters.

The identification of the cancer cell of origin is a fundamental issue in cancer biology. We used fluorescent lineage tracing of independent mouse skin stem cell populations, single-cell transcriptomics, and duplex sequencing to identify the origin of chemically induced skin tumors. Tumors arose predominantly from Lgr6+ and/or Lrig1+ stem cells of the upper hair follicle but only very rarely from the Lgr5+ and Krt19+ hair follicle bulge. Lgr6+ stem cells initiated by dimethylbenzanthracene responded to tumor promoter treatment, resulting in clonal expansion of initiated cells carrying the canonical Hras Q61L mutation. Spontaneous mutations in Kras also clonally expanded but did not generate tumors unless the Hras gene was deleted, thus revealing a competitive interaction between the Hras and Kras pathways that influences clonal selection.

Historical studies performed nearly a century ago using mouse skin models identified two key steps in cancer evolution: initiation, a likely mutational event, and promotion, driven by inflammation and cell proliferation. Initiation was proposed to be permanent, with promotion as the critical rate-limiting step for cancer development. In this study, we carried out whole-genome sequencing to demonstrate that initiated cells with thousands of mutagen-induced mutations can persist for long periods and are not removed by cell competition or by immune intervention, thus mimicking the persistence of cells with cancer driver mutations in normal human tissues. In the mouse, these cells do not give rise to tumors unless exposed to the tumor promoter 12-O-tetradecanoylphorbol-13-acetate (TPA). Tissue damage and regenerative proliferation, but not normal cell turnover, consistently trigger tumor formation. Wounding, promoter treatment, and obesity enhance promotion without increasing mutational burden, supporting the possibility of future cancer prevention efforts directed at promotional risk factors. SIGNIFICANCE:Using historical skin cancer models, we reveal that mutated cells can persist without tumor formation and give rise to cancer upon exposure to tumor promoters, underscoring the importance of tumor promotion over initiation as the rate-limiting step in carcinogenesis and the need for cancer prevention strategies targeting promotional factors. See related commentary by Lim and DeGregori, p. 1090.

While the causal role of mutagenic carcinogens in tumor development is well established, the relative contribution of environmental tumor promoting factors, wounding, and chronic inflammation is still unclear. Recent sequencing studies have suggested that most environmental carcinogens act as promoters rather than through mechanisms that involve direct induction of point mutations, but whether cancer risk factors such as obesity, chronic inflammation, wounding, or tumor promoters contribute directly or indirectly to mutation burden, or induce novel signatures, has not been investigated. Here, we present WGS analysis of over 100 mouse skin tumors to compare the effects of exposure to mutagens, the tumor promoter TPA, chronic wounding, obesity, or chemotherapy, on mutational burden and cancer risk. All tumors initiated by the carcinogen Dimethylbenzanthracene (DMBA) show a very strong A>T mutational signature (SBS.DMBA) attributable to a single exposure to this carcinogen. The number of SBS.DMBA mutations also showed a strong correlation with the “clock” signature SBS5, suggesting that one treatment with this mutagen can induce mutational signatures attributed to endogenous processes. No specific signatures could be attributed to obesity, high fat diet, wounding, or TPA. Cells carrying thousands of mutations persist over very long periods without inducing tumors or causing pathological changes but can give rise to tumors after short term exposure to TPA. Furthermore, normal cell turnover and proliferation during fetal and adult growth, is not sufficient for promotion, but tissue damage followed by regenerative proliferation seems to be required for tumor development. We conclude that tumor promoters, chronic inflammation, wounding, and obesity do not contribute significantly to tumor mutational burden, and that the rate-limiting determinant of tumor growth is exposure to a tumor promoter rather than the nature or number of genomic point mutations. These data are highly relevant to the recent demonstration of persistent oncogenic mutations in histologically normal human tissues during ageing.### Competing Interest StatementThe authors have declared no competing interest.

Stem cells play a critical role in cancer development by contributing to cell heterogeneity, lineage plasticity, and drug resistance. We created gene expression networks from hundreds of mouse tissue samples (both normal and tumor) and integrated these with lineage tracing and single-cell RNA-seq, to identify convergence of cell states in premalignant tumor cells expressing markers of lineage plasticity and drug resistance. Two of these cell states representing multilineage plasticity or proliferation were inversely correlated, suggesting a mutually exclusive relationship. Treatment of carcinomas in vivo with chemotherapy repressed the proliferative state and activated multilineage plasticity whereas inhibition of differentiation repressed plasticity and potentiated responses to cell cycle inhibitors. Manipulation of this cell state transition point may provide a source of potential combinatorial targets for cancer therapy.

Adult mammalian stem cells play critical roles in normal tissue homeostasis, as well as in tumor development, by contributing to cell heterogeneity, plasticity, and development of drug resistance. The relationship between different types of normal and cancer stem cells is highly controversial and poorly understood. Here, we carried out gene expression network analysis of normal and tumor samples from genetically heterogeneous mice to create network metagenes for visualization of stem-cell networks, rather than individual stem-cell markers, at the single-cell level during multistage carcinogenesis. We combined this approach with lineage tracing and single-cell RNASeq of stem cells and their progeny, identifying a previously unrecognized hierarchy in which Lgr6+ stem cells from tumors generate progeny that express a range of other stem-cell markers including Sox2, Pitx1, Foxa1, Klf5, and Cd44. Our data identify a convergence of multiple stem-cell and tumor-suppressor pathways in benign tumor cells expressing markers of lineage plasticity and oxidative stress. This same single-cell population expresses network metagenes corresponding to markers of cancer drug resistance in human tumors of the skin, lung and prostate. Treatment of mouse squamous carcinomas in vivo with the chemotherapeutic cis-platin resulted in elevated expression of the genes that mark this cell population. Our data have allowed us to create a simplified model of multistage carcinogenesis that identifies distinct stem-cell states at different stages of tumor progression, thereby identifying networks involved in lineage plasticity, drug resistance, and immune surveillance, providing a rich source of potential targets for cancer therapy.

An estimated 40% of all human cancers are suspected to be a result of modifiable risk factors such as obesity, high fat diet and chronic inflammation. Whole genome sequencing (WGS) of thousands of human tumors have revealed “mutational signatures” that provide a molecular footprint of cancer origins. Whether such signatures exist for modifiable cancer risk factors remains unclear. We studied the impact of lifestyle risk factors using a compendium of 107 mouse tumors that model obesity, high fat diet, wounding, chronic inflammation, or chemotherapy. We used a well-established 2-step skin carcinogenesis model composed of exposure to mutagen DMBA followed by tumor promoter TPA, generating squamous carcinomas that were analysed by WGS for identification of mutational signatures. In addition to recapitulating many COSMIC human signatures, we identified a novel SNV signature induced by a single treatment with DMBA (SBS.DMBA) which explains the majority of all detected mutations. While a single exposure of normal skin to DMBA induces a highly variable number of carcinogen-specific mutations, a very high mutation burden is insufficient for tumorigenesis. SNVs attributable to reactive oxygen species (ROS) are broadly found in about 25% of all mouse tumors, but are most prominent in tumors from mice that are exposed to DMBA in utero, suggesting that the developmental age of mutagen exposure may impact the repair of ROS generated mutations and that the timing of exposure is a poorly understudied component of carcinogenesis. In mouse tumor models of genetic and dietary obesity, the total mutational load and mutational signatures in tumors from obese mice were indistinguishable from those of lean mice despite dramatic differences in tumor latency and progression as well transcriptomic differences in immune activation. We found an enrichment in deleterious Tgfrb2 mutations in tumors from mice with low compared to high body mass index (BMI) (p<0.016). Using a conditionally activated RAS mouse model, we also show that a non-mutagenic inflammatory signal such as a chronic wound can act as the rate-limiting step for full tumor development. Surprisingly these tumors can be evoked even in somatic genomes with very few mutations apart from the initiating Ras driver. Finally, together with the Riva et al study, these chemically induced mouse tumors recapitulate >50% of established human cancer driver genes. DMBA caused 91% of all Hras/Kras mutations, but only a minority of other recurrent driver mutations in genes such as Trp53 and Tert, suggesting that these occur later during the process of carcinogenesis. Taken together, we have analyzed the largest compendium of WGS data from nearly 300 mouse tumors, showing that while exogenous promotional factors do not increase mutation burden or induce novel mutational patterns, they have a major rate-limiting role in determining cancer risk. Citation Format: Yun Rose Li, Kyle Halliwill, Eve Kandyba, Reyno Delrosario, Quan Tran, Nora Bayani, Di Wu, Olga Mirzoeva, Melissa McCreery Reeves, Mishu Islam, Laura Riva, Eric Bergstrom, John Digiovanni, Ludmil Alexandrov, Allan Balmain. The impact of carcinogens, obesity, and chronic inflammatory processes on mutational signatures and cancer risk in mouse tumor models [abstract]. In: Proceedings of the American Association for Cancer Research Annual Meeting 2022; 2022 Apr 8-13. Philadelphia (PA): AACR; Cancer Res 2022;82(12_Suppl):Abstract nr 2198.

Tumour cells are subjected to evolutionary selection pressures during progression from initiation to metastasis. We analysed the clonal evolution of squamous skin carcinomas induced by DMBA/TPA treatment using the K5CreER-Confetti mouse and stage-specific lineage tracing. We show that benign tumours are polyclonal, but only one population contains the Hras driver mutation. Thus, benign papillomas are monoclonal in origin but recruit neighbouring epithelial cells during growth. Papillomas that never progress to malignancy retain several distinct clones, whereas progression to carcinoma is associated with a clonal sweep. Newly generated clones within carcinomas demonstrate intratumoural invasion and clonal intermixing, often giving rise to metastases containing two or more distinct clones derived from the matched primary tumour. These data demonstrate that late-stage tumour progression and dissemination are governed by evolutionary selection pressures that operate at a multicellular level and, therefore, differ from the clonal events that drive initiation and the benign-malignant transition.

Muscle damage in Duchenne muscular dystrophy (DMD), the incurable diseases caused by the lack of dystrophin, is strongly linked to inflammatory reactions. Accordingly, modulation of inflammation can ameliorate the disease progression. In our recent studies we have generated the double knockout (dKO) mdx mice lacking additionally either the active heme oxygenase-1 gene (HO-1, Hmox1) or miR-378a coding sequence. The HO-1 is the crucial anti-inflammatory protein affecting also muscle stem (satellite) cells differentiation, while miR-378a is involved in regulation of muscle differentiation and cellular metabolism. Lack of HO-1 aggravates the severity of DMD in a mdx mice and Hmox1-/-mdx dKO animals demonstrate strongly impaired exercise capacity. Interestingly, mdx satellite cells (SCs) show disturbed and enhanced differentiation, which is further aggravated by Hmox1 deficiency. Interestingly, mdx SCs expressed lower level of Hmox1, and RNA sequencing demonstrated the significant changes in mdx SCs transcriptome. Importantly, restoration of Hmox1 expression normalized the SCs differentiation. Surprisingly, mice lacking miR-378a were able to run longer distance than their wild-type (WT) in a treadmill test. Lack of miR-378a ameliorated the impairment characteristic for mdx mice, as 3-months old dKO animals were stronger than mdx mice and ran similar distance to WT. The muscles of dKO animals showed less infiltration of inflammatory cells. RNA-sequencing revealed the fibroblast growth factor 1 (Fgf1) among the most downregulated genes in animals lacking miR-378-/and in dKO mice, what might be of relevance to lower fibrosis observed in dKO than mdx mice In summary, HO-1 and miR-378a can be considered as targets to alleviate muscle injury and DMD severity.

The cutaneous wound healing program is a product of a complex interplay among diverse cell types within the skin. One fundamental process mediated by these reciprocal interactions is the mobilization of local stem cells pools to promote tissue regeneration and repair. Using the epidermis specific ablation of caspase-8 as a model of mouse wound healing, we analyzed the signaling components responsible for epithelial stem cell proliferation. We found that IL-1α and IL-7 secreted from keratinocytes function synergistically to expand the activated population of resident epidermal γδ T-cells. A downstream effect of activated γδ T-cells is the preferential proliferation of hair follicle stem cells. On the other hand, IL-1α dependent stimulation of dermal fibroblasts optimally stimulates epidermal stem cell proliferation. These findings provide new mechanistic insights into the regulation and function of epidermal-immune cell interactions and how components classically associated with inflammation can differentially influence distinct stem cell niches within a tissue.

The cutaneous wound-healing program is a product of a complex interplay among diverse cell types within the skin. One fundamental process that is mediated by these reciprocal interactions is the mobilization of local stem cell pools to promote tissue regeneration and repair. Using the ablation of epidermal caspase-8 as a model of wound healing in Mus musculus, we analyzed the signaling components responsible for epithelial stem cell proliferation. We found that IL-1α and IL-7 secreted from keratinocytes work in tandem to expand the activated population of resident epidermal γδT-cells. A downstream effect of activated γδT-cells is the preferential proliferation of hair follicle stem cells. By contrast, IL-1α-dependent stimulation of dermal fibroblasts optimally stimulates epidermal stem cell proliferation. These findings provide new mechanistic insights into the regulation and function of epidermal cell–immune cell interactions and into how components that are classically associated with inflammation can differentially influence distinct stem cell niches within a tissue.

Article Figures and data Abstract eLife digest Introduction Results Discussion Materials and methods Appendix 1 Appendix 2 Data availability References Decision letter Author response Article and author information Metrics Abstract The control principles behind robust cyclic regeneration of hair follicles (HFs) remain unclear. Using multi-scale modeling, we show that coupling inhibitors and activators with physical growth of HFs is sufficient to drive periodicity and excitability of hair regeneration. Model simulations and experimental data reveal that mouse skin behaves as a heterogeneous regenerative field, composed of anatomical domains where HFs have distinct cycling dynamics. Interactions between fast-cycling chin and ventral HFs and slow-cycling dorsal HFs produce bilaterally symmetric patterns. Ear skin behaves as a hyper-refractory domain with HFs in extended rest phase. Such hyper-refractivity relates to high levels of BMP ligands and WNT antagonists, in part expressed by ear-specific cartilage and muscle. Hair growth stops at the boundaries with hyper-refractory ears and anatomically discontinuous eyelids, generating wave-breaking effects. We posit that similar mechanisms for coupled regeneration with dominant activator, hyper-refractory, and wave-breaker regions can operate in other actively renewing organs. https://doi.org/10.7554/eLife.22772.001 eLife digest Skin includes hundreds of thousands of hair follicles that cycle through different stages of activity. Each follicle grows hair, sometimes (in the case of long hairs like human head hair and horse tail hairs) for several years, before losing it. The follicle then goes through a resting stage before starting to grow another hair. To achieve high hair density, the follicles need to coordinate their hair-making activities. If they all worked independently from one another, bald patches would inevitably form that would compromise how effectively the skin works. Groups of cells can communicate using a variety of chemical signals. It was not known whether cells in hair follicles from different regions of the skin rely on the same signals to communicate, and whether follicles in neighboring regions are able to ‘understand’ one another. Through a combination of mathematical modeling and experimental results from mice, Wang, Oh et al. now show that hair follicles across the body use a common signaling system. This system consists of a pair of signals: ‘activators’ that stimulate hair growth, and ‘inhibitors’ that prevent it. The balance between these two signals affects the pattern of hair growth. For example, higher levels of activators allow fur to grow thickly on the belly of the mouse, likely to protect against heat loss and injuries from the ground. By contract, higher levels of inhibitors make the hairs on the ear sparse, which may prevent them from interfering with hearing. There is little evidence that hair follicles on the scalp communicate in adult humans. Learning to activate and control communication between these follicles could provide a way to treat male pattern baldness and similar conditions. Understanding how hair follicles communicate may also help researchers to develop ways of regenerating other fast-renewing organs, such as the gut and bone marrow. https://doi.org/10.7554/eLife.22772.002 Introduction Featuring prominent growth cycles, the hair follicle (HF) is a model system of choice for studying tissue regeneration. At the level of cellular activities, the hair growth cycle consists of three consecutive phases: anagen, phase of active proliferation; catagen, apoptotic involution phase; and telogen, relative quiescence phase (Al-Nuaimi et al., 2010; Paus and Foitzik, 2004; Schneider et al., 2009; Stenn and Paus, 2001). Cyclic regeneration is sustained by the bulge stem cells, located at the base of the permanent HF portion (Cotsarelis et al., 1990). During anagen initiation, signals from the niche, including the dermal papilla (DP), stimulate bulge stem cells and adjacent hair germ (HG) progenitors to proliferate (Enshell-Seijffers et al., 2010; Greco et al., 2009; Legrand et al., 2016). Activated progenitors generate all lower HF structures, including the outer root sheath (ORS) and hair matrix. During catagen, a widespread apoptotic program remodels the HF back toward a telogen state (Botchkarev et al., 2001b; Fessing et al., 2006; Foitzik et al., 2000; Lindner et al., 1997; Mesa et al., 2015). Conceptually, since the bulge produces downward migrating progeny (Hsu et al., 2011), it effectively serves as a progenitor source, while the matrix functions as a sink, and the ORS as a channel for progenitors transiting between them. The signaling mechanisms that time these coordinated cellular activities during hair regeneration remain incompletely understood (Al-Nuaimi et al., 2014; Bernard, 2012; Lin et al., 2009; Paus et al., 1999). The putative ‘hair cycle clock’ is thought to be composed of one or several activator/inhibitor pairs acting to time key cycle phase transitions at set thresholds of their activities. Accordingly, cycle pace will depend on the speed at which activators and inhibitors reach their respective thresholds (Chen et al., 2015). Importantly, HFs exist as large populations and at least in the dorsal skin they interact to coordinate growth cycles (Hodgson et al., 2014; Plikus et al., 2011, Plikus and Chuong, 2008a, Plikus et al., 2008b). Such coordination implies that at least some of the activators and inhibitors should be present between HFs, in the so-called skin macro-environment. Previous work on dorsal skin indicates that BMP and WNT pathways constitute important components of the hair cycle clock. Indeed, defects in either of these pathways can dramatically change hair cycle progression (Botchkarev et al., 2001a; Botchkarev and Sharov, 2004; Choi et al., 2013; Enshell-Seijffers et al., 2010; Kandyba and Kobielak, 2014; Kandyba et al., 2013; Kobielak et al., 2003, 2007; Sharov et al., 2005, 2006), and ligands and antagonists for both pathways mediate macro-environmental coordination between HFs (Chen et al., 2014; Plikus et al., 2011; Plikus and Chuong, 2014). Additionally, FGF, PDGF, TGFβ, TNFα and other pathways can modulate hair cycle (Chen et al., 2015; Festa et al., 2011; Higgins et al., 2014; Ito et al., 2003; Kimura-Ueki et al., 2012; Leishman et al., 2013; Oshimori and Fuchs, 2012; Plikus, 2012; Rivera-Gonzalez et al., 2016). Importantly, the combined signaling activities for the above pathways partition the hair cycle in the dorsal skin into four functional phases, each with its distinct activator/inhibitor profile: propagating anagen (P), autonomous anagen (A), refractory telogen (R) and competent telogen (C) (Hodgson et al., 2014; Plikus and Chuong, 2014; Plikus et al., 2008b). Interactions between HFs enable hair regeneration across dorsal skin to self-organize into dynamic patterns. Critical for this self-organization are the following HF-to-HF interactions: P-phase HFs can induce neighboring C-phase HFs to enter anagen via diffusible activators, leading to hair growth coupling, while A-phase anagen or R-phase telogen HFs cannot couple due to high levels of inhibitors (Murray et al., 2012; Plikus et al., 2011; Plikus and Chuong, 2014; Plikus et al., 2008b). It remains unknown, however, whether this self-organization mechanism and its underlying WNT/BMP signaling activities is a general feature of all body skin or a special case for dorsal skin only. Integrative understanding of large-scale hair regeneration requires a systems biology approach. Previous modeling on HFs include cellular automaton models (Halloy et al., 2000; Plikus et al., 2011), feedback-control model (Al-Nuaimi et al., 2012) and the FitzHugh-Nagumo (FHN) excitable medium model (Murray et al., 2012). Here, we present a unified three-dimensional and stochastic modeling framework for the HF that captures: (i) activator/inhibitor signaling dynamics in a single HF, (ii) cyclic growth of a single HF, and (iii) coupling between multiple HFs through diffusive signals. Using this model, we reveal that skin as a whole behaves as a heterogeneous regenerative field, where: (a) dominant hair cycle waves start in the ventrum, (b) propagate dorsally in a bilateral pattern, (c) stop at the boundary with hyper-refractory ear skin, and (d) break at non-propagating anatomical landmarks, such as eyelids and ears. We also show that WNT and BMP serve as a universal activator/inhibitor signaling pair, whose varying activities underlie distinct hair regeneration dynamics in all anatomical locations studied. These results provide new understanding of how the entire skin of the animal manages all of its hair regeneration. Results A multi-scale model recapitulates a single growing HF, as well as HF-to-HF communication First, in modeling the geometry of a single HF, we considered four key expression sites for activator/inhibitor ligands, antagonists, and receptors along the HF axis: bulge, HG, matrix, and DP (Figure 1A). During the cycle, the bulge (assigned as Region I) remains relatively static, whereas the DP moves up and down along the HF axis. Also dynamic are the HG and matrix. The former only exists during telogen, while the latter only exists during anagen. The HG grows down to make matrix during anagen onset, whereas during catagen, the matrix collapses, and a new HG reforms. Simplistically, cyclic HG→matrix→HG dynamics are coordinated with the DP; thus, in the model we identify them jointly as Region II. Next, we considered that both regions produce signaling factors. Although a biological simplification, we assumed that Region I does so at a rather constant rate, while Region II shows distinct temporal dynamics (Appendix 2—table 4). We also assumed that Region II is essential for sending hair cycle-promoting signal(s), while Region I is the primary signal target. In short, we hypothesized that the essential temporal molecular dynamics in the HF operate as follows: Region II generates a signaling ligand (L) gradient; Region I detects it and transmits it into ligand-bound receptors (LR) that then, through a series of intermediate signaling steps not captured in the model directly (such as activities of the downstream signaling pathways and involvement of additional cell populations), regulates cyclic HF growth (Figure 1A). The molecular signaling events, either activating or inhibitory, can be summarized as: (1)∂∂t[L] = Diffusion+Production+Reaction of L and R(2)∂∂t[LR] = Reaction of L and R +Degradation +Extra Source Figure 1 Download asset Open asset Model recapitulates hair cycling and its associated activator and inhibitor signaling dynamics. (A) Schematic depiction of HF growth dynamics during telogen and anagen. Telogen and anagen HFs are shown on the left and in the center, respectively. In both hair cycle phases, Region I (purple) represents bulge and Region II (orange) represents DP with HG during telogen phase, and DP with matrix during anagen phase. On the right, schematic drawing of diffusive activator (Act. L. in green) and inhibitor (Inh. L. in red) interactions with their corresponding receptors (Act. R and Inh. R, not depicted) that form ligand-bond-receptors (Act. LR and Inh. LR) and their coupling with physical growth of the HF (blue) is shown. (B) Typical noise-free dynamics of the activator (green) and inhibitor (red) and cyclic HF growth (blue) are shown. X-axis is time in simulated days. Y-axis for activator and inhibitor shows simulated signaling levels, and for HF growth – simulated length of the HF. Grey area demarcates one modeled hair growth cycle. (C) The duration of ~anagen and ~telogen phases as the function of inhibitor signaling strengths. X-axis shows modeled inhibitor levels with ‘0’ being an arbitrary baseline levels. Y-axis shows time in simulated days. Upon stronger inhibitory signaling (high Inh. L level) ~anagen shortens (yellow) and ~telogen lengthens (purple). The entire cycle (blue) becomes longer either with stronger or weaker inhibitory signaling. When inhibitory signaling becomes either very strong or very weak, the excitability of the system breaks down and HFs equilibrate in one state (grey regions). Also see Appendix 2—tables 1, 2 and 4. (D–E’’) A total of 236 putative activator genes (green) and 122 putative inhibitor genes (red) available from a whole skin microarray dataset were identified to recapitulate temporal dynamics of the simulated activator (D) and inhibitor (E), respectively. Multiple WNT pathway members are in the putative activator gene set (D’, D’’), while BMP pathway members are among the putative inhibitor genes (E’, E’’). See gene list in Dataset 1. For all genes log-transformed, zero-mean expression profile values were calculated using colorimetric ratio-scale algorithm as reported in (Lin et al., 2009). https://doi.org/10.7554/eLife.22772.003 where L, R and LR stand for ligands, receptors, and ligand-bound receptors, respectively. In the dynamics of LR (Equation 2), the ‘Extra Source’ describes stochastic signaling effects due to noise, and potential signaling contributions from Region I (Appendix 2-Governing equations for activators and inhibitors). As Equations 1 and 2 show, ligand-receptor interactions in the model take place only for the same signaling pathway, and no direct pathway cross-talk is set to occur. This, again, is a biological simplification. Recently, evidence for pathway interactions have emerged (Kandyba et al., 2013), and its effect is explored in Appendix 2-Possible interactions between the activator and inhibitor pathways do not qualitatively alter the HF dynamics. Our model integrates key signaling features of the hair growth cycle: strong activator signals enhance HF growth, while strong inhibitor signals prevent it. We modeled HF growth through the spatial average of LR concentration differences between the levels of activator and inhibitor in Region I (Equation 7 in Appendix 2-Modeling HF phases by concentration difference). We assumed the hair cycle has two critical ‘checkpoints’: (i) the event in late competent telogen, when production of activator starts to increase (Chen et al., 2014; Greco et al., 2009; Oshimori and Fuchs, 2012; Plikus et al., 2008b), and (ii) the event of anagen termination, when the HF starts to involute. Thus, our model recognizes two phases determined by these checkpoints: ~anagen, starting from the moment of activator amplification until anagen termination, and ~telogen, lasting until the next activator amplification event. In the context of the conventional hair growth cycle, ~anagen incorporates the late portion of competent telogen and the entire anagen, while ~telogen includes catagen, refractory telogen and the remainder of competent telogen (Plikus et al., 2011; Plikus and Chuong, 2014; Plikus et al., 2008b) (Appendix 2-Modeling HF phases by concentration difference; Appendix 2—figure 2). Model simulations produce several emergent behaviors. The cycle becomes autonomous – that is, it displays stable periodicity and excitability emerges naturally without a built-in ‘clock’ (Figure 1B). Cycling is maintained within a range of parameter values, allowing testing for various intrinsic and extrinsic signaling scenarios (Figure 1C). Associated with these dynamics are periodic changes in the system’s geometry – the signaling source in Region II moves cyclically. Simulations indicate that the moving HF geometry in the model is critical, greatly contributing to the regulation of the cycle. In a single HF model, activator/inhibitor diffusion occurs only along the HF axis. When a HF population is modeled, hair-to-hair communication emerges naturally as ligand diffusion from neighbors supplements intrinsic HF ligand levels. As such, hair cycle pace depends on interactive signaling between neighboring HFs – a feature that we explore below. HF cycling emerges from the growth-mediated coupling of activator and inhibitor Our model predicts that HF cycling occurs only within a certain range of signal strengths, that is the excitable regime (Figure 1C, white region). Within this regime, activator and inhibitor are predicted to inversely modulate duration of both ~telogen and ~anagen phases. At certain, either too high or too low signal strengths, the excitability is predicted to break down and the HF is expected to enter a non-cycling state of equilibrium (Figure 1C, grey regions). For example, when inhibitor levels are very high, the HF is predicted to equilibrate in an extended telogen (Appendix 2—figure 5A), while extended anagen is predicted for the opposite signaling condition (Appendix 2—figure 5B). Next, we used bioinformatic and experimental approaches to validate the model’s key prediction that the same activator or inhibitor pathway can inversely modulate telogen and anagen phase duration. Considering the established roles for BMP and WNT as respective inhibitor and activator pathways regulating telogen duration in the dorsal skin, we explored if they can also regulate anagen duration in the same skin region in a model-predicted fashion. First, we found that model-predicted temporal dynamics for inhibitor and activator during ~anagen (Figure 1D and 1E) match the actual anagen expression dynamics for multiple BMP and WNT pathway members established on a highly temporally resolved whole-tissue dorsal skin microarray dataset (Lin et al., 2009) (Figure 1D–E’’; Appendix 1-Identifying model predicted hair cycle activators and inhibitors). We also show that perturbing BMP (for details see Appendix 1-Validating model-predicted roles for BMP signaling in hair cycle control) and WNT in transgenic mice (for details see Appendix 1-Validating model-predicted roles for WNT signaling in hair cycle control) alters dorsal anagen phase duration and leads to hair length defects in a way that is consistent with the model’s predictions. Overall, this data shows that our model generates biologically meaningful outcomes and that its predictive power is robust. Model reveals skin is a heterogeneous regenerative field Next, we set out to explore novel aspects of hair regeneration at the population level. For this purpose, we modeled a linear array of HFs (i.e. two-dimensional organization; Appendix 2—figure 4A) and a grid of HFs (i.e. three-dimensional organization; Appendix 2—figure 4B). In both cases, the diffusion of activators and inhibitors accompanying each HF during growth naturally led to HF coupling (Appendix 1-Validating model-predicted roles for BMP signaling in hair cycle control) and emergence of several known features of collective hair growth behavior, including spontaneous anagen initiation and anagen wave spreading (Appendix 2—figures 11, 12). We then focused on the phenomenon of bilaterally symmetric hair growth that is prominent in young mice (Plikus et al., 2008b) yet remains unexplained. Conventionally, first anagen in the dorsal skin of newborn mice is considered synchronous. On the other hand, adult mice display fully asynchronous and asymmetric dorsal hair growth patterns (Chen et al., 2014; Plikus and Chuong, 2008a; Plikus et al., 2009). This, however, is preceded by prominent bilateral symmetry, which often persists into the fourth hair cycle (Plikus and Chuong, 2008a). We now show that in the three-dimensional model where all HFs are assumed to be identical, full asynchrony evolves within just one cycle, and bilateral symmetry cannot be achieved (Appendix 2-Dorsal and ventral HF patterns; Appendix 2—figures 11, 12; Appendix 2—video 1). Therefore, we hypothesized that first anagen is inherently asynchronous as a result of spatially patterned HF development. Indeed, spatial distribution of early anagen HFs in the dorsal skin of newborn mice (Figure 2A–D) reveals head-to-tail and subtle lateral-to-medial asynchronies. We modeled the impact of these asynchronies on hair growth pattern evolution. Simulations reproduced head-to-tail asynchrony (Appendix 2-Dorsal and ventral HF patterns; Appendix 2—figures 14, 15; Appendix 2—video 3); however, it persisted for at least 10 cycles, which is far more than the 3–4 cycles observable in mice. Moreover, prominent bilateral symmetry failed to form. Figure 2 Download asset Open asset Spatiotemporal patterning of early hair cycles. (A–D) Analysis of the whole mount dorsal skin samples from P0 (n = 3) (A) and P1 WT mice (n = 3) (A’) reveals subtle head-to-tail and lateral-to-medial hair cycle asynchronies. Asynchronies were inferred from examining the size of pigmented HFs. Larger HFs result from earlier anagen onset. (B) Heatmaps of skin samples from A and A’ built based on black pixel density (reflecting pigmented anagen HFs). (C) Quantification of anagen HFs at different phases confirms head-to-tail pattern asynchrony. Morphological definition of anagen phases used for this analysis is provided at the bottom on the panel. (D) Analysis of the whole mount dorsal skin samples from P0 220bpMsx2-hsplacZ mice, where lacZ reporter activates in anagen HFs starting from phase IIIb, confirms head-to-tail asynchrony. (E, F) Modeling rapid hair growth pattern evolution in the context of two heterogeneous domains. (E) Schematic depiction of the modeling conditions with central Dorsal domain flanked by two lateral Ventral sub-domains with coupling between Dorsal and Ventral HFs. (F) Compared to Dorsal domain HFs, Ventral domain HFs were assigned with higher levels of total available activator and inhibitor receptors, allowing shorter ~anagen and ~telogen duration. Furthermore, hair cycle asynchrony was introduced into Dorsal domain to model the initial head-to-tail asynchrony. In simulations, interactions between HFs across domain boundaries result in bilateral symmetry during the second cycle (simulated time t68-78.5). Also, initial asynchrony breaks down in the cycle 3 (t130), and partial bilateral symmetry maintains into the late cycles (see Appendix 2—video 4). Scale bars: A, A’, D – 5 mm. Images on A, A’ and D are composites. https://doi.org/10.7554/eLife.22772.004 We note that the above and previous simulations (Murray et al., 2012; Plikus et al., 2011) were performed on homogenous HF populations, where all HFs are assumed to be identical. We then considered that novel patterns might develop upon interaction of two or more HF populations, whose activator/inhibitor signaling levels are inherently different. In principle, dorsal skin HFs can interact with HFs from other body regions, such as ventral skin, where hair cycle dynamics are potentially distinct. Because all skin is continuous and forms an approximation of a cylinder, we modeled it as an unrolled sheet, where two Ventral sub-domains flank a rectangular Dorsal domain (Figure 2E). For initial modeling conditions (Appendix 2-Dorsal and ventral HF patterns), we considered that: (i) the first cycle on the dorsal skin has built-in head-to-tail asynchrony, and that (ii) ventral HFs develop with a 3- to 4-day delay relative to dorsal HFs (Appendix 1—figures 7–9). Because ventral HFs are known to produce distinctly shorter hairs (Candille et al., 2004), in the model we assumed that they have faster cycle dynamics compared to dorsal HFs (Appendix 2-Changes in the total amount of activator and inhibitor receptors results in different sensitivity of ~anagen and ~telogen lengths to signaling changes, Appendix 2-Dorsal and ventral HF patterns; Appendix 2—figures 8, 9). Indeed, in this configuration, our model readily reproduced patterns with aspects of bilateral symmetry already in the second cycle as the result of dominant waves spreading from the Ventral to the Dorsal domain (Figure 2F, t68-78.5). Importantly, after the second cycle, the effect of the initial built-in head-to-tail asynchrony started to disappear. Instead, the interaction between Ventral and Dorsal HFs continued to produce prominent bilateral symmetry in the third (Figures 2F, t130-145) and later cycles (Appendix 2—video 4). Taken together, the model predicts that rapid hair growth pattern evolution requires interaction of two or more skin domains with distinct hair cycle parameters. Ventral-dorsal interactions produce bilaterally symmetric hair growth patterns Next, we imaged Flash mice, whose luciferase reporter produces skin-specific WNT activity signal and allows to sensitively and non-invasively determine the location and percentage of anagen HFs across the entire body (Hodgson et al., 2014) (Figure 3A–C). Luminescence levels were measured both dorsally and ventrally and mice were followed up until day P119, encompassing up to five hair growth cycles. Combined analysis from multiple mice reveals prominent phase advancement in ventral over dorsal anagen, specifically during the second, third, and fourth hair cycles (Figure 3B, blue area). Additionally, the spatial luminescence signal mapping reveals distinct ventral-to-dorsal anagen propagation with features of bilateral symmetry during second (Figure 3C; Appendix 1—figure 6) and third cycles (Figure 3C’), supporting the patterning mechanism predicted by the model (Figure 2F). We also mapped body-wide hair growth patterns on the basis of anagen HF pigmentation between days P0-P55 (Figure 3D–G; Appendix 1—figures 7–12). This analysis confirms ventral over dorsal anagen phase advancement starting from the second cycle and provides the following additional insights: Figure 3 Download asset Open asset Dorsal-ventral HF interactions produce bilateral symmetry. (A) Time-lapse bioluminescence in dorsal and ventral skin of the representative Flash mouse between days P5-P48. Bioluminescent signal is color-coded according to the colorimetric scale on the right. (B) Combined temporal dynamics (from 6 Flash mice) of the bioluminescent signal-based anagen measurements over four hair cycles (days P5-P119). Dorsal skin dynamics are in brown and ventral skin dynamics are in blue. Prominent temporal advancement of ventral over dorsal anagen initiation can be seen during second, third and fourth cycles (light blue areas). This advancement is accompanied by dominant ventral-to-dorsal anagen wave spreading. (C, C’) Mapping of Flash-based anagen reveals ventral-to-dorsal hair growth wave propagation and bilateral pattern symmetry. New anagen areas for each time point are color-coded. Second anagen initiation is shown on panel C, and third anagen initiation on panel C’. Also see Appendix 1—figure 6. (D–G) Hair growth distribution patterns on P17 (D), P21 (E), P39 (F) and P55 (G). Three mice were analyzed at each time point. Inverted whole mount skin samples from representative mice are shown. Schematic pattern maps are provided with color-coded anagen (green), catagen (yellow) and telogen (red) regions. Also see Appendix 1—figures 7–12. (H, I) HF cycling dynamics in chin skin grafts remain faster compared to dorsal skin grafts. After transplantation, first anagen initiated similarly in both chin and dorsal skin grafts, however, second anagen initiated significantly faster in chin grafts. Representative chin and dorsal skin grafts are shown on (H). Combined temporal dynamics of skin grafts in anagen and telogen are shown on (I). Dorsal graft dynamics are in brown and chin graft dynamics are in blue. Temporal advancement of chin over dorsal second anagen initiation is highlighted with light blue color. Also see Appendix 1—figure 13. https://doi.org/10.7554/eLife.22772.005 (i) Ventrally, anagen phase is the shortest in the ‘chin domain’, ending around P10. It is longer in the ‘ventral domain’ proper, ending in the genital area around P14 and in the chest area around P17 (Figure 3D). (ii) Dorsally, anagen is the shortest in the ‘cranial domain’, ending around P14. In the ‘dorsal domain’ proper it ends as a head-to-tail wave between P15-P20 (Figure 3D; Appendix 1—figures 10, 11). (iii) First ventral telogen is shorter than dorsal telogen. Second anagen initiates in the chin and ventral domains already between P21-24 and then spreads toward ventral-dorsal boundaries in form of two bilaterally symmetric waves (Figure 3E; Appendix 1—figure 12). Second anagen also ends faster in the ventral skin, maintaining ventral-dorsal asynchrony and bilateral symmetry (Figure 3F). (iv) Third anagen initiates the fastest in the chin domain, as early as P42 (Figures 3G and 4H). When transplanted onto the back of pigmented SCID mice, chin skin grafts (n = 8) showed faster cycling compared to dorsal skin grafts (n = 8). While first post-transplantation anagen started with similar timing in both chin and dorsal grafts, consecutive anagen started significantly faster in chin grafts (Figure 3H and I; Appendix 1—figure 13). Furthermore, in many instances, grafts induced anagen in the surrounding dorsal host skin. Taken together, these data support that dominant ventral-to-dorsal hair wave spreading drives rapid hair growth pattern evolution and bilateral symmetry. Underlying this behavior are faster hair growth cycle dynamics in chin and ventral HFs, a property that is partially maintained upon skin grafting. Next, we asked if faster hair cycle dynamics in chin and ventral domains correlate with distinct molecular dynamics in putative activators and inhibitors. We performed RNA-seq profiling of whole skin from chin, ventral and dorsal domains at six hair cycle time points: first (aka competent) telogen, early anagen, mid-anagen, late anagen, catagen and early second (aka refractory) telogen. Analysis revealed non-overlapping transcriptomic trajectories of the hair cycle between the three domains (Figure 4A–B’’) and domain-specific expression patterns for multiple putative ac

Recently, we revealed how the underlying molecular mechanism of Bone Morphogenetic Protein (BMP) signaling governs the homeostasis of hair follicle stem cells (hfSCs) in vivo. Here, we focused further on identifying the key interactions of Smads in canonical BMP signaling in hfSC regulation, using our recently developed inducible constitutively active Bmpr1a transgenic system. This allowed us to perform ChipSeq analysis in vivo in hfSCsand receive the predicted target lists for Smads for both p22 and p62 telogen hfSCs. Comparison of our previously published global genes expression where BMP signaling has been inhibited in hfSCs to all overlapped in vivo ChipSeq data between p22 v.s p62 revealed that approximately 72% of genes which previously has been 2 fold up or down regulated after BMP ablation in hfSCs were found directly regulated by BMP/Smads pathway. Previously, we identified Id2 gene as a target of BMP signaling in hfSCs with transcriptional downregulation following targeted BMP inhibition. Here, we confirmed that Id2 gene has been consistently present in our recent ChipSeq analysis both at p22 and p62, therefore, we decided to further test the functional role of the BMP effector, Id2 in hfSCs using an in vivo Id2 gain of function approach. Our data demonstrated that Id2 overexpression in hfSCs results in prolonged telogen and a delay in anagen activation, maintaining stem cells quiescence. Thus, we confirmed the functional role of Id2 gene as a direct target and effector of BMP pathway in hfSCs in vivo.

The G-protein-coupled receptors LGR4, LGR5 and LGR6 are Wnt signaling mediators, but their functions in squamous cell carcinoma (SCC) are unclear. Using lineage tracing in Lgr5-EGFP-CreERT2/Rosa26-Tomato and Lgr6-EGFP-CreERT2/Rosa26-Tomato reporter mice, we demonstrate that Lgr6, but not Lgr5, acts as an epithelial stem cell marker in SCCs in vivo. We identify, by single-molecule in situ hybridization and cell sorting, rare cells positive for Lgr6 expression in immortalized keratinocytes and show that their frequency increases in advanced SCCs. Lgr6 expression is enriched in cells with stem cell characteristics, and Lgr6 downregulation in vivo causes increased epidermal proliferation with expanded lineage tracing from epidermal stem cells positive for Lgr6 expression. Surprisingly, mice with germline knockout of Lgr6 are predisposed to SCC development, through a mechanism that includes compensatory upregulation of Lgr5. These data provide a model for human patients with germline loss-of-function mutations in Wnt pathway genes, including RSPO1 or LGR4, who show increased susceptibility to squamous tumor development.

The control principles behind robust cyclic regeneration of hair follicles (HFs) remain unclear. Using multi-scale modeling, we show that coupling inhibitors and activators with physical growth of HFs is sufficient to drive periodicity and excitability of hair regeneration. Model simulations and experimental data reveal that mouse skin behaves as a heterogeneous regenerative field, composed of anatomical domains where HFs have distinct cycling dynamics. Interactions between fast-cycling chin and ventral HFs and slow-cycling dorsal HFs produce bilaterally symmetric patterns. Ear skin behaves as a hyper-refractory domain with HFs in extended rest phase. Such hyper-refractivity relates to high levels of BMP ligands and WNT antagonists, in part expressed by ear-specific cartilage and muscle. Hair growth stops at the boundaries with hyper-refractory ears and anatomically discontinuous eyelids, generating wave-breaking effects. We posit that similar mechanisms for coupled regeneration with dominant activator, hyper-refractory, and wave-breaker regions can operate in other actively renewing organs.

Inherited germline polymorphisms can cause gene expression levels in normal tissues to differ substantially between individuals. We present an analysis of the genetic architecture of normal adult skin from 470 genetically unique mice, demonstrating the effect of germline variants, skin tissue location, and perturbation by exogenous inflammation or tumorigenesis on gene signaling pathways. Gene networks related to specific cell types and signaling pathways, including sonic hedgehog (Shh), Wnt, Lgr family stem cell markers, and keratins, differed at these tissue sites, suggesting mechanisms for the differential susceptibility of dorsal and tail skin to development of skin diseases and tumorigenesis. The Pten tumor suppressor gene network is rewired in premalignant tumors compared to normal tissue, but this response to perturbation is lost during malignant progression. We present a software package for expression quantitative trait loci (eQTL) network analysis and demonstrate how network analysis of whole tissues provides insights into interactions between cell compartments and signaling molecules.

In this chapter we summarize recent progress in the characterization of different skin stem cell populations and their roles in epidermis, hair follicle, and sweat gland regeneration during homeostasis and injury. We discuss new discoveries in the differentiation capabilities of both embryonic stem cells and induced pluripotent stem cells into functional keratinocytes, with their potential to regenerate skin. We also review how different skin stem cell populations integrate their regulation by intrinsic and extrinsic mechanisms. Furthermore, we discuss a new intrinsic layer of feedback that regulates stem cells in a cell-autonomous manner that, together with newly discovered stem cell niche components, creates a unique environment to host and maintain proper skin stem cell characteristics. These exciting advances in skin stem cells and stem cell niches have helped unravel new concepts in basic stem cell biology and have great translational potential to considerably impact regenerative medicine. Significantly, fully understanding the niche composition and the signals within the microenvironment will help accelerate the design of effective therapeutic strategies to promote stem cell self-renewal and stem cell–based therapies for the regeneration of functional skin with associated skin appendages after injury or in human inherited skin diseases.