Single-cell analysis defines the divergence between the innate lymphoid cell lineage and lymphoid tissue-inducer cell lineage.

The precise lineage relationship between innate lymphoid cells (ILCs) and lymphoid tissue-inducer (LTi) cells is poorly understood. Using single-cell multiplex transcriptional analysis of 100 lymphoid genes and single-cell cultures of fetal liver precursor cells, we identified the common proximal precursor to these lineages and found that its bifurcation was marked by differential induction of the transcription factors PLZF and TCF1. Acquisition of individual effector programs specific to the ILC subsets ILC1, ILC2 and ILC3 was initiated later, at the common ILC precursor stage, by transient expression of mixed ILC1, ILC2 and ILC3 transcriptional patterns, whereas, in contrast, the development of LTi cells did not go through multilineage priming. Our findings provide insight into the divergent mechanisms of the differentiation of the ILC lineage and LTi cell lineage and establish a high-resolution 'blueprint' of their development.

Introduction

Innate lymphocytes lack B or T cell receptors and exert effector functions at mucosal barriers [1, 2]. These populations segregate into three general groups based on the expression of the transcription factors T-bet, GATA-3 and RORγt,. However, there is considerable heterogeneity among T-bet-expressing group 1 lymphocytes which comprise conventional (or classical) NK cells (cNKs), ILC1s, and tissue-resident NK cells, and in RORγt-expressing group 3 lymphocytes, which comprises CCR6+ lymphoid tissue inducer (LTi) cells and CCR6− ILC3s. In addition, some plasticity has been reported among CCR6− ILC3s which can upregulate T-bet and acquire group 1 properties [3], and among some populations of ILC2s which can acquire group 3 properties [4].

Lineage tracing and cell transfers have suggested that ILC1s, ILC2s and ILC3s, but not LTi cells or cNKs, were derived from a common dedicated precursor, the ILCP, characterized by expression of the transcription factor PLZF [5]. Similar to the LTi precursor (LTiP), the ILCP originates from an α4β7+ lymphoid precursor which was itself derived from the common lymphoid precursor (CLP). The Id2hi fraction of α4β7+ lymphoid precursors, termed the common helper innate lymphoid precursor (CHILP), is a heterogeneous population containing the PLZF-expressing ILCP as well as precursors to LTi cells [6], but it was not determined whether the CHILP population contained a common precursor to both ILCs and LTis, or separate precursors to these two lineages. A study has suggested that cNKs might originate from an earlier Id2loCXCR6+ fraction of α4β7-expressing lymphoid precursors (αLPs) [7]. Thus, the developmental relationships between these lineages remain incompletely established.

Several transcription factor genes including Nfil3, Tox, Id2, Gata3, Tcf7 and Zbtb16 (encoding PLZF) [ are required for the development of all or several of these innate lineages, suggesting an impact at a common precursor stage. However, partial rather than complete defects were often reported in mice lacking these transcription factors, suggesting significant redundancy and complexity within this early transcriptional network. Other transcription factor genes were found to selectively impact individual ILC lineages, such as Rora, Bcl11b and Gfi1 for ILC2 [17-19], suggesting more distal effects in the ILC differentiation pathway. A precise understanding of the general hierarchy of expression of these factors is missing, however, limiting the design and interpretation of mechanistic studies aiming at dissecting their interplay.

Here, we used cultures of single cells purified from the fetal livers of a Zbtb16-GFP reporter strain to precisely define the differentiation stages between CLP and ILC and to identify the stage of divergence between ILC and LTi cells. We also performed multiplex quantitative single-cell transcriptional analysis of 100 lymphoid genes to characterize the complexity of molecular events associated with this differentiation pathway. We derived a high-resolution map and inferred a precise ordering of transcription factor induction, defined new stages of development before and after PLZF expression, and identified the stage of bifurcation between ILC and LTi cell lineages. Notably, transcriptional priming for the different cytokine effector programs occurred at the ILCP stage itself through multilinage priming. In contrast, the LTiP proceeded to directly acquire its type 3 program without undergoing mixed transcriptional priming. Altogether, these findings further define the dichotomy between ILCs and LTis and provide new insights into the stages and mechanisms of their development, and the interplay of transcription factors that direct their differentiation.

Results

Bifurcation of αLP into ILCP and LTiP

Lin−IL-7Rα+ fetal liver cells (where Lin is a cocktail of antibodies against CD3ε, TCRβ, CD19, CD11c, GR-1, Ter119 and NK1.1), include the CLP, which is identified by a Flt3+ α4β7− profile, and a cell population expressing α4β7, which contains the precursors to innate lymphoid lineage and is thought to arise from the CLP (). This α4β7-expressing population was originally termed αLP, but the previous identification of ILCP and LTiP within αLP prompted us to limit the designation of αLP to the α4β7-expressing cells that are neither ILCP nor LTiP, but include their precursors [2, 20]. We further subdivided this αLP into Flt3+ and Flt3− subsets to identify early and late precursors, respectively. These populations are shown in the flow cytometry analysis of E15 fetal liver cells from the Zbtb16-GFP reporter strain (). Among the Lin−IL-7Rα+α4β7+ population, the ILCP was identified by GFP expression and the LTiP by its GFP−CXCR5+ profile [5, 21]. The LTiP was clearly distinguishable from the other subsets in this staining, although a small fraction of it seemed to express low amounts of GFP (). As expected, only the LTiP co-expressed CCR6 and CD4 ().

Both LTiPs and ILCPs could already be observed at low frequencies among Lin−IL-7Rα+ α4β7+ fetal liver cells at E12, the earliest time point in our analysis, and their absolute numbers increased 5-10 fold by E15 (). To establish the lineage relationships between αLP, ILCP and LTiP, we performed single-cell cultures of αLP subsets (Flt3+ and Flt3−) on OP9 stromal cells in the presence of IL-7 and stem cell factor (SCF), as described [5]. We scored ILC1, ILC2 and ILC3 colonies by high expression of NK1.1, ICOS and α4β7 respectively, as previously reported [5]. We distinguished LTi cells from ILC3 colonies by expression of CD4, which is found in nearly half the LTi cells but not in ILC3s (), although this method underestimates the real frequency of LTi cell colonies by approximately two-fold. We could not distinguish cNKs from ILC1s by eomesodermin expression, because this transcription factor was induced in ILC1s in our culture conditions (data not shown). However, because fetal progenitors do not produce cNKs [22], we counted all NK1.1+ colonies as ILC1s.

ILCPs gave rise nearly exclusively to single and mixed colonies of ILC1s, ILC2s and ILC3s, consistent with this cell-type being a common precursor to these lineages, as previously reported [5] (). In contrast, Flt3+ and Flt3− αLPs also generated a sizeable proportion of LTi cells, indicating that they contained precursors to the LTi cell lineage. Importantly, many wells containing LTi cells also included ILC1s or ILC2s (the presence of ILC3s could not be ascertained in wells containing LTis), indicating that a single Flt3+ or Flt3−αLP precursor could generate both LTi cells and ILC lineages. Thus, while the previously reported CHILP [6] contained a heterogeneous mixture of ILCP and other precursors, our observations identify the Flt3−αLP as the common proximal precursor, at the single cell level, to both ILCP and LTiP.

Notably, when the cultures were performed with OP9 stromal cells lacking Notch ligand, many fewer ILC2 colonies were derived from Flt3+ αLPs or Flt3− αLPs than from ILCPs. (, bottom). However, this defect was largely corrected in OP9-DL1 cultures (, bottom). This finding was consistent with the notion that a Notch signal is essential for ILC2 differentiation [14, 17], and further established that it must be delivered early at the αLP stage, rather than late at the ILCP stage.

Altogether, these experiments demonstrated that the fetal liver αLP pathway was dedicated to the formation of ILC and LTi lineages. They also identified the late Flt3−αLP as the common proximal precursor to the ILCP and LTiP.

Single-cell multiplex qRT-PCR analysis of ILC precursors

We performed transcriptional analysis of single cells along the fetal pathway linking αLP to ILCP and LTiP with the Biomark Fluidigm multiplex qRT-PCR system, using a panel of 100 lymphoid development factors including transcription factors, cytokine receptors and chemokine receptors. Data from 339 single cells, including 157 αLPs, 168 ILCPs, and 14 LTiPs were compiled from two independent single-cell sorting experiments using pooled E15 fetal livers. After filtering by housekeeping gene expression, we considered the transcriptional profiles from 299 single cells for downstream analysis.

Unsupervised hierarchical clustering analysis of these transcriptional profiles defined groups of cells that appeared to be in distinct developmental stages (). Consistent with the notion that acquisition of PLZF signified a distinctive developmental transition for ILC, we found strong, visually evident, separation between transcriptional profiles of the αLP (left clusters) and ILCP populations (right clusters). We used contiguous branches of the hierarchical clustering dendrogram to define 3 clusters of predominantly αLPs (AI-III), 1 cluster composed of mixed αLPs and ILCPs (B), and 4 clusters of predominantly ILCPs (CI-IV), based on our current understanding of ILC development. All clusters were significantly distinct from one another and were composed of a substantial number (>20) of similar cells (). One additional cluster of about 20 cells, all expressing conspicuously high levels of Il1rl1 encoding the IL-33 receptor chain IL-33Rα, was removed from the study because it was unrelated to the other clusters and, instead, seemed to represent contaminating mast cell precursors expressing low amounts of α4β7 and PLZF ().

Thus, this analysis identified further heterogeneity amongst precursors and generated a blueprint of their temporal sequence during ILC development.

Early developmental transitions prior to PLZF expression

To facilitate the examination of clusters, we generated a condensed heat map of all 299 single cells, limited to a set of 20 genes selected for their known function in innate lymphocyte differentiation (). Consistent with αLPs being early precursors to ILCPs and LTiPs, there was sparse expression of transcription factors and cytokines specific for these lineages in the A clusters. For example, Zbtb16, Tcf7, Gata3, Bcl11b, Cxcr5, Rorc, Rora and Tbx21 were not found in A clusters. In contrast, the A clusters expressed transcription factors that were implicated in early ILC and LTi development, including Id2, Tox, Nfil3, Sox4, Ets1, Runx1, Tcf7 and Zbtb16 (). This conclusion was confirmed by plots depicting the average mRNA expression per cell (), or the percentage of cells expressing these transcription factors within each cluster (). Notably, a clear temporal pattern of expression could be inferred from these graphs. Thus, cells in cluster AI expressed low amounts of Id2, but most lacked Tox and Nfil3. Nearly all cells in clusters AII-III, however, expressed Tox and Nfil3. Cells in clusters AI-II lacked expression of Ets1, Sox4 and Runx1, whereas a majority of cells in AIII expressed these factors. Finally, Tcf7 and Zbtb16 were not expressed in A (αLP) clusters, but were widely expressed in cluster B and the C (ILCP) clusters. We measured low Id2 expression across all A clusters, with a tendency towards more frequent and higher levels of expression in B and C clusters, in line with the suggestion that increased Id2 correlated with innate lineage commitment [6].

Thus, the Biomark analysis suggested that the temporal patterns of expression of these transcription factors were precisely regulated. Furthermore, unlike Tox, whose expression remained prevalent, Nfil3 expression was ultimately reduced in ILCP clusters, consistent with its temporally limited requirement as suggested by late gene ablation experiments [8].

Bifurcation between ILC and LTi branches

Cluster B comprised the most mixed representation of αLP and ILCP and appeared to be a developmental transition state linking early developmental events and ILC vs. LTi cell lineage specification. Notably, cells in cluster B expressed Tcf7 and Zbtb16, but largely lacked expression of ILC lineage specific factors. In fact, the levels of expression of Tcf7 and Zbtb16 in cluster B were lower than those observed in ILCP clusters, corroborating the conclusion that cluster B corresponded to a transitional stage (). A similar pattern was found for Id2, suggesting developmental continuity between αLP and ILCP clusters, with cluster B likely including the first cells to acquire Zbtb16 prior to ILC lineage specification. Although the induction of Tcf7 and Zbtb16 seemed nearly synchronous, a fraction of cluster B cells expressed Tcf7 without Zbtb16 ( and ), suggesting that they may include cells destined to the LTi cell lineage. In fact, detailed analysis of Tcf7 cluster B cells showed that most of these cells had a tendency to express factors associated with the LTi cell lineage, including Rorc, while conspicuously lacking Rora, which is expressed in more mature LTiPs (). In contrast, cluster B cells that expressed both Tcf7 and Zbtb16 rarely expressed Rorc at that stage, and instead displayed some expression of factors associated with ILC lineages, such as Gata3, while lacking Rora, suggesting that they were less mature than most ILCPs. Thus, we concluded that cluster B represented the stage of bifurcation of the αLP into the LTi cell and ILC lineages.

ILC lineage differentiation originates in ILCP

Cells in the ILCP clusters (CI-IV) showed signs of ILC maturation and expressed ILC lineage-defining transcription factors and cytokines (). These cells almost universally expressed high levels of Zbtb16, Tcf7, and Id2, and additional transcripts further distinguished them from the developmental intermediary cluster B. Cells in ILCP clusters almost all expressed Ets1 and Rora, and tended to have low Nfil3 expression. The ILCP clusters delineated by hierarchical clustering likely reflected sequential stages of ILC maturation more so than association with particular ILC lineages. Cells in clusters CI-II expressed high levels of Sell and Sox4 that diminished in clusters CIII-IV. Conversely, cells in clusters CIII-IV express higher levels of Ets1 and Il7r than those in clusters CI-II. Characteristically, cells in clusters CIII-IV more frequently expressed ILC lineage-defining transcription factors, including Tbx21, Bcl11b, and Rorc. The fact that we observed ILC lineage-defining transcription factors in multiple ILCP clusters indicated that we were likely capturing developmental stages of ILC differentiation. Further, cells that expressed ILC lineage-defining transcription factors in later clusters tended to express other lineage-associated cytokines and transcription factors and so appeared more differentiated.

ILCPs undergo multilineage transcriptional priming

To better illustrate the range of cells among the ILCP clusters that appeared to be differentiating toward ILC lineages, we compiled a subset of cells that expressed ILC lineage-associated factors (). Specifically, we identified cells with ILC1, ILC2, or ILC3 markers by their expression of Tbx21, Bcl11b, or two of Cxcr5, Rorc, Lta and Ltb, respectively. Among this population, we observed many cells co-expressing markers of different ILC lineages, for example Tbx21 and Bcl11b, or Tbx21 and Cxcr5 or Tbx21 and Rorc, and even cells expressing multiple markers of the three different lineages. In total we enumerated 60 cells expressing lineage differentiation markers out of 120 total cells in ILCP clusters, of which nearly one third (n=19) co-expressed markers of different lineages.

The multilineage priming of these cells was substantiated by the expression of other factors. In particular, the simultaneous expression of Bcl11b and Icos was associated with expression of canonical ILC2 genes, including Il1rl1 and Il13. We also observed increased expression of Gata3 in Bcl11b-expressing cells, although Gata3 was expressed throughout the ILCP without strict lineage association (as previously shown by flow cytometry [5]). Across cells in ILCP clusters, we observed relatively strong correlation between expression of Rxrg and Igf2 with Bcl11b. The expression of Tbx21 was associated with higher expression of Il2rb, and in some cases with expression of Eomes and Ncr1 (encoding NKp46). We also found that Irf8 appeared to preferentially associate with Tbx21 in differentiating ILCP. Finally, we found that expression of Rorc was associated with high levels of Cxcr5, as well as with the expression of Lta, Ltb, and Rankl. Generally, Il2rb and Cxcr5 expression was observed throughout the ILCP clusters, with higher levels and frequencies of expression associated with other ILC1 and ILC3 factors. There was substantial overlap between the expression of ILC1 and ILC3 factors, in particular between the expression of Tbx21 and Cxcr5, Lta, Ltb, and Rankl.

These results were in sharp contrast with the expression profile of LTiPs in the same fetal livers ( and , right panels). These cells showed no evidence of multilineage priming. In particular they conspicuously lacked expression of Gata3, Bcl11b or Tbx21, further supporting the notion that they followed a distinct lineage differentiation pathway.

Altogether, the expression of multiple lineage factors in ILCP clusters was suggestive of a multilineage potential. This further confirmed that cells navigating ILC lineage decisions were contained in the ILCP. Moreover, it suggested that during the initial lineage decisions, the relationships between distinct ILC lineage factors were more complicated than might have been anticipated from studies in T helper cells, for example.

Transcriptional programs correlate with lineage decisions

To functionally evaluate the significance of these transcriptional programs in ILCPs, we tested their lineage potential in single-cell cultures. We sorted fetal ILCP subsets based on their expression of CD122, ICOS or CXCR5 (), and studied their progeny in single-cell cultures. As shown in (top row), these subsets showed some bias towards the corresponding ILC1, ILC2 and ILC3 programs, respectively. Thus, the earliest precursors with a bias to ILC1 or ILC2 lineages could be discerned at the ILCP stage by their CD122+CXCR5− or ICOShi profiles, repectively. However, the biases were incomplete, illustrating the residual multipotency of these populations. For example, while ICOShi ILCPs generated mostly ILC2 colonies, they also gave rise to some single ILC1 and ILC3 colonies, as well as to dual and triple colonies. In contrast, CD122+CXCR5− ILCPs mostly generated ILC1 colonies, but also some single ILC3 colonies as well as dual colonies. Even the dual colony-producers maintained a corresponding bias as, for example, all dual colonies originating from ICOShi precursors included an ILC2 colony, and those originating from a CD122+CXCR5− precursor included an ILC1 colony (, bottom row). Together with the multilineage priming observed at the transcriptional level, these observations implied that lineage polarization was proceeding by stepwise restriction of alternative programs from multipotential precursors, ultimately leading to the canonical polarized lineages. Using intracellular staining for transcription factors, we further confirmed that a substantial fraction of fetal ILCPs coexpressed T-bet, RORγt and GATA3 or ICOS, indicative of multilineage priming, whereas mature ILCs were strictly limited to their lineage-specific transcription factors (). This scenario is different from current models of polarization in T helper cells, where acquisition of cytokine effector programs does not typically involve intermediates with mixed lineage patterns.

Discussion

Using the Zbtb16-GFP reporter system, we identified the stage where the αLP bifurcates into ILCP and LTiP during early fetal lymphopoiesis and used single-cell multiplex qPCR analysis to produce a high-resolution map of innate lymphocyte lineage development. It should be stressed that the proposed developmental progression, discussed below, is inferred based on the continuity of gene expression between clusters, but remains to be experimentally demonstrated.

We found that Zbtb16 and Tcf7 were simultaneously upregulated in cluster B at the bifurcation of ILC and LTi lineages, with Tcf7 marking both lineages and Zbtb16 identifying the ILC lineage. A fraction of these Tcf7+Zbtb16− precursors expressed Rorc and seemed closer to LTiP than to ILCP. Notably, cluster B cells lacked expression of Rora, which is characteristically induced later in development, further supporting the conclusion that cluster B cells navigate the bifurcation between LTi and ILC lineages.

The A clusters prior to Zbtb16 and Tcf7 induction must therefore represent earlier precursors, based on the expression of Id2, Tox, Nfil3, Runx1 [7–10, 16, 23–25] , which have been implicated in both lineages. Id2 appeared to be the earliest transcription factor expressed at substantial level at stage A-I, when Nfil3 and Tox transcripts were barely detected. Nfil3 and Tox levels rapidly ascended during A-II and reached maximum levels at A-III. These findings stand in apparent contrast with a recent report that Nfil3 could bind to and induce Id2 and that Nfil3 ablation could be complemented by Id2 [8], but it is consistent with the presence of Id2 transcripts reported in arrested precursors lacking Nfil3 or Tox [8, 10]. Nfil3 might therefore exert a positive feedback loop rather than a primary trigger for Id2 expression. Nfil3 and Tox also appeared nearly simultaneously induced. Thus, a report that Nfil3 could induce Tox and that Nfil3 ablation could be complemented by Tox [7] might also reflect a positive feedback mechanism.

Additional transcription factors that were induced early after Id2/Tox/Nfil3, but before Zbtb16/Tcf7, include Sox4, Runx1 and Ets1. Although their functions in ILC development are currently unknown, Sox4 was previously associated with the development of fetal γδ T cells in conjunction with Tcf7 [26], whereas Runx1 was previously shown to be important for NK cell and LTi cell development [23], and Ets1 was associated with NK cell development [27]. Thus, these factors are likely integral components of the early transcriptional network of innate lymphocytes.

The C clusters defined several stages of ILCP with C-I associated with the induction of Rora, which was maintained in all ILC lineages, although it only seems to be required for ILC2 [17, 28]. C-II was associated with the expression of Gata3, which was maintained at high level during the remaining ILCP stages and ultimately is downregulated in mature ILC1 and ILC3 [5]. Therefore, transient but high level of Gata3 expression is an intrinsic developmental event that may distinguish ILC3s from LTi cells. Clusters C-III and C-IV were marked by the emergence of lineage-specific factors. Intriguingly there was a clear and extensive pattern of coexpression of factors from the three different lineages in nearly a third of these cells, indicating that a phase of multilineage transcriptional priming preceded the polarization of mature lineages. In sharp contrast, the LTiP did not go through this singular process and directly expressed Rorc and other LTi cell lineage attributes without coexpression of alternative ILC lineage markers. The extent of multilineage priming revealed in our study goes well beyond the coexpression of low levels of T-bet and RORγt by GATA3+ ILC precursors, which was previously detected by FACS analysis in the fetal intestine [29]. Indeed, we found that it included an extended list of canonical markers of all three lineages. ILCPs sorted based on the expression of some of these lineage differentiation markers showed a relative bias in differentiation towards the corresponding lineages after single-cell culture, although these subsets also maintained some multipotency as shown by the generation of mixed colonies. These findings suggested that the terminal differentiation of ILCP into polarized ILC lineages was more complex than previously envisioned. Instead of direct acquisition of one of three programs, ILC precursors first activated multiple effector programs simultaneously and then progressively turned off alternative lineage priming. In the absence of exogenous polarizing cytokines, such developmental strategy can take advantage of the well-established antagonistic effects between these programs [30]. In ILCP, external or internal inputs that could influence lineage decision include retinoids, which were shown to favor type 3 innate lymphocytes [31, 32], or Notch signals which were needed for type 2 ILCs [14, 33]. Multilineage transcriptional priming has been proposed as a more general strategy for developmental decisions, such as the differentiation of various hematopoietic lineages from a common progenitor [34-37]. It has not been widely reported for the differentiation of Th1/Th2/Th3 programs in CD4 T-helper cells, which instead relies on polarizing cytokines released during infection or allergy [30], although mixed transcriptional patterns have been observed in cells cultured under mixed cytokine conditions [38] .

Our studies further emphasized that, while sharing many functional properties, ILC3s and LTi cells had distinct precursors and different developmental histories. Only ILC3 precursors transited through a PLZFhigh GATA3high stage with multilineage priming at the single cell level. These results should encourage studies aimed at further understanding the different functions of LTi cells and ILC3s, as illustrated, for example, by the report that ILC3s and LTi cells have a different topographical location in the lamina propria [39].

In summary, this study allowed the identification of several new developmental transitions and transcription factors sequentially associated with the lineage progression of innate lymphocytes. In particular, our results characterized in cellular and molecular details the previously poorly defined bifurcation of lymphoid precursors into LTiP and ILCP, and demonstrated differential mechanisms of acquisition of polarized effector programs by these lineages.

Materials and Methods

Mice

Zbtb16-GFP reporter mice were generated in our lab as previously described and contain an IRES after the last exon of Zbtb16, followed by an eGFP-Cre fusion gene 40. For fetal ontogeny experiments, the morning a vaginal plug was identified was counted as embryonic day 0 (E0). Mice were housed in a specific pathogen free environment at the University of Chicago and experiments were performed in accordance with the guidelines of the Institutional Animal Care and Use Committee.

Preparation of cell suspensions

Fetal livers were mechanically dissociated through a 70μm cell strainer and resuspended in HBSS (Gibco) containing 0.25% BSA (Sigma-Aldrich) and 5mM sodium azide (Sigma-Aldrich).

Flow Cytometry

Cell suspensions were pre-incubated with BD Fc Block for 10 minutes on ice. Fetal liver cells were pre-enriched for α4β7-expressing cells unless otherwise indicated. For pre-enrichment of α4β7+ cells, fetal liver or bone marrow cells were stained with allophycocyanin (APC)- conjugated antibody and subsequently bound to anti-APC microbeads (Miltenyi Biotec). The cells were then enriched using the autoMACS (Miltenyi Biotec), positive selection double sensitive program. Lineage (CD3ε, CD11c, CD19, NK1.1, TCRβ, Ter119 and GR-1) depletion of fetal liver or bone marrow cells was accomplished by incubation of cells with biotinylated Lineage antibodies and then bound to streptavidin conjugated microbeads (Miltenyi Biotec) and separated using the autoMACS depletion sensitive program. Fluorochrome or biotin conjugated monoclonal antibodies (clone in parentheses) against mouse α4β7 (DATK32), CCR6 (29-2L17), CD3ε (145-2C11), CD4 (GK1.5), CD8α (53-6.7), CD11c (N418), CD19 (6D5), CD27 (LG.3A10), CD45.1 (A20), CD45.2 (104), CD49a (Ha31/8), CD90.2/Thy1.2 (53-2.1), CD122 (5H4), CD127/IL-7Rα (A7R34), CXCR5(L138D7), FcεRIα (MAR-1), GR-1 (RB6-8C5), ICOS (C398.4A), IL-33Rα/ST2 (DIH9), NK1.1 (PIK136), NKp46 (29A1.4), Sca-1 (D7), TCRb (H57-597), Ter119 (TER-119), GATA3 (TWAJ), RORγt (Q31-378), T-bet (4B10), mouse IgG1/κ (MG1-45), mouse IgG2a/κ (MG2a-53), rat IgG2b/κ (RTK45-30) were purchased from Biolegend, eBioscience, or BD Biosciences unless otherwise noted. The D-9 antibody against PLZF was conjugated to Alexa-Fluor 647 using the labeling kit from Molecular Probes Life technologies. For intracellular staining, cells were fixed and permeabilized using the Foxp3 Transcription Factor Staining Buffer Set (eBioscience). Cells were then blocked with unlabeled isotype control before addition of fluorochrome-conjugated anti-transcription factor antibodies. As negative control, a 20-fold excess of unlabeled anti-transcription factor antibody was added prior to the conjugated antibody (“cold” competition). Data was acquired on a LSRII (BD Biosciences) or sorted using a FACS Aria II (BD Biosciences) and analyzed using FlowJo software (Tree Star).

Single-cell cultures

Stocks of OP9 and OP9-DL1 stromal cells were a gift from J. C. Zúñiga-Pflücker. All experiments were performed in Opti-MEM with GlutaMAX (Gibco) containing 10% FCS (Gibco), 1% penicillin/streptomycin (Gibco), and 60mM 2-mercaptoethanol (Sigma-Aldrich) and maintained in a 37°C incubator (Thermo Scientific) with 5% CO2. Stromal cells were plated at 70% confluency. Before addition of lymphocytes, stromal cells were irradiated (1,500 rad) and culture media supplemented with murine IL-7 (25 ng ml−1; R&D Systems) and SCF (25 ng ml−1; BioLegend). Fetal liver lymphocytes were enriched for α4β7 by MACS and single-cell sorted onto 96 well plates containing stromal cells and cytokines as described above. Cultures were analyzed after 6 or 10 days of culture and only colonies with more than 10 CD45.2+ cells were considered.

Biomark

Cells were sorted in 96-well qPCR plates in 10µl of CellsDirect™ One-Step qRT-PCR Kit (Life Technologies), containing mixtures of diluted primers (0.05X final concentration, see Supplementary Table 1 and 2). Preamplified cDNA was obtained after reverse transcription (15’ at 40°C, 15’ at 50°C and 15’ at 60°C), and preamplification (22 cycles : 15’’ at 95°C, 4’ at 60°C), and diluted 1 : 5 in TE pH8 Buffer (Ambion). Sample mix was as follows : diluted cDNA (2.9 µl), Sample Loading Reagent (0.29 µl, Fluidigm), TaqmanUniversal PCR Master Mix (3.3µl, Applied Biosystem) or Solaris qPCR Low ROX Master Mix (3.3 µl, GE Dharmacon). Assay mix was as follows: Assay Loading Reagent (2.5µl, Fluidigm), Taqman (2.5µl, Applied Biosystem) or Solaris (2.5µl, GE Dharmacon). 48.48 or 96.96 dynamic array integrated fluidic circuit (IFC, Fluidigm) was primed with control line fluid, and the chip was loaded with assays (either Taqman or Solaris) and samples using an HX IFC controller (Fluidigm). The experiments were run on a BiomarkHD (Fluidigm) for amplification and detection (2’ at 50°C, 10’ for Taqman reagents or 15’ for Solaris reagents at 95°C, 40 cycles : 15’’ at 95°C, 60’’ at 60°C).

Analysis of single-cell multiplexed qPCR data

Independent single-cell qPCR experiments were performed for αLP, ILCP, and LTiP from two different pools of PLZF-GFPCre E15 fetal livers. Cells not expressing detectable levels of all three housekeeping genes Actb, Gapdh, and Hprt were removed from further downstream analysis. To appropriately compare measurements of distinct cells, expression levels of each gene for a given cell were normalized with respect to the average housekeeping gene expression level for that cell. Specifically the cycle threshold value for each gene, Ctg, was converted to the measure −ΔCtg = − Ctg, the difference between the average housekeeping gene cycle threshold and the gene cycle threshold. For quantification purposes, genes that were not detected were given a −ΔCt value of −30, close to the minimum value detected. Hierarchical clustering was performed using the Euclidian distance metric with complete-linkage agglomeration.

Visual inspection of the intercellular distances used for hierarchical clustering corroborates the strong distinction between αLP and ILCP clusters and the similarities between adjacent clusters in the prescribed developmental order (). We performed permutation analysis to empirically evaluate the significance of our clustering assignments. Specifically, we calculated the sum of square Euclidean distances from the group mean for each cluster and used the total sum of squares within groups (SSW) value to compare our clustering assignment to 10,000 randomly permuted clustering assignments. We found the SSW of our clustering assignment to be substantially lower than those of all the random permutations, implying P < 1e-4. We further directly compared all pairs of clusters using the same approach and similarly found that all clusters were significantly distinguishable from one another with P < 1e-4. Finally, we performed hierarchical clustering with each gene iteratively removed from our complete dataset to evaluate the stability of our clustering method. We found that hierarchical intercellular relationships were widely preserved under these perturbations. For instance, in ~85% of all cell pairs, the number of branches connecting these pairs of cells to their common branch point on the clustering dendrogram were changed by less than 2 on average. Hierarchical clustering and analysis of qPCR data was performed with custom scripts using the base packages in R (v3.1.2) and heatmap displays were generated using the NMF package (v0.20.6).

Statistical analysis

Two-tailed Student’s t-test was performed using Prism (GraphPad Software) *p<0.05; **p<0.001; ***p< 0.0001