Geographic clonal tracking in macaques provides insights into HSPC migration and differentiation

Adult hematopoietic stem and progenitor cells (HSPCs) reside in bone marrow (BM) microenvironmental niches that were first proposed theoretically, then identified functionally and quantitatively as a limited resource controlling HSPC numbers and behavior, and recently characterized at increasingly complex cellular and molecular levels (Schofield, 1978; Morrison and Scadden, 2014). During fetal development, hematopoiesis occurs sequentially in several anatomical sites, including the yolk sac, aorta-gonado-mesonephros region, and liver, before finally localizing in the BM before birth (Costa et al., 2012). HSPCs have been shown to move from site to site via the vasculature during this developmental process. This migratory ability has been taken advantage of clinically, via harvest of HSPCs from the BM or after pharmacologic mobilization of HSPCs into the peripheral blood (PB) and transplantation via simple infusion of HSPCs into the bloodstream, followed by engraftment of transplanted HSPCs in BM niches (Körbling and Freireich, 2011).

At steady state, a very small number of HSPCs can be found circulating in the blood in model animals and humans (Goodman and Hodgson, 1962; Papayannopoulou and Scadden, 2008). The number of HSPCs in blood increases with stress, during recovery from myelosuppression, and in various pathological HSPC conditions such as myeloproliferative disorders (Richman et al., 1976; Hoggatt et al., 2013). The physiological role of these cells is unclear. Various theories regarding control of self-renewal versus commitment, for both hematopoietic and nonhematopoietic stem cell types, suggest that with cell division, one daughter stem cell remains in a niche and the other differentiates; dies upon movement out of the niche, because of loss of niche signals required to retain stemness; or migrates to an open niche (Yamashita et al., 2007; Morrison and Scadden, 2014). Parabiosis as well as nonablative transplantation experiments in mice demonstrate extremely slow mixing of HSPCs in the BM, suggesting that, at least in mice, available niches are full at steady state and that exit from the BM is primarily a death pathway (Abkowitz et al., 2003; Chen et al., 2006; Czechowicz et al., 2007). Mobilization of endogenous HSPCs out of the BM with cytokines or antibodies interrupting the interaction of HSPC receptors with niche factors can enhance engraftment of exogenous transplanted HSPCs (Chen et al., 2006).

After myeloablative transplantation, HSPCs home to BM niches and proliferate very rapidly, regenerating the long-term repopulating stem cell pool within months in mice (Pawliuk et al., 1996). During the recovery process, very rapid HSPC proliferation after initial niche engraftment combined with the inflammatory stress related to conditioning might be expected to result in release of daughter HSPCs into the circulation and reseeding into new distant niches, predicting rapid homogenization of the progeny of individual HSPCs throughout an organism’s entire BM space. However, a recent mouse study demonstrated differences in chimerism levels between bones after competitive transplantation (Rundberg Nilsson et al., 2015). Further insights into the geographic process of clonal HSPC spread requires methodology able to identify and localize the output of individual clones in vivo.

We transplanted mouse HSPCs transduced with multihued lentiviral gene ontology (LeGO) fluorescent lentiviral vectors, allowing discrimination of at least 50 different clones simultaneously via confocal imaging of the BM, and were surprised to find geographically restricted macroscopic output from individual HSPCs as late as 4 mo posttransplantation, long after counts normalized, and at a time the BM had reached 100% cellularity (Malide et al., 2012). Another group reported asymmetric distribution of mouse HSPC clones across the skeleton and reequilibration upon a single G-CSF challenge (Verovskaya et al., 2014). Luciferase imaging of immunodeficient mice after transplantation of human HSPCs revealed initial local foci of hematopoiesis, followed by eventual spread to distant locations; however, the resolution of this approach is low, and the output of individual HSPCs cannot be tracked in this manner (Wang et al., 2003).

Given the differences between human and mouse hematopoiesis, it is unclear whether insights regarding niche interactions and trafficking can be extrapolated from mouse models (Doulatov et al., 2012). Mouse BM remains almost 100% cellular throughout life, in marked contrast to the patchy cellularity of adult human BM. Even after transplantation or recovery from a cytotoxic event, human BM biopsies are characterized by macroscopic foci of active hematopoiesis, separated by large areas of fatty BM. It is not clear whether areas without active hematopoiesis exist because of loss of niche activity or a HSPC life cycle in humans that does not rely on constant vascular release and rehoming.

We have previously used rhesus macaque HSPC autologous transplantation combined with lentiviral cellular barcoding to study hematopoiesis at a clonal level in a model with great relevance to human hematopoiesis (Donahue and Dunbar, 2001; Wu et al., 2014; Koelle et al., 2017). Macaques and humans have prolonged lifespans, homologous BM and secondary lymphoid organ morphology and cellularity, and similar HPSC cycling and dynamics (Shepherd et al., 2007). Here, we applied this approach to visualize hematopoiesis geographically at a clonal level and directly analyze the location of lineage-specific differentiation pathways. Even in the context of rapid rhesus macaque HSPC expansion after fully myeloablative transplantation, hematopoiesis remained highly geographically localized within the BM at a clonal level, with very slow mixing across different BM sites for a period of a year or more. Lineage analyses revealed a novel BM-restricted CD3⁺ T cell developmental pathway and further insights into natural killer (NK) cell production and maturation.

As shown in Fig. 1 A, autologous CD34⁺ HSPCs were tagged with a high-diversity lentiviral barcoded library, labeling each HSPC and all its progeny with a unique cellular DNA barcode after vector transduction and integration, allowing tracking of in vivo output from individual labeled HSPCs as described previously (Donahue et al., 2005; Wu et al., 2014; Koelle et al., 2017). We applied this approach to track the clonal geographic distribution of hematopoiesis after myeloablative transplantation. Transplantation and lentiviral transduction parameters for the seven animals included in this study are summarized in Fig. 1 A and Table 1. We obtained concurrent PB and large (10–15-ml) BM aspirates from both left and right iliac crests or ischial pads at various time points from 3.5 to 46 mo after transplantation. The left versus right iliac crests or ischial pads are separated by joints or ligamentous connections and not connected via BM cavities; thus HSPC clones present on both sides must have spread hematogenously. Granulocytes (Gr), monocytes (Mono), and T, B, and NK cells were sorted to high purity from all samples, and CD34⁺ HSPCs were purified from each BM sample (Fig. S1, A and B). The percentage of GFP⁺CD34⁺ cells, and thus lentivirally barcoded CD34⁺ HSPCs in left and right BM samples over time, are plotted in Fig. 1 B and in various PB populations in Fig. S1 C. GFP levels in BM CD34⁺ HSPCs varied between animals; however, in each individual animal, the levels between the left and right sides correlated closely (r = 0.991). Lymph nodes (LNs) were also collected in three animals at selected time points to allow comparison of the clonal composition of LNs to PB. Quantitative retrieval of barcodes from each cell population was performed via low-cycle PCR, Illumina HiSeq sequencing, application of custom algorithms to identify valid barcodes, and application of thresholds to account for sequencing errors and sampling constraints, as described and validated in detail in our previous publications (Wu et al., 2014; Koelle et al., 2017).

To be able to separate BM tissue-resident cells from circulating cells in these highly vascular tissues, we adopted an intravascular staining approach used previously to discriminate tissue-resident versus circulating lymphocytes (Anderson et al., 2012). An anti-CD45 antibody conjugated with the fluorescent moiety Alexa Fluor 647 was administrated intravenously to one barcoded macaque 4 mo posttransplantation, and BM and PB cells were collected concurrently 4 h postinfusion. Hematopoietic cells of all lineages other than maturing erythroid cells are CD45 positive, and those cells within the circulation bind this antibody rapidly. Cells stained with CD45–Alexa Fluor 647 represent PB-resident circulating hematopoietic CD45⁺ cells, defined here as intravascular (IVas)⁺ cells, whereas CD45-Alexa Fluor 647–negative (IVas⁻) cells in BM or other tissue sample represent tissue-resident hematopoietic CD45⁺ cells, either locally produced or homed back to the tissue from the blood before antibody infusion. 99.6% of PB CD45⁺ mononuclear cells were IVas⁺ at this time point, demonstrating that antibody dosing was sufficient to stain virtually all intravascular cells. In BM aspirates, 20–24% of mononuclear CD45⁺ cells were IVas⁺. Virtually all BM CD34⁺ HSPCs were IVas⁻ (Fig. 1 C), given the very low level of circulating CD34⁺ cells at steady state, supporting the use of anti-CD45 intravenous infusion followed by tissue sampling as a technology to discriminate tissue-resident from circulating cells in hematopoietic tissues.

We have previously reported that, by 3–6 mo posttransplantation, clonal hematopoietic patterns in the blood of rhesus macaques stabilize, with Gr, Mono, B, and T cells primarily derived from thousands of multilineage long-term repopulating HSPC clones (Wu et al., 2014; Koelle et al., 2017). Blood counts in these animals recover by 1 mo posttransplantation, but the earliest time points that sufficient spicules and CD34⁺ HSPCs can be reproducibly aspirated from the BM for barcode analysis is 3–7 mo. In large-volume BM samples collected at these initial time points, clonal contributions within CD34⁺ HSPCs obtained from left versus right skeletal BM sites were almost completely distinct in five of seven animals, visualized at a clonal level for all clones via clonal dot plots and clonal histograms plotting the degree of left–right bias (Fig. 2 A), and via Pearson correlations between the populations of clones contributing to CD34⁺ HSPCs from each geographically distinct BM sample (Fig. 2 B). In the remaining two animals, one at 5 mo (ZJ38) and another at 7 mo (ZG66), some clones were shared between left and right BM samples; however, many clones were still completely segregated between the BM sites. Both large and small clones were highly segregated, as can be seen from the position of clones along the x and y axes in the dot plots and the variable clone sizes in the histogram stacks in Fig. 2 A. The median Pearson correlation between left and right BM sites for all seven animals sampled between 3 and 7 mo was only 0.002 (range, 0–0.541).

We compared clonal patterns in CD34⁺ HSPCs at each BM site to Gr purified from the same BM sample and to circulating Gr. At time points ranging from 3 to 7 mo, the clonal composition of BM Gr correlated very closely with CD34⁺ HSPCs from the same BM location (Fig. 2 B), supporting local production of these myeloid cells from CD34⁺ HSPCs. Heat maps analyzing contributions from the top 10 clones in each sample, mapped over all samples and grouped by clustering of clones with similar patterns of contribution, provide additional information and allow visualization of the contributions from individual clones or groups of clones (Fig. 2 C). The top 10 clones from each sample were chosen and mapped for clarity; however, heat maps of contributions from the top 100 clones from each sample (Fig. S2) show the same patterns. A focus on larger clones further obviates sampling constraints and provides information on clones with the greatest impact on hematopoiesis. The major clones in each BM site were found in both CD34⁺ HSPCs and Gr samples from the same site, but not to CD34⁺ HSPCs or Gr from the other BM location sampled concurrently in ZJ31, ZK22, ZL40, and ZH33. A few major clones had begun to mix in ZJ38 and ZH19, and there was already more extensive mixing of CD34⁺ HSPC clones in ZG66 at 7 mo. Even the highest contributing Gr clones in the PB were not found in the BM CD34⁺ HSPCs or Gr from either side in some animals (ZK22, ZL40, ZH19), suggesting that these particular PB clones were not being produced at either of these BM sites, and that the largest clones contributing to ongoing hematopoiesis were not necessarily mixed early or evenly across the entire BM. Major PB clones in some animals (ZJ31, ZH33) were found primarily in Gr but not CD34⁺ HSPCs on one side or the other, suggesting greater PB contamination of some BM aspirates, with some heterogeneity between each site and animal (Fig. 2 C). In ZJ38, a few major PB clones were shared with Gr and CD34⁺ HSPCs primarily on the right side, indicating that the site of production of these clones happened to have been sampled. In ZJ38 and ZG66, some of the largest PB Gr clones were found in Gr and in CD34⁺ HSPCs on both sides of BM, as expected, given some clonal mixing between CD34⁺ HSPCs on the two sides at these time points in these two animals.

In summary, even long after massive proliferation of transplanted engrafting HSPCs and normalization of blood counts, HSPCs spread locally, rather than immediately mixing via release into the blood and reseeding of distant BM niches. These findings suggest that even regenerating hematopoiesis spreads to adjacent niches in the BM for months posttransplantation.

We observed gradual increases in correlations between clonal contributions to CD34⁺ HSPCs at different sites when BM was serially sampled at later times posttransplantation, visualized at both population and clonal levels (Fig. 3, A–C). The fraction of clones highly biased toward contributions to left versus right BM demonstrated a gradual loss of geographic clonal segregation over time (Fig. 3 D). The pace and degree of clonal mixing varied between animals, but by 1 y, the majority of hematopoiesis was generated by clones that were found to be present and contributing at more than one location (Fig. 3 D). After mixing had occurred, clonal patterns were quite stable in BM CD34⁺ cells, mirroring the clonal stability we recently reported in circulating PB cells (Koelle et al., 2017). At longest follow-up, once CD34⁺ HSPC clone distribution homogenized between different BM sites, the top contributing clones in the PB Gr were now also generally the top contributing clones in both sampled BM sites (Fig. 3, E and F). Over time, the Shannon diversity index of CD34⁺ HSPCs in the BM increased and approached that of the blood, reflecting more complete clonal mixing and homogeneous distribution (Fig. S3).

We administered the HSPC-mobilizing cytokines G-CSF and stem cell factor (SCF) to ZH33 6.5 mo after transplantation to collect larger numbers of HSPCs for ongoing studies. Samples collected subsequent to cytokine mobilization showed an increase in clones shared between CD34⁺ HSPCs at left and right BM sites, but the pace of mixing seemed no faster than the changes that occurred in other animals between 6 mo and 1 y (Fig. 3), suggesting that most of the mobilized HSPCs did not home back to a new niche and commence contributing to hematopoiesis; however, the results in a single animal are not definitive.

In a subset of animals, we performed flow cytometric sorting of cells from each BM site and concurrent PB to gain further insights into sites of production for each hematopoietic lineage. We used lineage-defining phenotypic markers (Fig. S1 and Table S1) for purification of T, B, Mono, Gr, and CD16⁺ and/or CD56⁺ NK cells. In ZJ31 at 3.5 mo, Fig. 4 A demonstrates groups of CD34⁺ HSPC clones present in only the left or only the right BM and contributing to local production of Gr, Mono, and B cells, as documented by the presence of geographically restricted clones in both CD34⁺ HSPCs and highly purified cells of these lineages. We defined highly biased clones as those showing at least a 10-fold higher fractional contribution to CD34⁺ HSPCs on one side compared with the other side, and in ZJ31 at 3.5 mo, ∼90% of CD34⁺ HSPCs were highly biased on each side. Highly biased clones accounted for the majority of Gr, Mono, and B cell production when analyzing the entire population of clones (Fig. 4 B). A similar pattern was seen for samples obtained from ZH19 at 6 mo, although more mixing between BM sites had already occurred in this animal by this time point (Fig. 4, C and D).

Some Mono, Gr, and B cells in the left and right BM were produced from unbiased clones in both animals. These shared clones were also the dominant clones found in the PB lineages at the same time points, and some clones were also present in a small fraction of left and right BM CD34⁺ HSPCs. These unbiased clones could represent contaminating PB cells present in the BM aspirates, mature cells that migrated from a different BM production site and then homed back to the both the left and right BM sites, or local output of CD34⁺ HSPC clones that had already migrated to both BM sites in the case of cells matching CD34⁺ HSPCs from all samples.

NK cells in the rhesus macaque are defined as CD3⁻CD20⁻CD14⁻ mononuclear cells expressing either CD56 or CD16 (Webster and Johnson, 2005; Wu et al., 2014). Macaque CD56⁺CD16⁻ NK cells (here referred to as CD56⁺) are equivalent in phenotypic characteristics and functions such as cytokine secretion to human CD56^bright NK cells, and are hypothesized to be less mature NK cells. Macaque CD56⁻CD16⁺ NK cells (designated CD16⁺) are equivalent to human CD56^dimCD16⁺ cytotoxic cells and are hypothesized to be fully mature NK cells. At 3.5 mo in ZJ31, the vast majority of BM NK cells were CD56⁺, in contrast to the predominant CD16⁺ subset in the PB (Fig. S1), similar to the distribution of the two main human NK subsets between BM and PB and a prior study in macaques (Webster and Johnson, 2005). These NK subsets were sorted from the PB and from left and right BM of ZJ31 at 3.5 mo and ZH19 at 6 mo. The clonal pattern in BM CD56⁺ NK cells showed marked geographic segregation (Fig. 4), with locally produced CD56⁺ NK clones found only in the left or right BM. Of note, the subset of local CD34⁺ HSPCs producing BM CD56⁺ NK cells seemed to differ somewhat in relative contributions to CD56⁺ NK cells versus local HSPC-produced Gr, Mono, and B cells found at the same site, suggesting some clonal bias of progenitors producing this NK subset in the BM versus B and myeloid cells.

As we have previously reported (Wu et al., 2014), the majority of circulating NK cells, specifically those within the CD16⁺ subset, were clonally distinct from circulating Gr, Mono, B, and T cells, which we demonstrated were produced from long-term stable multipotent clones beginning 3–4 mo posttransplantation. The heat maps show that the BM CD16⁺ NK clones contribute at similar levels on left and right sides, do not match local CD34⁺ HSPCs, and are thus not produced from local CD34⁺ HSPCs at 3.5 mo (ZJ31) and 6 mo (ZH19). Instead, their clonal distribution closely matches that of CD16⁺ PB NK cells (Fig. 4). Similar results were found in ZJ31 at 8.5 mo (Fig. S4). These results directly demonstrate that CD56⁺ NK cells are produced in the BM, in contrast to CD16⁺ NK, which are produced outside the BM. Of note, CD56⁺ and CD16⁺ NK cells in the PB, reflecting total-body production of NK cells, are clonally distinct from each other, with CD56⁺ clonality more closely matching myeloid and other lymphoid lineages, confirming our prior study, and suggesting that CD16⁺ NK cells are not produced in an ongoing manner from the same subset of CD34⁺ HSPCs that generate other hematopoietic lineages (Wu et al., 2014).

In monkey ZL40, intravascular staining confirmed these results, revealing that virtually 100% of CD16⁺ NK cells found in both left and right BM aspirates were IVas⁺, indicating intravascular localization in the hours before BM aspiration, and not local BM production or homing back and long-term residence in the BM after production elsewhere (Fig. 5 A). These BM CD16⁺ NK cells were likely present intravascularly within the BM space, and thus represent PB contamination. In contrast, BM CD56⁺ NK cells were primarily IVas⁻, suggesting local production or rehoming of these cells outside the vasculature (Fig. 5 A), with the clonal mapping supporting local production, given shared clonality with local HSPCs.

CD3 surface expression is associated with T cell maturation and gene rearrangement followed by surface expression of the TCR complex, and is thought to occur almost exclusively in the thymus. Barcoded GFP⁺ T cells, derived from transplanted CD34⁺ HSPCs, regenerate more slowly than other lineages after transplantation (Fig. S1 C), and T cell production from engrafted CD34⁺ cells after ablative transplantation lags behind that of other lineages, for reasons thought to be related to thymic damage from irradiation. We expected that all CD3⁺ T cells present in the BM and derived from transplanted CD34⁺ HSPCs would carry barcodes shared between T cells in both the left and right BM samples and match the barcodes in PB T cells, reflecting homing of T cells back to the BM via the PB, after maturation in the thymus. Indeed, in ZJ31 at 3.5 mo, there was a cluster of T clones highly represented in PB and in both left and right BM, as expected (Fig. 4 A). Surprisingly, we also detected a prominent group of CD3⁺ T clones present in only the left or only the right BM, clonally related to the localized CD34⁺ HSPCs and to locally produced B and myeloid cells (Fig. 4, A and B), even with very stringent gating and sorting purities of >98%. Less abundant local clonally restricted T clones were found in ZH19 at 6 mo (Fig. 4, C and D). By 8.5 mo in ZJ31, these locally produced T cells were no longer detectable (Fig. S4), despite some degree of continued geographic segregation of clones contributing to CD34⁺ HSPCs and B and myeloid cells in ZJ31. These results suggest that a population of CD3⁺ T cells is produced and matured locally in the BM early after transplantation, before robust thymic recovery, independent of thymic homing and maturation.

We examined this surprising finding more intensively in ZL40, after intravenous administration of the anti-CD45–Alexa Fluor 647 antibody. 2–3% of CD3⁺ cells in the BM were IVas⁻, suggesting either local production or homing of these CD3⁺ T cells back to the BM (Fig. 5 B). At the same time point, PB CD3⁺ T cells were 99.9% IVas⁺ (Fig. 5 B). The GFP⁺ fractions of resident BM CD3⁺IVas⁻ cells were 46.0% (left BM) and 55.5% (right BM), similar to the BM CD34⁺ HSPC GFP⁺ fraction (53.9% and 51.4% in left and right BM, respectively), and considerably higher than the GFP⁺ fraction in BM CD3⁺IVas⁺ T and circulating PB CD3⁺IVas⁺ T cells (19.9% in left BM, 18.9% in right BM, and 20.3% in PB; Fig. 5 B), also supporting the concept of two subsets of T cells, one locally produced in the BM early after transplantation, and the other requiring thymic production. The CD3⁺ IVas⁻ cells expressed CD8 (22.6%), CD4 (10.3%), and CD56 (83.8%), and 38% were positive for the C-C chemokine receptor type 5 (CCR5) often expressed on tissue-resident T cells (Mattapallil et al., 2005). The vast majority expressed TCR-αβ (89.5%), with smaller populations expressing TCR-γδ (2.5%) and TCR-Vα24 (4.52%), a stereotypical TCR expressed on rhesus macaque NKT cells (Rout et al., 2010, 2012), as shown in Fig. 5 (C and D).

These findings suggest that BM CD3⁺ tissue-resident T cells may mature directly from local BM HSPCs, without transit of T precursors to the thymus to complete maturation. We further analyzed barcodes present in sorted CD3⁺IVas⁻ and CD3⁺IVas⁺ cells. Barcodes present in the CD3⁺IVas⁻ cells sorted from the left BM matched the left BM CD34⁺ barcodes, providing two lines of evidence documenting local BM CD3⁺ T cell production from CD34⁺ HSPCs, rather than homing back from the thymus. This was true for both CD4⁺ and CD8⁺ T cells (Fig. 5 E). In contrast, the CD3⁺IVas⁺ BM T cells had barcodes shared with PB CD3⁺IVas⁺ T cells, again in both CD4 and CD8 subsets (Fig. 5 E).

We compared the clonal composition of concurrent PB and LN T and B cells in ZG66 (3 mo) and ZH33 (9.5 mo) and found very strong positive correlations between clonal contributions found in PB versus whole LNs (Fig. 6 A). The exception was a single massively expanded T cell clone present in the PB of ZG66, accounting for ∼20% of marked cells at this time point, but not showing clonal dominance in the concurrently harvested node, suggesting that this clone was expanding at another location besides the single LN sampled.

We also harvested left and right LNs and PB concurrently from ZG66 at 6.5 mo and sorted T, B, and NK cells and NK subsets (Fig. 6 B). The clonal composition of the T and B cells in the left LN highly correlated with the right LN and PB T or B cells (Fig. 6 C). Interestingly, the right and left LN whole NK populations also strongly correlated with each other, but not with PB NK cells (Fig. 6 D). We further separated the PB NK cells into CD16⁺ and CD56⁺ NK subsets and found the left and right LN NK cells, which were primarily CD56⁺ cells, to be correlated with the PB CD56⁺ NK subset, and not with the circulating CD16⁺ NK cell subset (Fig. 6 E). At least early after myeloablative transplantation, it appears that LNs are repopulated quite homogeneously by circulating T, B, and CD56⁺ NK cells, and that specific expansion of clones in the node in response to antigens has not yet occurred.

In this study, we used clonal tracking and separation of BM-resident versus circulating cells to document prolonged geographic segregation of HSPCs and their progeny after myeloablative transplantation in a rhesus macaque model. Even during the intense proliferative period after HSPC homing, engraftment, and hematopoietic recovery, daughter HSPCs generally spread locally, rather than exiting the BM and traveling via the circulation to a distant niche. This pattern of hematopoietic recovery in macaques is particularly relevant in light of the extremely heterogeneous BM cellularity observed in patients recovering from HSPC transplantation, chemotherapy, or aplastic anemia, with macroscopic areas of 100% cellularity interspersed between areas containing only adipocytes and stromal cells. This geographically restricted pattern of hematopoiesis suggests an actual limitation in the number of functional niches, perhaps resulting from irradiation, or physiology favoring local spread of HSPCs.

Continued localization of proliferating and differentiating HSPCs after initial engraftment of individual cells has been observed using live imaging of the mouse BM (Lo Celso et al., 2009). BM structures termed “hemospheres,” consisting of several niche cell types, have been shown to encase large proliferating clusters of endogenous or transplanted mouse HSPCs, with or without prior ablative irradiation (Wang et al., 2013). There is evidence that several types of HSPC-derived cells are important for maintenance of quiescent HSCs and/or differentiation of cycling progenitors, including monocyte/macrophages, megakaryocytes, neutrophils, and T reg cells (Winkler et al., 2010; Zhao et al., 2014; Sánchez-Aguilera and Méndez-Ferrer, 2017), often via secretion of cytokines. Thus, once an HSPC has engrafted, it can remodel the BM neighborhood to be more conducive to hematopoiesis. Conversely, adipocytes inhibit proliferation and survival of HSPCs, and thus may prevent facile spread of initially engrafted cells or their daughter cells to distant sites (Naveiras et al., 2009; Chitteti et al., 2010). Our results strongly support a model of local modulation of niche function, with a highly potent positive feedback loop during rapid BM regeneration. Daughter cells resulting from initial HSPC divisions and their downstream progeny spread contiguously, before final maturation and release into the bloodstream. HSPCs that do enter the circulation, at least after the initial factors driving rapid proliferation abate, may either die before finding a new niche or remain quiescent once they reach a new niche in a hematopoietically inactive site.

HSPC mobilization accelerates clonal mixing, both in our rhesus macaque model, as suggested in one animal (ZH33) after G-CSF and SCF treatment at 6.5 mo posttransplantation, and in a prior mouse study (Verovskaya et al., 2014). Whether mixing involves settlement of HSPCs in previously unoccupied niches or mixing of HSPCs within already occupied and conditioned niches is unclear, and an interesting question to address in future studies.

In ZJ38, which received the lowest dose of CD34⁺ HSPCs (Table 1), we observed the earliest HSPC mixing in different bone locations, with overall clonal diversity lower than in any of the other barcoded animals. It is possible that in the setting of low HSPC doses, proliferation and continued hematopoietic stress overwhelm the ability of cells to spread contiguously, leading to faster migration and mixing. It is interesting to note that in a xenograft model tracking clonal behavior of barcoded human acute lymphoblastic leukemia cells, lower leukemic cell doses resulted in more marked diffusion of clones throughout the recipient mice compared with higher cell doses (Elder et al., 2017).

HSPCs differentiate into a variety of daughter lineages via a series of developmental steps. Full maturation of neutrophils, Mono, and B cells is thought to occur in the BM (Jacobsen and Osmond, 1990; Kondo, 2010; Melchers, 2015). In this study, we provide direct evidence for local production of Gr (predominantly neutrophils), Mono, and B cells via demonstration of shared clonality between local CD34⁺ HSPCs and these cell types.

Knowledge regarding the development and ontogeny of NK cells is more limited than for most other lineages. Human NK cells are not produced robustly in xenograft models, and mouse and human NK cells differ in many characteristics. NK cell differentiation from HSPCs has been inferred from patterns of acquisition of cell surface markers at different sites in vivo and during in vitro maturation (Freud et al., 2014). NK cells are hypothesized to develop in the BM from HSPCs via a common lymphoid progenitor, completing BM maturation at the CD56^brightCD16⁻ stage in humans, before entering the circulation and further maturing to CD56^dimCD16⁺ NK cells in secondary lymphoid tissues (Freud et al., 2014). However, we previously reported that the major circulating phenotypic subpopulation of rhesus macaque cytotoxic CD16⁺ NK cells are clonally distinct from T, B, and myeloid cells in barcoded rhesus macaques (Wu et al., 2014; Koelle et al., 2017). In our current study, BM CD16⁺ NK cells were clonally identical to PB CD16⁺ NK cells and clonally distinct from local CD34⁺ HSPCs, suggesting ongoing production of CD16⁺ NK cells independent of BM HSPCs, outside the BM. BM CD16⁺ NK cells almost 100% stained with the intravascular anti-CD45 antibody, suggesting that BM CD16⁺ NK cells were largely intravascular cells present in blood flowing through the BM, rather than NK cells that had homed back to the BM from the PB more than 4 h previously. Ongoing studies are focused on localizing the site of production of this clonally unique NK subset in tissues and organs with substantial populations of NK cells, including tonsils, LNs, liver, spleen, salivary glands, and uterus (Renoux et al., 2015; Björkström et al., 2016). In contrast, CD56⁺ NK cells were clonally related to local CD34⁺ HSPCs in the BM and were almost 100% BM resident by anti-CD45 intravascular staining, implying local production.

Perhaps the most unexpected finding our clonal tracking model provided was evidence for local T cell production in the BM early after transplantation. In the classic T cell developmental model, T cell precursors migrate from the BM to the thymus, where T cell receptor rearrangement and commitment to T cell functional subsets occur, as well as deletion of autoreactive cells (Koch and Radtke, 2011). Mature T cells then exit the thymus, circulate in PB, or take up residence in secondary lymphoid tissues (Shah and Zúñiga-Pflücker, 2014). However, evidence for extrathymic T cell generation or final maturation has been reported, suggesting that these extrathymic T cells may have different functions from canonical thymic-dependent T cells (Tsark et al., 2001; Wang et al., 2001; García-Ojeda et al., 2005; Nonaka et al., 2005).

In our model, we have evidence for local production of mature CD3⁺ T cells from CD34⁺ HSPCs during the earliest stages of T recovery after myeloablative transplantation. At these time points, the overall output of new GFP⁺ barcoded T cells from transplanted HSPCs was very low, likely because of thymic dysfunction resulting from irradiation. By later time points, these BM locally produced T cells almost disappeared, implying they arose in the context of hematopoietic stress immediately after transplantation, perhaps before any degree of thymic recovery. These locally produced T cells included populations that were either CD4 or CD8 positive, expressed tissue-resident T cell markers such as CCR5, and were primarily TCRαβ positive. Many of these cells coexpressed CD56, which can occur under situations such as chronic inflammation or aging (Michel et al., 2007; Lemster et al., 2008). The function, diversity, and specificity of locally produced BM T cells deserve further investigation.

Production of homogeneously BM-distributed, presumably thymic-dependent, T cells increased over time, and barcode diversity and distribution patterns in circulating T cells mirror those in other hematopoietic lineages, as shown in the current study and our prior study (Koelle et al., 2017). This polyclonality is in contrast to a recent study tracking clonal diversity of human T cells produced from barcoded cord blood engrafted in immunodeficient mice, reporting highly restricted clonal patterns in thymic precursors and mature human T cells compared with HSPCs, suggesting a thymic niche bottleneck (Brugman et al., 2015). These differences may be caused by xenograft-specific factors.

In conclusion, clonal tracking of HSPCs and their progeny both geographically and temporally, across hematopoietic lineages, and in a model organism with great relevance to human biology has provided new insights, some with clinical relevance. Our findings help explain the highly heterogeneous cellularity observed even in individual human BM biopsies, particularly after transplantation or recovery from cytotoxic therapies. Barcoded human primary leukemia samples transplanted into immunodeficient mice also demonstrated asymmetric engraftment levels across different BM locations (Belderbos et al., 2017). Sampling of BM samples for acquired somatic mutations, chimerism, or cytogenetic abnormalities may be misleading, and any assays that can instead be performed on PB may provide a better assessment of overall allele frequencies or disease burden. It seems possible that leukemic or MDS stem cells and their progeny may also spread contiguously, and would be an interesting area for further investigation.

The barcoded lentiviral vector consists of the backbone pCDH (Systems Biosciences) expressing the CopGFP marker gene, followed by a 6-bp library identifier and a 27- or 35-bp highly diverse DNA barcode sequence (Lu et al., 2011). The vectors were produced with the χHIV packaging system, modified from an HIV packaging system, allowing efficient transduction of rhesus macaque cells (Uchida et al., 2009, 2012). Diversity was validated using Monte Carlo simulation to determine the number of CD34⁺ target cells able to be transduced with a library, resulting in a greater than 95% probability that more than 95% of barcodes represent single transduction cells (Wu et al., 2014; Koelle et al., 2017).

Animal protocols were approved by the National Heart, Lung, and Blood Institute Animal Care and Use Committee. Rhesus macaque PB stem cell mobilization, collection, and CD34 immunoselection were performed as described (Donahue et al., 2005; Wu et al., 2014). CD34⁺ HSPCs were transduced once with the high-diversity barcoded libraries at a multiplicity of infection of 25, targeting a transduction efficiency of ∼30%, ensuring that the majority of transduced HSPCs contain only one barcode per cell (Wu et al., 2014; Koelle et al., 2017). Cells were reinfused into the irradiated (500 cGy, 2×) autologous macaque 24 h posttransduction. Table 1 gives the transduction and transplantation parameters for all rhesus macaques included in the current article. Details of transplantation and the clonal patterns in PB from animals ZH33, ZG66, ZJ31, and ZH19 have been previously reported (Wu et al., 2014; Koelle et al., 2017).

The rhesus macaque–specific anti-CD45 monoclonal antibody MB4-6D6 was produced by Miltenyi for intravenous administration, using an azide-free, formaldehyde-free endotoxin <0.5 U/mg protein. This antibody was conjugated to Alexa Fluor 647. The anti-CD45–Alexa Fluor 647 compound was administered intravenously at dose of 30 µg/kg 4 h before tissue harvesting.

10–15-ml BM aspirates were obtained from the posterior iliac crest or ischial tuberosities. Cells from LNs were isolated by gentle mechanical disruption in the presence of 1 mg/ml collagenase D (Roche; Seggewiss et al., 2007). BM and PB samples were separated into Mono and Gr fractions via Ficoll gradient centrifugation. The Gr pellet was collected after red cell lysis with ACK buffer (Quality Biological), with Gr purities assessed by cytospin of >95%. BM mononuclear cells were passed over an immunoselection column to purify CD34⁺ HSPCs as described (Donahue et al., 2005). BM CD34⁻ cells flowing through the column, PB mononuclear cells, or LN cells were stained with lineage-specific antibodies and sorted for CD3⁺ T cells, CD3⁺CD4⁺ and/or CD8⁺ T cells, CD20⁺ B cells, CD3⁻CD20⁻CD14⁺ Mono, and CD3⁻CD20⁻CD14⁻CD56⁺/CD16⁺ NK cells (Fig. S1). A list of monoclonal antibodies used can be found in Table S1.

Genomic DNA was extracted with the DNeasy Blood & Tissue kit (Qiagen) and quantified by NanoDrop. 200 ng (ZH33, ZG66, ZJ31, and ZL40) or 500 ng (ZH17, ZH19, and ZJ38) DNA underwent PCR with primers bracketing the barcode and multiplex sequencing as described (Wu et al., 2014). Equal amounts of gel-purified product from individual samples were pooled together for multiplex sequencing (Illumina HiSeq 2000).

Sequencing output was processed using custom Python code supplied with our prior publications (Wu et al., 2014; Koelle et al., 2017) and available on Github (www.github.com/dunbarlabNIH/). Data analysis and plot generation was done using R code, ggplot2, corrplot, and reshape2 software, and Python code and associated R functions used for analyses are also freely available on Github.

Fig. S1 shows the purification strategies for HSPCs and specific hematopoietic lineages and the GFP marking for PB lineages. Fig. S2 shows clonal patterns of the top 100 high-contributing clones from left versus right BM CD34⁺ HSPCs and Gr from BM and PB at 3.5–7 mo posttransplantation. Fig. S3 shows clonal diversity of BM CD34⁺ HSPCs and Gr from BM and PB over time after transplantation. Fig. S4 shows clonal heat map of geographic localization of top 10 contributing clones in animal ZJ31 at 8.5 mo. Table S1 lists antibodies used in this study.

We thank Naoya Uchida for χHIV plasmid, and Keyvan Keyvanfar (National Heart, Lung, and Blood Institute [NHLBI] FACS Core and DNA Sequencing and Genomics Core, Bethesda, MD) for technical assistance. We thank Allen Krouse, Mark Metzger, Sandra Price, Aylin Bonifacino, Barrington Thompson, and the NHLBI Primate Facility for animal care and sample procurement.

This research was supported by the intramural research programs of the National Heart, Lung, and Blood Institute and the Vaccine Research Program of the National Institute of Allergy and Infectious Diseases. D. Yang was funded by the scientific research training program for young talents from China.

The authors declare no competing financial interests.

Author contributions: Conceptualization: C. Wu, M. Roederer, C.E. Dunbar; investigation: C. Wu, D.A. Espinoza, S.J. Koelle, B. Li, D. Yang, X. Fan, R.E. Donahue; formal analysis and software; D.A. Espinoza, S.J. Koelle, B. Li, R. Lu; resources: R. Lu, M. Roederer; methodology: E.L. Potter, R. Lu, M. Roederer; visualization: C. Wu, D.A. Espinoza; writing: C. Wu, C.E. Dunbar.