Adult hematopoiesis has been studied in terms of progenitor differentiation potentials, whereas its kinetics in vivo is poorly understood. We combined inducible lineage tracing of endogenous adult hematopoietic stem cells (HSCs) with flow cytometry and single-cell RNA sequencing to characterize early steps of hematopoietic differentiation in the steady-state. Labeled cells, comprising primarily long-term HSCs and some short-term HSCs, produced megakaryocytic lineage progeny within 1 wk in a process that required only two to three cell divisions. Erythroid and myeloid progeny emerged simultaneously by 2 wk and included a progenitor population with expression features of both lineages. Myeloid progenitors at this stage showed diversification into granulocytic, monocytic, and dendritic cell types, and rare intermediate cell states could be detected. In contrast, lymphoid differentiation was virtually absent within the first 3 wk of tracing. These results show that continuous differentiation of HSCs rapidly produces major hematopoietic lineages and cell types and reveal fundamental kinetic differences between megakaryocytic, erythroid, myeloid, and lymphoid differentiation.

Introduction

Hematopoiesis is a continuous lifelong process whereby billions of new blood cells are generated every day to maintain essential functions such as oxygen transport (erythrocytes), coagulation (platelets), and immune defense (myeloid cells and lymphocytes). Adult hematopoiesis in mammals occurs primarily in the bone marrow (BM), which comprises a heterogeneous mixture of blood cell types at different stages of differentiation. At the top of the differentiation hierarchy is the hematopoietic stem cell (HSC), a multipotent cell type that can regenerate and sustain multilineage hematopoiesis when transplanted into myeloablated recipients (Eaves, 2015). This unique capacity of HSCs enables BM transplantation, a life-saving procedure that is widely used to treat cancer and other disorders of the blood (Copelan, 2006). On the other hand, aberrant activity of HSCs is thought to contribute to aging-associated abnormalities, anemia, and leukemogenesis (Elias et al., 2014; Adams et al., 2015).

Hematopoiesis is thought to proceed through a hierarchy of stem and progenitor cells with progressively restricted lineage potentials (Shizuru et al., 2005). Thus, true HSCs with long-term reconstitution capacity are thought to give rise to short-term HSCs (ST-HSCs) and/or multipotent progenitors (MPPs), which in turn produce lineage-committed progenitors such as common myeloid and common lymphoid progenitors (CMPs and CLPs, respectively) and finally, cell type–specific progenitors such as granulocyte/monocyte progenitors (GMPs) or megakaryocyte progenitors (MkPs). This HSC-driven hierarchical scheme of hematopoiesis has been established primarily in the transplantation settings, and its relevance to endogenous steady-state hematopoiesis has become a subject of controversy. In particular, it has been argued that HSCs barely contribute to myeloid cells (Sun et al., 2014) or provide a relatively infrequent contribution to hematopoiesis (Busch et al., 2015), emphasizing the putative role of downstream progenitors such as ST-HSCs. In contrast, other recent studies suggested a major sustained contribution of HSCs to steady-state hematopoiesis in mice (Sawai et al., 2016; Yu et al., 2016; Chapple et al., 2018) and humans (Biasco et al., 2016).

Similarly, the precise hierarchy of lineage branching points and the stages of lineage commitment are being hotly debated. For example, the bifurcation of erythroid/megakaryocytic/myeloid versus lymphoid cell fates was originally proposed as the earliest major branching point (Shizuru et al., 2005), as supported recently by the observed clonal divergence of lymphoid and myeloid development in the steady-state (Pei et al., 2017). On the other hand, evidence has been provided for early divergence of megakaryocytic and/or erythroid lineages (Notta et al., 2016; Rodriguez-Fraticelli et al., 2018) and the existence of a common lymphoid-primed MPP (Adolfsson et al., 2005). Furthermore, clonal analyses of stem/progenitor cell output during transplantations or in culture suggested that lineage commitment may occur before the lineage-specific progenitor stages, e.g., in HSCs or MPPs (Naik et al., 2013; Yamamoto et al., 2013; Perié et al., 2015; Lee et al., 2017; Carrelha et al., 2018). This notion has been supported by single-cell RNA sequencing (scRNA-Seq), which revealed preestablished lineage-specific signatures in phenotypically defined CMPs (Paul et al., 2015). On the other hand, progenitor populations with multilineage transcriptional signatures have been detected, consistent with their multipotent nature and ongoing lineage commitment (Drissen et al., 2016; Olsson et al., 2016; Tusi et al., 2018). Collectively, these studies provided fundamental insights into HSC/progenitor differentiation by analyzing its long-term outcomes and/or the static composition of progenitor populations. In contrast, little is known about the sequence of lineage development and the emergence of progenitor populations from HSCs on a real-time scale. Such kinetic information, however, would be critical for the understanding of adult hematopoiesis and of its hierarchical structure.

Recently, we generated a system for inducible genetic labeling of HSCs in vivo, based on the expression of tamoxifen-regulated Cre recombinase-estrogen receptor fusion (CreER) from an HSC-specific transgene. Using this system for long-term lineage tracing, we demonstrated an extensive contribution of adult HSCs to all major hematopoietic lineages except certain embryo-derived cells such as tissue macrophages (Sawai et al., 2016). Here we combined this system with high-dimensional single-cell analysis to characterize the early stages of HSC differentiation. The results provide an unbiased kinetic roadmap of hematopoietic differentiation and reveal major differences in the speed of HSC contribution to different lineages. In particular, they suggest that the early emergence of megakaryopoiesis and the subsequent divergence of erythroid and myeloid development from the relatively delayed lymphopoiesis represent fundamental kinetic hallmarks of steady-state hematopoiesis.

Results

HSCs contribute sequentially to platelets, myeloid cells, and lymphocytes

We have previously described Pdzk1ip1-CreER mice crossed with a Cre-inducible R26tdTomato reporter, in which a single dose of tamoxifen (1 mg) induced the expression of RFP tdTomato (Tom) in a fraction of HSCs (Sawai et al., 2016). The labeled HSCs had a particularly undifferentiated phenotype, showed preferential label retention, and gave rise to other cells within the HSC population, placing them at the top of hematopoietic hierarchy. We further characterized the specificity of HSC labeling induced by an even lower dose of tamoxifen (0.5 mg). Fig. S1 (A and B) shows that median 18% of the phenotypically defined Lin c-Kit+ Sca-1+ (LSK) CD150+ CD48 HSC population became Tom+ 3 d after tamoxifen. A much lower labeling was observed in other immature LSK subpopulations (Kiel et al., 2005; Cabezas-Wallscheid et al., 2014; Pietras et al., 2015), including CD150 CD48 ST-HSCs (median 4%) and CD150+ CD48+ MPP2 (median 2%, with no labeling in half of the mice). No labeling was observed in MPP3/4 or any erythroid or myeloid c-Kit+ Sca-1 progenitors in the majority of animals. Conversely, among the labeled Tom+ cells, median 60% were phenotypic HSCs, 17% were ST-HSCs, and the remaining cells comprised different progenitors that varied between different animals (Fig. S1 B).

To confirm the functionality of labeled cells, we sorted Tom+ and Tom fractions of the LSK CD150+ CD48 HSC population and transferred 100 cells into irradiated recipients, along with 2 × 105 wild-type BM cells. We found that the Tom+ HSCs produced robust multilineage reconstitution, whereas the Tom HSCs produced lower reconstitution of myeloid cells and only poor and variable reconstitution of lymphocytes (Fig. S1, C and D). Thus, a single low dose of tamoxifen allowed labeling of functional HSCs, with the total labeled population comprising predominantly phenotypic HSCs and some ST-HSCs.

Previous long-term tracing of Tom+ HSCs revealed a sustained label accrual in mature hematopoietic cells so that >60% of platelets and myeloid cells became Tom+ within 6–9 mo. Nearly identical kinetics and magnitudes of labeling were observed after animal rederivation into a different facility with a more restricted pathogen-free status (data not shown). Examination of labeling frequencies in the peripheral blood (PB) at earlier time points (4–12 wk) revealed significant differences between lineages (Fig. 1 A). Platelets showed the highest labeling at 4–8 wk. Labeling in granulocytes and monocytes followed and equalized with that of platelets by 12 wk. In contrast, the labeling of all lymphocytes, including natural killer (NK) and B cells, was barely detectable at 4–8 wk and was still very low at 12 wk (Fig. 1 A).

To confirm the differential kinetics of lineage emergence observed in the PB, we analyzed the BM at weekly endpoints after HSC labeling (counting day 3 after tamoxifen administration as a “start”). We designed a panel of antibodies to 17 surface markers, capturing major hematopoietic cell populations including immature erythroid cells (but not megakaryocytes [Mk’s] or platelets). Labeled cells comprised predominantly c-Kithi HSCs/progenitors at the start and 1 wk, with lineage marker-positive (Lin+) c-Kit cells emerging by 2 wk (Fig. 1 B). The stained total BM cells were down-sampled, concatenated, and analyzed using visualization of t-distributed stochastic neighbor embedding (viSNE) dimensionality reduction algorithm (Amir et al., 2013), allowing simultaneous visualization of multiple cell populations. Labeled granulocytes and monocytes emerged at 2 wk and expanded further at 3 wk and beyond (Fig. 1, C and D). Immature CD71+ erythroid cells became detectable at the same time, with further analysis of erythroid cells precluded by their low expression of Tom. In contrast, lymphoid cells including immature B cells were virtually absent at 2–3 wk and became detectable only by 6 wk (Fig. 1, C and D). To confirm that the observed low contribution to lymphocytes was not due to insufficient sampling, we analyzed BM cells without down-sampling (5–6 × 106 total). As shown in Fig. 1 E, the labeling of c-Kit+ progenitors and CD11b+ myeloid cells was evident by 2 wk, and the labeling of B220+ Ly6c+ plasmacytoid dendritic cells (pDCs) became detectable by 3 wk. In contrast, the more abundant B220+ Ly6c B cells (∼1.2–1.5 × 106 cells captured) contained hardly any Tom+ cells at 2–3 wk, but became detectably labeled by 6 wk. These data suggest major differences in the kinetics of lineage specification from the labeled HSC compartment: an initial emergence of platelets is followed by myeloid and erythroid differentiation, whereas lymphoid differentiation is substantially delayed.

HSCs contribute sequentially to megakaryocytic, erythroid and myeloid, and lymphoid progenitors

Having established the differential kinetics of HSC contribution to major lineages, we focused on hematopoietic progenitors within the first 1–3 wk after HSC labeling. The Lin compartment of the BM samples stained with a multicolor panel (Fig. 1, C–E) revealed the presence of megakaryocytic, myeloid, and erythroid progenitors (EryPs) at 2–3 wk, whereas very few lymphoid progenitors were present at these stages (Fig. 2 A). Next, we analyzed tamoxifen-induced Pdzk1ip1-CreER R26tdTomato reporter mice continuously by BM biopsy from the starting point (day 3 after tamoxifen) at weekly intervals. The overall kinetics of labeling was similar to the endpoint analysis, further confirming that BM biopsy does not affect hematopoietic differentiation (Sawai et al., 2016). At 1 wk, the initially labeled HSC population gave rise to Sca-1 c-Kit+ progenitors that were CD150+ CD48lo, encompassing megakaryocytic and EryPs (Fig. S2 A). At 2 wk and even more so at 3 wk, progenitor populations within the LSK compartment, as well as Sca-1 c-Kit+ CD150 CD48hi myeloid progenitors (MyPs), became apparent. To further test the early labeling of megakaryocytic progenitors including CD150+ CD41+ MkPs (Fig. 2 B), we analyzed the BM at the endpoint of 1 wk after labeling. We found that MkPs were labeled, whereas CD150+ CD41 megakaryocytic/erythroid progenitors (MEPs) and CD150 MyPs manifested minimal or no labeling at that point (Fig. 2 C). Notably, MPP2 showed substantial labeling at 1 wk, but MPP3/4 did not, consistent with the proposed role of MPP2 as an intermediate between HSCs and MkPs (Pietras et al., 2015). The preferential labeling of MkPs compared with that of MEPs/MyPs continued up to 4 wk (Fig. 2 C).

To test whether the observed kinetic differences could be caused by tamoxifen itself, we induced the CreER transgene using another estrogen analogue, fulvestrant (ICI 182,780). Whereas tamoxifen acts as an estrogen receptor agonist in hematopoietic cells (Sánchez-Aguilera et al., 2014), fulvestrant is a pure antagonist. As expected from the low affinity for fulvestrant of the CreER T2 variant used in the transgene (Feil et al., 1997), the efficiency of labeling was extremely low. Nevertheless, the overall pattern of labeling at 2–3 wk was similar to that induced by tamoxifen, with the preferential labeling of CD150+ MkPs and MEPs compared with MyPs (Fig. S2 B). Although short of proof, these results suggest that the observed kinetic differences are not caused by tamoxifen activity on the endogenous estrogen receptor.

Next, we extended the analysis to lymphoid progenitors, including Kitlo Flt3+ CD127+ CLPs and their Ly6D+ fraction of committed B cell progenitors (Fig. 2 B; Inlay et al., 2009). While total unlabeled CLP populations were easily detectable at all time points, Tom+ CLPs were essentially absent at 2 wk, but were labeled comparably with other progenitors at 7 wk (Fig. 2, B and D). Also of note, the preferential labeling of MkPs at 2 wk was largely absent at 7 wk, confirming that it reflects kinetic differences rather than a lineage bias.

To correlate the phenotypic analysis with cell function in vitro, we cultured cells in serum-free medium with the cytokines stem cell factor (SCF) and thrombopoietin (TPO), which support HSC/progenitor expansion and/or differentiation into Mk’s (Umemoto et al., 2012). About half of sorted LSK CD150+ CD48 HSCs gave rise to undifferentiated growing colonies (U) that produced a mixture of undifferentiated cells and Mk’s after 2 wk (U+Mk; Fig. 2 E and Fig. S2 C). Other HSCs differentiated into a single Mk or produced multiple Mk’s, as observed previously (Roch et al., 2015). In contrast, LSK CD150 cells produced only U colonies that frequently collapsed by 2 wk, whereas MkPs produced only Mk’s, but no U or U+Mk colonies (Fig. 2 E and Fig. S2 C). We then sorted single Lin Tom+ cells from the BM at the start of labeling (day 3 after tamoxifen) or 1 or 2 wk after labeling and cultured them with SCF + TPO. At the start, the majority of positive wells produced U+Mk colonies, as well as some Mk (Fig. 2 F). At 1 wk, the majority of wells produced Mk’s, with very few U+Mk colonies at the end of culture. At 2 wk, the number of resulting Mk’s was lower, whereas the number of U+Mk colonies was still reduced compared with the start. The megakaryocytic differentiation in these conditions was comparable to that in a full Mk-inducing medium including SCF, TPO, IL-6, and IL-11 (Nishikii et al., 2015), in which bulk cultures of Tom+ cells from all time points produced Mk’s (Fig. S2 D and data not shown). Nearly all U or U+Mk colonies from Tom+ cells at the start of labeling showed efficient replating in liquid medium (Fig. 2 G) and produced colonies of different types in semisolid medium with differentiation-inducing cytokines (data not shown). The replating capacity declined in colonies from 1 wk and was almost absent in colonies from 2 wk after labeling (Fig. 2 G). These results suggest that (1) labeled Tom+ cells at the start, but not at subsequent time points, are comparable to phenotypically defined HSCs in their capacity for expansion and differentiation; and (2) labeled cells at 1 wk show a transient increase in megakaryocytic production. Collectively, the analysis of progenitors is consistent with the tracing of mature lineages and suggests that the labeled HSC population rapidly gives rise to MkPs followed by Ery and MyPs, whereas lymphoid differentiation lags by several weeks.

Initial proliferation of HSC progeny is coupled to megakaryocytic differentiation

We asked how the observed rapid differentiation of HSCs is related to their proliferation; in other words, how many HSC divisions are required to produce lineage-committed progeny. To this end, we crossed Pdzk1ip1-CreER mice to the “Switch” (R26SW) reporter, in which Cre recombination switches the expression of Tom to that of GFP. Treatment of the resulting Pdzk1ip1-CreER R26SW animals with a single dose of tamoxifen (5 mg) yielded the labeling of HSCs that was comparable to the R26Tom reporter (Fig. S3 A). After 1 wk, the label expanded to downstream progenitors, with MkPs showing the highest level of labeling (Fig. S3 A), in agreement with the single-color Tom reporter (Fig. 2 C). Because of its extreme stability in mammalian cells, the fluorescent Tom protein persists after Cre recombination in GFP+ cells until diluted by cell division (Muzumdar et al., 2007), thus representing an endogenous proliferation “tracer.” Indeed, BM biopsy of Pdzk1ip1-CreER R26SW animals at the start of labeling revealed GFP+ cells that were uniformly Tomhi (Fig. 3 A). In the same animals assessed 1 wk later, GFP+ cells diversified into cells with high (Tomhi), low (Tomlo), or negative (Tomneg) Tom expression (Fig. 3 A). Tomhi cells retained the original LSK phenotype, Tomlo cells showed the down-regulation of Sca-1, and Tomneg cells largely lost Sca-1 to acquire the Sca-1 c-Kit+ progenitor phenotype (Fig. 3 B). Notably, these GFP+ Tom progenitors were predominantly CD150+ CD48+/low, and the majority of them represented CD41+ MkPs (Fig. 3 C). Accordingly, whole-mount bone sections revealed multiple GFP+ Tom Mk’s, whose frequency (4.5 ± 1.5, n = 5) was comparable to that of GFP+ MkP (Fig. S3 B). Thus, rapid dilution of Tom “tracer” in the progeny of GFP-labeled HSCs directly correlated with differentiation, which was predominantly directed toward the megakaryocytic lineage.

A detailed analysis of GFP+ cells at 1 wk after HSC labeling showed discreet peaks of Tom fluorescence whose mean fluorescence intensity differed twofold (Fig. 3 D). This was reminiscent of fluorescent dye dilution in proliferation assays, although the resolution was lower due to the low numbers (300–500) of GFP+ cells. Nevertheless, in all analyzed samples (n = 5), we could discern cells with minimal or no Tom dilution (peak 0), up to four peaks with progressive Tom dilution (1–4), and a population that largely diluted Tom (5+; Fig. 3, D and E). While peak 0 remained within the LSK gate, subsequent peaks progressively lost Sca-1 expression (Fig. 3, F and G). Conversely, CD150+ progenitors became prominent starting at peak 2, with the fraction of CD150+ CD41+ MkP progressively increasing. A small fraction of CD150 MyPs was present in peaks 2–5+, and c-Kit progenitors were detectable in peak 5+. Although it is possible that labeled cells other than HSCs (mainly ST-HSCs; Fig. S3 A) contribute to the observed rapid proliferation, this is not compatible with the limited megakaryocytic potential of ST-HSCs (Pietras et al., 2015). These results are consistent with the observed differential kinetics of HSC contribution to different lineages and suggest that (1) some cells within the labeled HSC population undergo differentiation into MkP that takes only two to three cell divisions, and (2) the differentiation of labeled HSCs into c-Kit+ MyPs trails by one to two divisions, whereas c-Kit progenitors require extensive proliferation (five or more divisions).

scRNA-Seq reveals stepwise HSC differentiation in vivo

To characterize the molecular events during the differentiation of endogenous HSCs, we analyzed the labeled Tom+ HSCs and their progeny that appeared over time by global mRNA-Seq. Initially, we performed RNA-Seq on bulk populations of 1–2 × 103 Tom+ cells sorted from individual animals at the start of labeling or after 2 wk (data not shown). The genes up-regulated at 2 wk mainly comprised markers of erythroid (Klf1 and Tfrc/CD71) and myeloid differentiation, including monocytes/macrophages (F13a1 and Ly6c2), neutrophilic granulocytes (Hp and Alas1), and mast cells/basophils (Prss34 and Mcpt8). In contrast, lymphoid genes (e.g., Rag1, Rag2, Dntt, Blk, and Il7r/CD127) were virtually undetectable in all samples. We then performed scRNA-Seq on Lin Tom+ cells sorted from individual animals at the start of labeling or 1 or 2 wk after labeling. By phenotype, labeled cells at the start comprised primarily HSCs; at 1 wk, HSCs plus MkPs and MEPs; and at 2 wk, a variety of Kit+ and Kit progenitors (Fig. S4 A). In total, 188, 107, and 283 cells from the respective time points passed quality control and were pooled and analyzed by unsupervised clustering and dimensionality reduction using the Seurat software (Satija et al., 2015). The resulting cell pool segregated into five broad clusters (Fig. 4 A) that expressed genes characteristic of either HSCs (cluster 0), MkPs (cluster 1), MPPs (cluster 2), MyPs (cluster 3) or EryPs (cluster 4; Fig. 4 B and Fig. S4 B). The expression of genes encoding key transcription factors (TFs; Fig. 4 C) and the genes encoding surface markers (Fig. S4 C) was also consistent with the proposed identities of the clusters. As expected, the expression of proliferation-related genes was observed mainly in lineage-restricted progenitors (clusters 1, 3, and 4; Fig. S4 D). Of note, virtually no lymphoid lineage-specific markers or TFs could be detected (Fig. 4 C and Fig. S4 B).

Importantly, cell populations from different time points were composed of different clusters (Fig. 4 D). As expected, >80% of cells at the start of labeling were from the HSC cluster. Cells from week 1 contained predominantly the HSC and MkP clusters, and cells from week 2 also contained large fractions of MPPs and MyPs. Conversely, the MPP, EryP, and MyP clusters were composed largely of cells from the 2-wk time point (Fig. 4 E). These data collectively suggest that (1) the contribution of labeled HSC population to MkPs occurs rapidly (within 1 wk) and before the emergence of MPPs defined by gene expression, as well as of MPP3/4 defined by phenotypic analysis (Fig. 2 C and Fig. S3 A); (2) HSCs yield MPPs, as well as definitive erythroid and myeloid progeny, within 2 wk; and (3) lymphoid progeny are virtually undetectable at these stages, confirming a major delay of lymphoid relative to myeloid differentiation.

HSCs rapidly generate a broad spectrum of progenitors

The data described above showed the emergence of multipotent and lineage-restricted progenitors from labeled HSCs after 2 wk of tracing. To better resolve this earliest stage of multilineage differentiation, we sequenced single Lin Tom+ cells from two additional animals at 2 wk after labeling. The results were merged with the same time point from the first experiment, and the pooled 647 cells were analyzed using Seurat as above. As expected, this analysis revealed four broad cell clusters corresponding to MyP (cluster 0), HSC + MPP (cluster 1), MkP (cluster 2), and EryP (cluster 3; Fig. 5, A and B). These clusters contained cells from both experiments, suggesting that the clustering by transcriptional profiles overcame batch effects (although the latter were detectable within the clusters, particularly clusters 1 and 2; Fig. S5 A). All clusters showed prominent expression of the expected lineage-specific genes (Fig. 5 C) and surface markers (Fig. S5 B). Even with the increased number of cells, little evidence of lymphoid development could be detected (note only two cells expressing Rag2 in Fig. 5 C).

Given the simultaneous emergence of EryPs and MyPs at this stage, we sought a potential progenitor with expression features of both lineages. The vast majority of cells in the EryP cluster (cluster 3) expressed erythroid-specific genes such as Tspo2, Blvrb, and Car1; however, only a fraction of EryP cells showed hallmarks of advanced erythroid differentiation such as the expression of glycophorin A (Gypa) and hemoglobin (Hbb) genes (Fig. 6 A). Indeed, the EryP cluster could be divided into two distinct subclusters, EryP-A and EryP-B (Fig. 6 B), each of which contained cells from both merged experiments (Fig. S5 A). The EryP-B subcluster expressed erythroid-specific TF (Klf1 and Gata1), marker genes (Cldn13, RHD, Ctse, and the Hbb genes), and proliferation signature genes (Fig. 6, C–E), likely representing committed progenitors undergoing erythroid differentiation. In contrast, EryP-A cells expressed not only erythroid TFs, but also TFs typical of other cell types including progenitors (Hlf and Meis1), Mk’s (Gata2), and myeloid cells (Spi1). Furthermore, EryP-A cells lacked the above-mentioned erythroid genes, but expressed multiple genes characteristic of granulocytes (Csf2ra and Clec4d), particularly basophils and mast cells (Cpa3 and Csrp3; Fig. 6 D and data not shown). Finally, EryP-A cells lacked the proliferation signature genes expressed by EryP-B cells (Fig. 6 E). Thus, the expression profile of these cells has features of multiple lineages including erythrocytes and basophils, resembling a previously described dual progenitor of erythrocytes and basophils/eosinophils (Drissen et al., 2016). Indeed, rare cells with the same expression profile could be detected in the previous static analyses of steady-state progenitor populations (Fig. S5 C and data not shown). These data suggest that the earliest steps of HSC differentiation into erythroid and myeloid lineages yield cells with expression features of either or both lineages.

Rapid specification of myeloid cell types during HSC differentiation

Finally, we analyzed the most abundant myeloid cluster (cluster 0) for potential heterogeneity and emergence of distinct myeloid cell type expression programs. Although most cells in this cluster expressed general markers of myeloid differentiation (e.g., Plac8 and Mpo), the expression of more cell type–specific myeloid genes was dichotomous (e.g., monocyte-enriched Elane and Ly6c2) or very rare (e.g., basophil-enriched Slc18a2 and Csrp3; Fig. 7 A). Accordingly, the MyP cluster was heterogeneous and could be separated into at least three distinct subpopulations containing cells from both experiments (Fig. 7, B–D). The more abundant MyP-A subcluster featured genes characteristic of monocytes and/or macrophages (Ly86 and Csf1r); the MyP-B subcluster featured granulocyte-specific genes such as Hdc and S100a8; and the small and tight MyP-C subcluster expressed multiple basophil/mast cell–specific genes (Mcpt8, Prss34, Csrp3, and Cpa3; Fig. 7 B). The expression of lineage markers was consistent with that of TF, with MyP-A expressing monocyte/dendritic cell regulators Irf8 and Klf4, MyP-B expressing granulocytic TFs Cebpe and Gfi1, and MyP-C expressing Gata2 that is required for basophil development (Fig. 7 C). Thus, the earliest steps of myeloid differentiation comprise cells with features of mononuclear, granulocytic neutrophilic, and granulocytic basophilic cell fates (Fig. 7 D).

It is an open question whether myeloid differentiation proceeds directly from HSCs/MPPs into predetermined cell fates or involves intermediate MyPs undergoing cell fate choice. One example of the latter is monocytic versus granulocytic lineage bifurcation, which was shown to involve rare MyPs that coexpress the key TFs IRF8 and GFI1 (Olsson et al., 2016). Despite the overall segregation of Irf8 and Gfi1 expression between MyP-A and MyP-B subclusters, rare cells expressing both factors could be detected (Fig. 7 E), confirming the early emergence of intermediate Irf8+ Gfi1+ progenitors during hematopoietic differentiation.

Another example of myeloid cell fate decision is the specification of interferon-producing pDCs and antigen-presenting classical dendritic cells (cDCs) from common progenitors (Merad et al., 2013; Murphy et al., 2016). Both pDC and a distinct cross-presenting cDC subset termed cDC1 require IRF8, and Irf8 was indeed abundantly expressed by nearly all MyP-A cells (Fig. 7 C). The pDC specification is regulated by the E protein TF TCF4 (E2-2), which is expressed in HSCs and all progenitors and is further up-regulated during pDC differentiation (Grajkowska et al., 2017). TCF4 represses the expression of ID2, an E protein antagonist that blocks pDC and facilitates cDC1 development (Ghosh et al., 2010). The MyP-A subcluster contained a distinct group of cells that expressed high level of Tcf4, as well as of another pDC-promoting TF Bcl11a (Ippolito et al., 2014) and of pDC-specific markers SiglecH (Blasius et al., 2006) and CD300c (Fig. 7 F; Kaitani et al., 2018). Conversely, multiple cells adjacent to this Tcf4hi group expressed Id2, as well as a common dendritic cell progenitor marker Cx3cr1 (Fig. 7 F). Only few cells expressed the cDC1-specific marker Clec9a and TF Batf3, consistent with incipient cDC differentiation at this stage. Notably, cells that coexpress Tcf4 and Id2 could be readily identified within this group (Fig. 7 G). Collectively, these data reveal myeloid cells expressing genes for mutually antagonistic TFs that drive specification of alternative cell fates such as IRF8/GFI1 (monocytic/granulocytic) and TCF4/ID2 (pDC/cDC1). The emergence of these intermediate cell states at the onset of myeloid differentiation is consistent with ongoing cell fate choices in early MyPs.

Discussion

We took advantage of the efficiency and precise time frame of inducible HSC labeling in Pdzk1ip1-CreER animals to characterize the initial steps of hematopoietic differentiation in the BM of unperturbed adult animals. The majority of labeled cells comprised phenotypically defined HSCs, although some ST-HSCs and a small variable number of progenitors were labeled as well. Such an extent of phenotypic heterogeneity can be expected, given that the HSCs are heterogeneous and their functions may not strictly segregate with phenotypes defined by nonfunctional markers such as CD150. Notably, the labeled cells were fully functional in the transplantation assay (more so than the unlabeled HSC fraction), resembled phenotypically defined HSCs in the growth and differentiation capacity in vitro, and sustained multilineage hematopoiesis in the steady-state. Thus, notwithstanding the labeling of some phenotypic non-HSCs, our system labeled the cell population at the top of hematopoietic hierarchy. Another caveat applying to all CreER-based inducible tracing systems is the potential effect of tamoxifen on hematopoietic differentiation. Notably, we used the single dose of tamoxifen (0.5 mg for the Tomato reporter) that is substantially lower than the cumulative doses used in other CreER-based HSC-labeling studies (5 mg in Busch et al., 2015; 12.5 mg in Chapple et al., 2018) or shown to affect hematopoiesis (7–28 mg; Sánchez-Aguilera et al., 2014). Furthermore, high doses of tamoxifen preferentially impair progenitors such as MPPs (Sánchez-Aguilera et al., 2014), which are not initially labeled in our system. Finally, a similar pattern of differentiation was observed with a functionally different estrogen analogue, suggesting that the effects of tamoxifen are likely minimal.

We found that labeled HSCs rapidly gave rise to the differentiated progeny including a broad spectrum of megakaryocytic, erythroid, and myeloid differentiation within 2 wk after labeling. These data strongly support a major ongoing contribution of HSCs to multilineage hematopoiesis in the steady-state, with the bulk of HSCs in the body contributing during the lifetime and only a minor HSC fraction kept in “reserve” (Bernitz et al., 2016; Sawai et al., 2016). Using single-cell approaches to analyze the contribution of HSCs on a real-time scale, we observed a fundamentally different kinetics whereby different lineages emerge from HSCs. The megakaryocytic lineage including MkPs and their products, Mk’s and platelets, emerged before other lineages within 1 wk of HSC labeling. This rapid Mk contribution is unlikely to be caused by labeled non-HSCs such as ST-HSCs, which have a limited Mk differentiation potential in vivo and in vitro (Pietras et al., 2015). These observations are consistent with the early divergence of Mk lineage potential in the hematopoietic hierarchy, as supported by (1) the transcriptional “priming” of HSCs toward the Mk lineage (Kent et al., 2009; Guo et al., 2013; Sanjuan-Pla et al., 2013); (2) an Mk lineage-specific contribution of HSCs upon transplantation (Carrelha et al., 2018), in inflammatory conditions (Haas et al., 2015), or in human HSC/progenitor culture (Notta et al., 2016); and (3) long-term clonal divergence of Mk from other lineages in steady-state hematopoiesis (Rodriguez-Fraticelli et al., 2018). Using an endogenous proliferation tracer, we estimated that HSCs can produce phenotypic Mk progeny after only two to three cell divisions, whereas MyPs require more extensive proliferation. This rapid proliferation-coupled differentiation may be mediated by a fraction of bona fide HSCs that preferentially contributes to the Mk lineage while maintaining multipotency (Sanjuan-Pla et al., 2013; Carrelha et al., 2018). Alternatively, a population of committed MkP may be present within the phenotypic HSCs (Roch et al., 2015); this possibility could not be ruled out by our polyclonal labeling method and remains to be tested in future studies. However, ∼50% of HSCs showed no or minimal label dilution at 1 wk, and their slow division is likely to maintain the HSC pool. Thus, the megakaryocytic lineage priming of some cells within the HSC population may enable rapid megakaryocytic differentiation without major transcriptional and epigenetic changes, whereas the latter require multiple cell divisions in the case of other lineages.

The erythroid and myeloid progeny of HSCs was virtually undetectable at 1 wk, but became prominent at 2 wk after labeling. This simultaneous emergence raises the question of whether it proceeds through a classically defined CMP with erythroid potential or diverges at the higher hierarchical levels, as argued recently (Paul et al., 2015; Perié et al., 2015). These scenarios are not mutually exclusive, and indeed, our scRNA-Seq data reveal a continuum of progenitors with lineage-specific and intermediate expression profiles. The observed continuous spectrum of differentiation is consistent with other recent studies of HSCs and progenitors using scRNA-Seq (Nestorowa et al., 2016; Rodriguez-Fraticelli et al., 2018; Tusi et al., 2018). In particular, cells with an erythroid expression program and features of certain myeloid programs, particularly those of basophils, were readily detectable. Consistent with their definition as uncommitted progenitors, these cells lacked additional markers of either lineage, but expressed genes and TFs characteristic of HSCs and MkP. The expression profile of these cells resembles that of erythroid/mast cell progenitors described by Drissen et al. (2016), and similar cells are also detectable in other scRNA-Seq experiments (Nestorowa et al., 2016; Olsson et al., 2016). Notably, these cells appear to represent a relatively minor fraction of steady-state progenitor populations (Drissen et al., 2016), whereas they were relatively abundant in our analysis of recently emerged progenitors. Its relative predominance at early stages of HSC differentiation suggests this progenitor population as an important intermediate step in erythroid and/or myeloid development.

Similar to the bifurcation of erythroid and myeloid lineages, within the latter, we observed progenitors with both myeloid cell type–specific and intermediate expression programs. Cells with expression features of monocytes/macrophages, neutrophils, basophils, and dendritic cells (particularly pDCs) were already predominant at the onset of myeloid differentiation, revealing a rapid HSC specification into these myeloid cell types. At the same time, cells with an intermediate profile expressing antagonistic TFs could be readily detected. These include the previously described progenitors coexpressing Gfi1 and Irf8, which were shown to undergo fate choice toward granulocytic or monocytic lineages (Olsson et al., 2016). Another example is the appearance of progenitors coexpressing Tcf4 and its antagonist Id2, which are expressed in a mutually exclusive manner in mature pDCs and cDCs. The existence of these cells was predicted based on the genetic basis of dendritic cell development (Grajkowska et al., 2017), but they have not been described previously. These results confirm that the diversification of myeloid cell types proceeds rapidly and, at least in some cases, involves clearly detectable intermediate states.

Despite our system’s ability to detect very rare cells by flow cytometry and scRNA-Seq (such as basophil/mast cell progenitors), lymphoid differentiation was virtually undetectable at 2–3 wk after HSC labeling. This delay does not reflect a lineage bias of labeled HSCs, because they eventually produce lymphocytes and the labeling of lymphoid progenitors approaches that of MyPs (Sawai et al., 2016). Neither does it reflect the lag period required to undergo V(D)J recombination, because even the earliest preceding steps such as the expression of CD127 (IL-7R), were undetectable. Thus, the delay of lymphopoiesis compared with other lineages appears to represent a built-in kinetic feature, the mechanism of which remains to be elucidated. One possibility is that lymphoid lineage commitment represents a rare (possibly stochastic) event that diverts an MPP from the default myeloid differentiation pathway. Another possibility is that lymphoid and erythroid/myeloid differentiations originate from different sets of HSC clones, e.g., those that demonstrate myeloid or lymphoid bias upon transplantation (Sieburg et al., 2006; Dykstra et al., 2007; Carrelha et al., 2018). This possibility is supported by the different clonal composition of adult lymphoid and myeloid compartments upon clonal cell tracing from the embryo (Pei et al., 2017). Irrespective of the mechanism(s), the fundamental kinetic dichotomy of erythroid/myeloid and lymphoid development should be incorporated into hierarchical models of hematopoiesis; as a minimum, it argues against a simple continuous divergence of myeloid and lymphoid development at a fixed cellular stage.

Collectively, our results suggest a preliminary model describing the kinetics of HSC-driven hematopoiesis in the steady-state (see the graphical abstract). The earliest emergence of Mk’s and platelets is consistent with their essential role, short lifespan (∼5 d for platelets), and ancient evolutionary origins that precede the immune system. Indeed, platelet reduction below a certain threshold represents a danger of the organism’s demise within minutes and must be prevented by hard-wired features of hematopoiesis such as rapid HSC-to-Mk differentiation. The contribution of HSCs to erythroid and myeloid cells is also rapid but delayed relative to megakaryopoiesis. Accordingly, the reduction of myeloid cells or erythrocytes would impair the organism’s survival on the scale of days, and its danger is further mitigated by myeloid cell redundancy (e.g., with tissue macrophages) and by the long lifespan of erythrocytes. Finally, the significantly delayed contribution to lymphocytes may reflect the role of the latter in long-term survival of the organism, as well as the recent evolutionary origin of the adaptive immune system. It is also consistent with long lifespans of peripheral lymphocytes and with their location in lymphoid organs, which precludes their rapid simultaneous loss. The proposed kinetic scheme of hematopoiesis may provide a useful starting point for the integration with cross-sectional single-cell analyses to yield a comprehensive molecular and cellular picture of steady-state hematopoiesis.

Materials and methods

Animals

All animal studies were performed according to the investigator’s protocol approved by the Institutional Animal Care and Use Committees of Columbia University or New York University School of Medicine. Wild-type C57BL/6 mice and C57BL/6 mice congenic for CD45.1 (B6.SJL-PtprcaPepcb/BoyCrl) were obtained from Taconic and Charles River, respectively. Pdzk1ip1-CreER transgenic mouse strain was previously described (Sawai et al., 2016) and was crossed to the Cre-inducible reporter strains R26Tom/Tom (B6;129S6-Gt(ROSA)26Sortm9(CAG-tdTomato)Hze/J; Madisen et al., 2010) or R26SW/SW (Gt(ROSA)26Sortm4(ACTB-tdTomato,-EGFP)Luo/J; Muzumdar et al., 2007; both from the Jackson Laboratories). All mice were crossed to C57BL/6 background for >10 generations. 8–12-wk-old mice of either sex that were either heterozygous or homozygous for the reporter alleles were used in the studies, as no major differences in the efficiency or specificity of labeling were observed between them. To induce recombination in the reporter mice, a single dose of 0.5–1 mg (for R26Tom) or 5 mg (for R26SW) tamoxifen was administered by gavage. To induce recombination using fulvestrant, a single dose of 0.5–1 mg fulvestrant (Sigma) was administered by i.p. injection. Day 3 after tamoxifen or fulvestrant administration was considered a starting point of labeling in all experiments. PB collection and BM biopsy were performed as previously described (Sawai et al., 2016). In brief, blood (∼0.1 ml) was collected by submandibular vein puncture with a sterile disposable lancet. BM biopsy was performed on isoflurane-anesthetized animals by inserting a 281/2-G insulin syringe needle into the joint surface of the femur through the patellar tendon and then into the bone cavity. Up to 20 µl of the BM suspension (∼106 cells) was collected.

Flow cytometry and cell sorting

For analysis of Tom expression in the PB, cells were stained with directly conjugated antibodies against CD45, B220, TCRb, CD11b, Gr-1, NK1.1, and Ter119. Platelets were identified in a separate staining using antibodies against CD41, CD150, and Ter119. The gating scheme for the identification of cell populations was described previously (Sawai et al., 2016). For stem/progenitor cell identification in the BM, cell suspensions were subjected to red blood cell lysis and stained with antibodies against CD150, CD48, c-Kit, Sca-1, and CD41 in addition to a cocktail of lineage markers (CD11b, CD19, B220, TCRb, CD4, CD8, DX5, Gr-1, and Ter119). Lymphoid progenitors were identified using antibodies against CD127, CD135, CCR9, Ly6D, c-Kit, Sca-1, and a mixture of lineage markers including B220, CD11b, CD19, TCRb, Gr-1, Ter119, and NK1.1. Samples were acquired on an Attune NxT flow cytometer (Thermo Fisher Scientific). For high-dimensional flow cytometry experiments, cells were stained with antibodies against CD150, CD48, CD41, c-Kit, Sca1, CD34, CD71, CD127, CD201, FcγR III/II, Flt3, CD45.2, Ly6C, CD11c, B220, CD11b, and Ter119. Sample acquisition was performed on ZE5 Cell Analyzer (Propel Labs). A detailed list of the antibodies used for various experiments can be found in Table S1.

For sorting, total BM cells were isolated from femurs, tibias, ilia, and humeri of Pdzk1ip1-CreER R26Tom animals on day 3, 10, and 17 following tamoxifen treatment. Following red blood cell lysis, cells were stained with biotinylated antibodies to lineage markers and negatively selected using streptavidin magnetic beads and MACS Cell Separation Columns (Miltenyi Biotech) according to the manufacturer’s protocol. The resulting lineage-depleted cells were stained for HSC/progenitor markers and sorted on FACS Aria IIu SORP cell sorter (BD Biosciences) using tdTomato fluorescence as the sole sorting criterion. For bulk RNA isolation and sequencing, ∼2,000 Tom+ cells were sorted from each sample. For single cell RNA sequencing, Tom+ single cells were sorted into 96-well plates and subsequently processed for RNA isolation and sequencing. For cell culture, Tom+ cells were sorted into 96-well plates as single cells (for HSC expansion) or at 1,000–2,000 cells/well (for megakaryocytic differentiation).

Flow cytometry analysis

Conventional flow cytometry analysis was performed using the FlowJo software (FlowJo, LLC). Cell populations in the PB were defined as previously described (Sawai et al., 2016). Cell populations in the BM were defined as follows: HSC, LSK CD150+ CD48; ST-HSC, LSK CD150 CD48; MPP2, LSK CD150+ CD48+; MPP3/4, LSK CD150 CD48+; MyP, Lin c-Kit+ Sca-1 CD150 CD41; MEP, Lin c-Kit+ Sca-1 CD150+ CD41; MkP, Lin c-Kit+ Sca-1 CD150+ CD41+; CLP, Lin Kitlo/− Sca-1 Flt3+ CD127+ Ly6D; and CLP', Lin Kitlo/− Sca-1 Flt3+ CD127+ Ly6D+.

High-dimensional flow cytometry data were analyzed using the viSNE algorithm (Amir et al., 2013) as implemented in FCS Express (De Novo Software). To facilitate the analysis of rare labeled cells, all Tom+ cells from each acquired BM sample (4–6 × 106 cells) were gated, and the entire sample was then down-sampled to 5 × 104 cells. The resulting entire Tom+ cell population and the down-sampled total cell population from each sample were concatenated, and the resulting file was analyzed by viSNE using hyperlog scaling of fluorescence parameters, 1,000 iterations, and Barnes-Hut approximation of 0.5. The cells were plotted using the resulting values of t-distributed stochastic neighbor embedding (tSNE) parameters. For gated cells (Fig. 1 D and Fig. 2 A), all events were sampled; for the entire file (Fig. 1 C), 50,000 events were sampled, and tSNE was estimated for the unsampled events. In the analysis of cell lineages (Fig. 1), all fluorescence parameters except tdTomato were used. In the analysis of progenitors (Fig. 2 A), the gating parameters (CD11b, Ter119, and B220) were excluded from viSNE.

The populations in Fig. 1 C were defined as follows: mature lineages, granulocytes (blue, CD11b+ CD48lo), monocytes (light blue, CD11b+ Ly6chi CD48hi), pDC (brown, B220+ Ly6c+), T/NK cells (orange, CD11b B220 CD45hi), B cells (B220+ Ly6c [including pro–B cells], light yellow, CD150lo CD127+; immature B cells, yellow, CD150lo CD127; mature B cells, dark yellow, CD150hi CD127+), and erythroid (green, CD71+ Ter119+; dark green, CD71 Ter119+); lineage (B220/CD11b/Ter119)-negative progenitors: LSK (pink, c-Kit+ Sca-1+), MkP (purple, c-Kit+ Sca-1 CD150+), myeloid (turquoise, c-Kit+ Sca-1 CD150 CD41), and erythroid (light green, CD71+ Ter119). The populations in Fig. 1 D were defined the same way, except erythroid (green) encompassed all CD71+ cells. The populations in Fig. 2 A were gated as lineage (B220/CD11b/Ter119) negative and defined as for conventional flow analysis described above, with the following additions: EryPs (CD71+) and monocyte/dendritic cell progenitors (Ly6c+).

Cell transplantations

Pdzk1ip1-CreER+ R26Tom animals were treated with 1 mg tamoxifen and sacrificed 3 d later, and lineage-depleted cells were prepared and stained as above. Tom+ and Tom fractions of the Sca1+ cKit+ CD48 CD150+ populations were sorted, and 100 sorted cells were injected i.v. into lethally irradiated CD45.1+ congenic B6.SJL animals. Total BM cells from untreated B6.SJL mice were injected in parallel (2 × 105 per recipient) to ensure recipient survival and provide competitor BM. The fraction of donor-derived CD45.2+ cells among the major cell types in the PB was measured monthly.

Cell culture

Sorted single Tom+ cells were cultured in 96-well plates containing serum-free medium StemSpan SFEM II (STEMCELL Technologies) supplemented with 50 ng/ml recombinant murine SCF and 50 ng/ml TPO (Peprotech). Cultures were examined on days 3–4, and wells with living cells were ascertained by brightfield and (where appropriate) fluorescent microscopy for Tom expression. Wells were scored throughout the culture period based on morphology as either U colonies, differentiating into Mk’s, or a mixture of undifferentiated cells and Mk’s (U+Mk). After 7 d in culture, an equal number of randomly selected U and U+Mk colonies was replated whereby half of the well was transferred into a 96-well plate containing a fresh batch of the medium described above and the other half was cultured in semisolid MethoCult GF M3434 (STEMCELL Technologies) supplemented with TPO to assay for CFU potential.

To assess megakaryocytic differentiation, bulk sorted Tom+ cells were transferred to 24-well plates at a density of 1,000–2,000 cells/well and cultured in the same liquid medium as above with the addition of 20 ng/ml IL-6 and 20 ng/ml IL-11 (Peprotech). Brightfield and fluorescent microphotographs were acquired using BZ-X800 microscope (Keyence) with manufacturer’s software.

Whole-mount immunostaining

Whole-mount tissue preparation of sternal bones was performed as described (Kunisaki et al., 2013). In brief, sternal bones from tamoxifen-induced Pdzk1ip1-CreER R26SW were isolated, transected into two to three fragments and sagittally bisected to expose the BM. The tissues were fixed in 4% paraformaldehyde for 1 h at 4°C, stained with biotinylated anti-CD48 (HM48-1; BioLegend), rat anti-mouse CD41 (MWReg30; BioLegend), and biotinylated lineage markers (anti-Ter119 [TER119], anti-B220 [6B2], anti-Gr1 [RB6-8C5], anti-CD3e [500A2], anti-CD11b [HM1/70]; eBiosciences). Secondary reagents included Alexa Fluor 647–conjugated goat anti-rat IgG and Alexa Fluor 488–streptavidin conjugate (Life Technologies). Tissues were imaged on a 35-mm glass-bottom dish using a Zeiss LSM 710 confocal microscope.

RNA isolation and library preparation

For bulk RNA-Seq, cells were lysed using RLT buffer (Qiagen) and DNA/RNA was separated from lysed cell debris using AMPure beads. Immediately after isolation, a modified Smart Seq 2 protocol was used to perform reverse transcription, cDNA synthesis and amplification (Bracken, 2018). In brief, samples were resuspended in a mixture containing 5× Maxima Transcription Buffer, dNTPs, and SUPERase RNase inhibitor (Thermo Fisher Scientific). A unique 12-bp barcode was added to each sample to allow for pooling and multiplexing. Reverse transcription and cDNA amplification were performed with Maxima H Minus Reverse transcription (Thermo Fisher Scientific) and KAPA HotStart Ready Mix (Kapa Biosystems), respectively. Sequencing libraries were made using Nextera XT sample preparation kit (Illumina). Before sequencing with Hiseq 2500, the libraries were quantified using Qubit and BioAnalyzer to check their quality.

For scRNA-Seq, cells were sorted into 96-well plates containing a mix with 5× Maxima Transcription Buffer, dNTPs, and NXGen RNase inhibitor (Lucigen). Once sorted in wells, plates were spun down and put on dry ice for 5 min to lyse the cells. A modified version of the SMARTseq2 protocol (Picelli et al., 2014) was performed using SuperScript II Reverse transcription and KAPA HotStart Ready Mix to generate and amplify the cDNA for each cell (Bracken, 2018). Each cell/transcript was assigned a unique 12-bp barcode that facilitated individual cell cDNA pooling at the Nextera XT library preparation step. QC analysis was done using Qubit and BioAnalyzer kits.

Gene expression analysis

For all sequencing results, raw sequencing reads were processed using DropSeq tools v1.12 (Macosko et al., 2015). To obtain a reads/cell/gene count table, reads were mapped to mouse GRCm38.84 reference genome. scRNA-Seq analysis was performed using Seurat (Satija et al., 2015). Raw gene counts were log-normalized to correct for the differences in sequencing depth between single cells. Cells with unique gene counts <1,400 or >9,000 and/or cells with >8% of mitochondrial genes were removed. Latent variable regression was performed, regressing out the following variables: the percentage of mitochondrial genes and the ratio of unique molecular identifier counts per gene counts. Variable genes were detected based on average expression and dispersion using the following parameters: x.low.cutoff = 0.16, x.high.cutoff = 3, and y.cutoff = 0.7. The resulting variable genes were subjected to principal component analysis, and top principal components were selected based on the inflection point of the elbow plot of P values and used for cell clustering. The dataset in Fig. 4 yielded 1,990 variable genes out of total 22,971 detected genes across 578 cells, and principal components 1–6 were used for clustering. The dataset in Fig. 5 yielded 2,073 variable genes out of total 34,615 genes across 647 cells, and principal components 1–9 were used for clustering. Further software-based subclustering was precluded by insufficient numbers of cells in each cluster. Therefore, the subclusters within the myeloid and erythroid clusters (Fig. 6 B and Fig. 7 D) were identified manually based on the heterogeneous expression of key subset-specific genes including TFs.

Data availability

scRNA-Seq data have been deposited in the Gene Expression Omnibus (GEO) repository under the accession no. GSE120239.

Statistical analysis

Normal distribution of data were not assumed and statistical significance of differences between experimental groups was determined by nonparametric Mann-Whitney test using the Prism software (GraphPad). Differences were considered significant for P values <0.05 (*), <0.01 (**), and <0.001 (***).

Online supplemental material

Supplemental material includes five figures and one table. Fig. S1 presents characterization of labeled cell population in induced Tom reporter mice. Fig. S2 provides additional characterization of HSC contribution to hematopoietic progenitors in Tom reporter mice. Fig. S3 describes the tracing of HSC differentiation using the Switch reporter mice. Fig. S4 provides additional analysis of HSC differentiation at 0–2 wk by scRNA-Seq. Fig. S5 describes additional analysis of differentiating HSC progeny at 2 wk after HSC labeling. Table S1 lists antibody conjugates used for cell staining.

Acknowledgments

This paper is supported by the National Institutes of Health grants AI072571, AG049074 and AI115382 (to B. Reizis), and AI100853 (A. Rashidfarrokhi).

The authors declare no competing financial interests.

Author contributions: S. Upadhaya, C.M. Sawai, E. Papalexi, A. Rashidfarrokhi, and G. Jang performed experiments, analyzed, and interpreted results. P. Chattopadhyay, R. Satija, and B. Reizis analyzed and interpreted results. S. Upadhaya and B. Reizis wrote the original manuscript. All authors reviewed and edited the manuscript.

References

References
Adams
,
P.D.
,
H.
Jasper
, and
K.L.
Rudolph
.
2015
.
Aging-Induced Stem Cell Mutations as Drivers for Disease and Cancer
.
Cell Stem Cell.
16
:
601
612
.
Adolfsson
,
J.
,
R.
Månsson
,
N.
Buza-Vidas
,
A.
Hultquist
,
K.
Liuba
,
C.T.
Jensen
,
D.
Bryder
,
L.
Yang
,
O.J.
Borge
,
L.A.
Thoren
, et al
2005
.
Identification of Flt3+ lympho-myeloid stem cells lacking erythro-megakaryocytic potential a revised road map for adult blood lineage commitment
.
Cell.
121
:
295
306
.
Amir
,
A.D.
,
K.L.
Davis
,
M.D.
Tadmor
,
E.F.
Simonds
,
J.H.
Levine
,
S.C.
Bendall
,
D.K.
Shenfeld
,
S.
Krishnaswamy
,
G.P.
Nolan
, and
D.
Pe’er
.
2013
.
viSNE enables visualization of high dimensional single-cell data and reveals phenotypic heterogeneity of leukemia
.
Nat. Biotechnol.
31
:
545
552
.
Bernitz
,
J.M.
,
H.S.
Kim
,
B.
MacArthur
,
H.
Sieburg
, and
K.
Moore
.
2016
.
Hematopoietic Stem Cells Count and Remember Self-Renewal Divisions
.
Cell.
167
:
1296
1309.e10
.
Biasco
,
L.
,
D.
Pellin
,
S.
Scala
,
F.
Dionisio
,
L.
Basso-Ricci
,
L.
Leonardelli
,
S.
Scaramuzza
,
C.
Baricordi
,
F.
Ferrua
,
M.P.
Cicalese
, et al
2016
.
In Vivo Tracking of Human Hematopoiesis Reveals Patterns of Clonal Dynamics during Early and Steady-State Reconstitution Phases
.
Cell Stem Cell.
19
:
107
119
.
Blasius
,
A.L.
,
M.
Cella
,
J.
Maldonado
,
T.
Takai
, and
M.
Colonna
.
2006
.
Siglec-H is an IPC-specific receptor that modulates type I IFN secretion through DAP12
.
Blood.
107
:
2474
2476
.
Bracken
,
B.
2018
. Barcoded Plate-Based Single Cell RNA-Seq. Available at: https://www.protocols.io/view/barcoded-plate-based-single-cell-rna-seq-nkgdctw (accessed March 23, 2018).
Busch
,
K.
,
K.
Klapproth
,
M.
Barile
,
M.
Flossdorf
,
T.
Holland-Letz
,
S.M.
Schlenner
,
M.
Reth
,
T.
Höfer
, and
H.R.
Rodewald
.
2015
.
Fundamental properties of unperturbed haematopoiesis from stem cells in vivo
.
Nature.
518
:
542
546
.
Cabezas-Wallscheid
,
N.
,
D.
Klimmeck
,
J.
Hansson
,
D.B.
Lipka
,
A.
Reyes
,
Q.
Wang
,
D.
Weichenhan
,
A.
Lier
,
L.
von Paleske
,
S.
Renders
, et al
2014
.
Identification of regulatory networks in HSCs and their immediate progeny via integrated proteome, transcriptome, and DNA methylome analysis
.
Cell Stem Cell.
15
:
507
522
.
Carrelha
,
J.
,
Y.
Meng
,
L.M.
Kettyle
,
T.C.
Luis
,
R.
Norfo
,
V.
Alcolea
,
H.
Boukarabila
,
F.
Grasso
,
A.
Gambardella
,
A.
Grover
, et al
2018
.
Hierarchically related lineage-restricted fates of multipotent haematopoietic stem cells
.
Nature.
554
:
106
111
.
Chapple
,
R.H.
,
Y.J.
Tseng
,
T.
Hu
,
A.
Kitano
,
M.
Takeichi
,
K.A.
Hoegenauer
, and
D.
Nakada
.
2018
.
Lineage tracing of murine adult hematopoietic stem cells reveals active contribution to steady-state hematopoiesis
.
Blood Adv.
2
:
1220
1228
.
Copelan
,
E.A.
2006
.
Hematopoietic stem-cell transplantation
.
N. Engl. J. Med.
354
:
1813
1826
.
Drissen
,
R.
,
N.
Buza-Vidas
,
P.
Woll
,
S.
Thongjuea
,
A.
Gambardella
,
A.
Giustacchini
,
E.
Mancini
,
A.
Zriwil
,
M.
Lutteropp
,
A.
Grover
, et al
2016
.
Distinct myeloid progenitor-differentiation pathways identified through single-cell RNA sequencing
.
Nat. Immunol.
17
:
666
676
.
Dykstra
,
B.
,
D.
Kent
,
M.
Bowie
,
L.
McCaffrey
,
M.
Hamilton
,
K.
Lyons
,
S.J.
Lee
,
R.
Brinkman
, and
C.
Eaves
.
2007
.
Long-term propagation of distinct hematopoietic differentiation programs in vivo
.
Cell Stem Cell.
1
:
218
229
.
Eaves
,
C.J.
2015
.
Hematopoietic stem cells: concepts, definitions, and the new reality
.
Blood.
125
:
2605
2613
.
Elias
,
H.K.
,
C.
Schinke
,
S.
Bhattacharyya
,
B.
Will
,
A.
Verma
, and
U.
Steidl
.
2014
.
Stem cell origin of myelodysplastic syndromes
.
Oncogene.
33
:
5139
5150
.
Feil
,
R.
,
J.
Wagner
,
D.
Metzger
, and
P.
Chambon
.
1997
.
Regulation of Cre recombinase activity by mutated estrogen receptor ligand-binding domains
.
Biochem. Biophys. Res. Commun.
237
:
752
757
.
Ghosh
,
H.S.
,
B.
Cisse
,
A.
Bunin
,
K.L.
Lewis
, and
B.
Reizis
.
2010
.
Continuous expression of the transcription factor e2-2 maintains the cell fate of mature plasmacytoid dendritic cells
.
Immunity.
33
:
905
916
.
Grajkowska
,
L.T.
,
M.
Ceribelli
,
C.M.
Lau
,
M.E.
Warren
,
I.
Tiniakou
,
S.
Nakandakari Higa
,
A.
Bunin
,
H.
Haecker
,
L.A.
Mirny
,
L.M.
Staudt
, and
B.
Reizis
.
2017
.
Isoform-Specific Expression and Feedback Regulation of E Protein TCF4 Control Dendritic Cell Lineage Specification
.
Immunity.
46
:
65
77
.
Guo
,
G.
,
S.
Luc
,
E.
Marco
,
T.W.
Lin
,
C.
Peng
,
M.A.
Kerenyi
,
S.
Beyaz
,
W.
Kim
,
J.
Xu
,
P.P.
Das
, et al
2013
.
Mapping cellular hierarchy by single-cell analysis of the cell surface repertoire
.
Cell Stem Cell.
13
:
492
505
.
Haas
,
S.
,
J.
Hansson
,
D.
Klimmeck
,
D.
Loeffler
,
L.
Velten
,
H.
Uckelmann
,
S.
Wurzer
,
A.M.
Prendergast
,
A.
Schnell
,
K.
Hexel
, et al
2015
.
Inflammation-Induced Emergency Megakaryopoiesis Driven by Hematopoietic Stem Cell-like Megakaryocyte Progenitors
.
Cell Stem Cell.
17
:
422
434
.
Inlay
,
M.A.
,
D.
Bhattacharya
,
D.
Sahoo
,
T.
Serwold
,
J.
Seita
,
H.
Karsunky
,
S.K.
Plevritis
,
D.L.
Dill
, and
I.L.
Weissman
.
2009
.
Ly6d marks the earliest stage of B-cell specification and identifies the branchpoint between B-cell and T-cell development
.
Genes Dev.
23
:
2376
2381
.
Ippolito
,
G.C.
,
J.D.
Dekker
,
Y.H.
Wang
,
B.K.
Lee
,
A.L.
Shaffer
III
,
J.
Lin
,
J.K.
Wall
,
B.S.
Lee
,
L.M.
Staudt
,
Y.J.
Liu
, et al
2014
.
Dendritic cell fate is determined by BCL11A
.
Proc. Natl. Acad. Sci. USA.
111
:
E998
E1006
.
Kaitani
,
A.
,
K.
Izawa
,
A.
Maehara
,
M.
Isobe
,
A.
Takamori
,
T.
Matsukawa
,
M.
Takahashi
,
Y.
Yamanishi
,
T.
Oki
,
H.
Yamada
, et al
2018
.
Leukocyte mono-immunoglobulin-like receptor 8 (LMIR8)/CLM-6 is an FcRγ-coupled receptor selectively expressed in mouse tissue plasmacytoid dendritic cells
.
Sci. Rep.
8
:
8259
.
Kent
,
D.G.
,
M.R.
Copley
,
C.
Benz
,
S.
Wöhrer
,
B.J.
Dykstra
,
E.
Ma
,
J.
Cheyne
,
Y.
Zhao
,
M.B.
Bowie
,
Y.
Zhao
, et al
2009
.
Prospective isolation and molecular characterization of hematopoietic stem cells with durable self-renewal potential
.
Blood.
113
:
6342
6350
.
Kiel
,
M.J.
,
O.H.
Yilmaz
,
T.
Iwashita
,
O.H.
Yilmaz
,
C.
Terhorst
, and
S.J.
Morrison
.
2005
.
SLAM family receptors distinguish hematopoietic stem and progenitor cells and reveal endothelial niches for stem cells
.
Cell.
121
:
1109
1121
.
Kunisaki
,
Y.
,
I.
Bruns
,
C.
Scheiermann
,
J.
Ahmed
,
S.
Pinho
,
D.
Zhang
,
T.
Mizoguchi
,
Q.
Wei
,
D.
Lucas
,
K.
Ito
, et al
2013
.
Arteriolar niches maintain haematopoietic stem cell quiescence
.
Nature.
502
:
637
643
.
Lee
,
J.
,
Y.J.
Zhou
,
W.
Ma
,
W.
Zhang
,
A.
Aljoufi
,
T.
Luh
,
K.
Lucero
,
D.
Liang
,
M.
Thomsen
,
G.
Bhagat
, et al
2017
.
Lineage specification of human dendritic cells is marked by IRF8 expression in hematopoietic stem cells and multipotent progenitors
.
Nat. Immunol.
18
:
877
888
.
Macosko
,
E.Z.
,
A.
Basu
,
R.
Satija
,
J.
Nemesh
,
K.
Shekhar
,
M.
Goldman
,
I.
Tirosh
,
A.R.
Bialas
,
N.
Kamitaki
,
E.M.
Martersteck
, et al
2015
.
Highly Parallel Genome-wide Expression Profiling of Individual Cells Using Nanoliter Droplets
.
Cell.
161
:
1202
1214
.
Madisen
,
L.
,
T.A.
Zwingman
,
S.M.
Sunkin
,
S.W.
Oh
,
H.A.
Zariwala
,
H.
Gu
,
L.L.
Ng
,
R.D.
Palmiter
,
M.J.
Hawrylycz
,
A.R.
Jones
, et al
2010
.
A robust and high-throughput Cre reporting and characterization system for the whole mouse brain
.
Nat. Neurosci.
13
:
133
140
.
Merad
,
M.
,
P.
Sathe
,
J.
Helft
,
J.
Miller
, and
A.
Mortha
.
2013
.
The dendritic cell lineage: ontogeny and function of dendritic cells and their subsets in the steady state and the inflamed setting
.
Annu. Rev. Immunol.
31
:
563
604
.
Murphy
,
T.L.
,
G.E.
Grajales-Reyes
,
X.
Wu
,
R.
Tussiwand
,
C.G.
Briseño
,
A.
Iwata
,
N.M.
Kretzer
,
V.
Durai
, and
K.M.
Murphy
.
2016
.
Transcriptional Control of Dendritic Cell Development
.
Annu. Rev. Immunol.
34
:
93
119
.
Muzumdar
,
M.D.
,
B.
Tasic
,
K.
Miyamichi
,
L.
Li
, and
L.
Luo
.
2007
.
A global double-fluorescent Cre reporter mouse
.
Genesis.
45
:
593
605
.
Naik
,
S.H.
,
L.
Perié
,
E.
Swart
,
C.
Gerlach
,
N.
van Rooij
,
R.J.
de Boer
, and
T.N.
Schumacher
.
2013
.
Diverse and heritable lineage imprinting of early haematopoietic progenitors
.
Nature.
496
:
229
232
.
Nestorowa
,
S.
,
F.K.
Hamey
,
B.
Pijuan Sala
,
E.
Diamanti
,
M.
Shepherd
,
E.
Laurenti
,
N.K.
Wilson
,
D.G.
Kent
, and
B.
Göttgens
.
2016
.
A single-cell resolution map of mouse hematopoietic stem and progenitor cell differentiation
.
Blood.
128
:
e20
e31
.
Nishikii
,
H.
,
Y.
Kanazawa
,
T.
Umemoto
,
Y.
Goltsev
,
Y.
Matsuzaki
,
K.
Matsushita
,
M.
Yamato
,
G.P.
Nolan
,
R.
Negrin
, and
S.
Chiba
.
2015
.
Unipotent Megakaryopoietic Pathway Bridging Hematopoietic Stem Cells and Mature Megakaryocytes
.
Stem Cells.
33
:
2196
2207
.
Notta
,
F.
,
S.
Zandi
,
N.
Takayama
,
S.
Dobson
,
O.I.
Gan
,
G.
Wilson
,
K.B.
Kaufmann
,
J.
McLeod
,
E.
Laurenti
,
C.F.
Dunant
, et al
2016
.
Distinct routes of lineage development reshape the human blood hierarchy across ontogeny
.
Science.
351
:
aab2116
.
Olsson
,
A.
,
M.
Venkatasubramanian
,
V.K.
Chaudhri
,
B.J.
Aronow
,
N.
Salomonis
,
H.
Singh
, and
H.L.
Grimes
.
2016
.
Single-cell analysis of mixed-lineage states leading to a binary cell fate choice
.
Nature.
537
:
698
702
.
Paul
,
F.
,
Y.
Arkin
,
A.
Giladi
,
D.A.
Jaitin
,
E.
Kenigsberg
,
H.
Keren-Shaul
,
D.
Winter
,
D.
Lara-Astiaso
,
M.
Gury
,
A.
Weiner
, et al
2015
.
Transcriptional Heterogeneity and Lineage Commitment in Myeloid Progenitors
.
Cell.
163
:
1663
1677
.
Pei
,
W.
,
T.B.
Feyerabend
,
J.
Rössler
,
X.
Wang
,
D.
Postrach
,
K.
Busch
,
I.
Rode
,
K.
Klapproth
,
N.
Dietlein
,
C.
Quedenau
, et al
2017
.
Polylox barcoding reveals haematopoietic stem cell fates realized in vivo
.
Nature.
548
:
456
460
.
Perié
,
L.
,
K.R.
Duffy
,
L.
Kok
,
R.J.
de Boer
, and
T.N.
Schumacher
.
2015
.
The Branching Point in Erythro-Myeloid Differentiation
.
Cell.
163
:
1655
1662
.
Picelli
,
S.
,
O.R.
Faridani
,
A.K.
Björklund
,
G.
Winberg
,
S.
Sagasser
, and
R.
Sandberg
.
2014
.
Full-length RNA-seq from single cells using Smart-seq2
.
Nat. Protoc.
9
:
171
181
.
Pietras
,
E.M.
,
D.
Reynaud
,
Y.A.
Kang
,
D.
Carlin
,
F.J.
Calero-Nieto
,
A.D.
Leavitt
,
J.M.
Stuart
,
B.
Göttgens
, and
E.
Passegué
.
2015
.
Functionally Distinct Subsets of Lineage-Biased Multipotent Progenitors Control Blood Production in Normal and Regenerative Conditions
.
Cell Stem Cell.
17
:
35
46
.
Roch
,
A.
,
V.
Trachsel
, and
M.P.
Lutolf
.
2015
.
Brief Report: Single-Cell Analysis Reveals Cell Division-Independent Emergence of Megakaryocytes From Phenotypic Hematopoietic Stem Cells
.
Stem Cells.
33
:
3152
3157
.
Rodriguez-Fraticelli
,
A.E.
,
S.L.
Wolock
,
C.S.
Weinreb
,
R.
Panero
,
S.H.
Patel
,
M.
Jankovic
,
J.
Sun
,
R.A.
Calogero
,
A.M.
Klein
, and
F.D.
Camargo
.
2018
.
Clonal analysis of lineage fate in native haematopoiesis
.
Nature.
553
:
212
216
.
Sánchez-Aguilera
,
A.
,
L.
Arranz
,
D.
Martín-Pérez
,
A.
García-García
,
V.
Stavropoulou
,
L.
Kubovcakova
,
J.
Isern
,
S.
Martín-Salamanca
,
X.
Langa
,
R.C.
Skoda
, et al
2014
.
Estrogen signaling selectively induces apoptosis of hematopoietic progenitors and myeloid neoplasms without harming steady-state hematopoiesis
.
Cell Stem Cell.
15
:
791
804
.
Sanjuan-Pla
,
A.
,
I.C.
Macaulay
,
C.T.
Jensen
,
P.S.
Woll
,
T.C.
Luis
,
A.
Mead
,
S.
Moore
,
C.
Carella
,
S.
Matsuoka
,
T.
Bouriez Jones
, et al
2013
.
Platelet-biased stem cells reside at the apex of the haematopoietic stem-cell hierarchy
.
Nature.
502
:
232
236
.
Satija
,
R.
,
J.A.
Farrell
,
D.
Gennert
,
A.F.
Schier
, and
A.
Regev
.
2015
.
Spatial reconstruction of single-cell gene expression data
.
Nat. Biotechnol.
33
:
495
502
.
Sawai
,
C.M.
,
S.
Babovic
,
S.
Upadhaya
,
D.J.H.F.
Knapp
,
Y.
Lavin
,
C.M.
Lau
,
A.
Goloborodko
,
J.
Feng
,
J.
Fujisaki
,
L.
Ding
, et al
2016
.
Hematopoietic Stem Cells Are the Major Source of Multilineage Hematopoiesis in Adult Animals
.
Immunity.
45
:
597
609
.
Shizuru
,
J.A.
,
R.S.
Negrin
, and
I.L.
Weissman
.
2005
.
Hematopoietic stem and progenitor cells: clinical and preclinical regeneration of the hematolymphoid system
.
Annu. Rev. Med.
56
:
509
538
.
Sieburg
,
H.B.
,
R.H.
Cho
,
B.
Dykstra
,
N.
Uchida
,
C.J.
Eaves
, and
C.E.
Muller-Sieburg
.
2006
.
The hematopoietic stem compartment consists of a limited number of discrete stem cell subsets
.
Blood.
107
:
2311
2316
.
Sun
,
J.
,
A.
Ramos
,
B.
Chapman
,
J.B.
Johnnidis
,
L.
Le
,
Y.J.
Ho
,
A.
Klein
,
O.
Hofmann
, and
F.D.
Camargo
.
2014
.
Clonal dynamics of native haematopoiesis
.
Nature.
514
:
322
327
.
Tusi
,
B.K.
,
S.L.
Wolock
,
C.
Weinreb
,
Y.
Hwang
,
D.
Hidalgo
,
R.
Zilionis
,
A.
Waisman
,
J.R.
Huh
,
A.M.
Klein
, and
M.
Socolovsky
.
2018
.
Population snapshots predict early haematopoietic and erythroid hierarchies
.
Nature.
555
:
54
60
.
Umemoto
,
T.
,
M.
Yamato
,
J.
Ishihara
,
Y.
Shiratsuchi
,
M.
Utsumi
,
Y.
Morita
,
H.
Tsukui
,
M.
Terasawa
,
T.
Shibata
,
K.
Nishida
, et al
2012
.
Integrin-αvβ3 regulates thrombopoietin-mediated maintenance of hematopoietic stem cells
.
Blood.
119
:
83
94
.
Yamamoto
,
R.
,
Y.
Morita
,
J.
Ooehara
,
S.
Hamanaka
,
M.
Onodera
,
K.L.
Rudolph
,
H.
Ema
, and
H.
Nakauchi
.
2013
.
Clonal analysis unveils self-renewing lineage-restricted progenitors generated directly from hematopoietic stem cells
.
Cell.
154
:
1112
1126
.
Yu
,
V.W.C.
,
R.Z.
Yusuf
,
T.
Oki
,
J.
Wu
,
B.
Saez
,
X.
Wang
,
C.
Cook
,
N.
Baryawno
,
M.J.
Ziller
,
E.
Lee
, et al
2016
.
Epigenetic Memory Underlies Cell-Autonomous Heterogeneous Behavior of Hematopoietic Stem Cells
.
Cell.
167
:
1310
1322.e17
.

Author notes

*

S. Upadhaya and C.M. Sawai contributed equally to this paper.

C.M. Sawai’s present address is Institut National de la Santé et de la Recherche Médicale Unit 1218 ACTION Laboratory, University of Bordeaux, Bergonié Cancer Institute, Bordeaux, France.

This article is distributed under the terms of an Attribution–Noncommercial–Share Alike–No Mirror Sites license for the first six months after the publication date (see http://www.rupress.org/terms/). After six months it is available under a Creative Commons License (Attribution–Noncommercial–Share Alike 4.0 International license, as described at https://creativecommons.org/licenses/by-nc-sa/4.0/).

Supplementary data