key: cord-0283384-c2l4a7mf
authors: Elmentaite, R; Kumasaka, N; King, HW; Roberts, K; Dabrowska, M; Pritchard, S; Bolt, L; Vieira, SF; Mamanova, L; Huang, N; Goh Kai’En, I; Stephenson, E; Engelbert, J; Botting, RA; Fleming, A; Dann, E; Lisgo, SN; Katan, M; Leonard, S; Oliver, TRW; Hook, CE; Nayak, K; Perrone, F; Campos, LS; Dominguez-Conde, C; Polanski, K; Van Dongen, S; Patel, M; Morgan, MD; Marioni, JC; Bayraktar, OA; Meyer, KB; Zilbauer, M; Uhlig, H; Clatworthy, MR; Mahbubani, KT; Saeb Parsy, K; Haniffa, M; James, KR; Teichmann, SA
title: Cells of the human intestinal tract mapped across space and time
date: 2021-04-07
journal: bioRxiv
DOI: 10.1101/2021.04.07.438755
sha: 2e647c3710880c59437ea61911edb54c734ac996
doc_id: 283384
cord_uid: c2l4a7mf

The cellular landscape of the human intestinal tract is dynamic throughout life, developing in utero and changing in response to functional requirements and environmental exposures. To comprehensively map cell lineages in the healthy developing, pediatric and adult human gut from ten distinct anatomical regions, as well as draining lymph nodes, we used singlecell RNA-seq and VDJ analysis of roughly one third of a million cells. This reveals the presence of BEST4+ absorptive cells throughout the human intestinal tract, demonstrating the existence of this cell type beyond the colon for the first time. Furthermore, we implicate IgG sensing as a novel function of intestinal tuft cells, and link these cells to the pathogenesis of inflammatory bowel disease. We define novel glial and neuronal cell populations in the developing enteric nervous system, and predict cell-type specific expression of Hirschsprung’s disease-associated genes. Finally, using a systems approach, we identify key cell players across multiple cell lineages driving secondary lymphoid tissue formation in early human development. We show that these programs are adopted in inflammatory bowel disease to recruit and retain immune cells at the site of inflammation. These data provide an unprecedented catalogue of intestinal cells, and new insights into cellular programs in development, homeostasis and disease.

Intestinal tract physiology relies on the integrated contribution of epithelial, mesenchymal, endothelial, mesothelial, innate and adaptive immune, and neuronal cell lineages, whose relative abundance and cell networking fluctuate from embryonic development to adulthood. Factors contributing to these dynamics include gut function, environmental challenges and disease states that vary at different life stages. Adding further complexity is that the intestinal tract is formed of distinct anatomical regions that develop at different rates and carry out diverse roles in digestion, nutrient absorption, metabolism, and immune regulation in adulthood.

Single-cell analyses of intestinal organoids from iPSCs and human primary tissue have provided a unique opportunity to study the earliest stages of gut development. Although focused mainly on epithelial cells, they have led to identification of factors driving cellular differentiation 1,2 and novel markers of rare cell types 3 . Analysis of rare fetal intestinal tissues have further resolved the formation of villi-crypt structures and seeding of immune cells into the gut environment 4, 5 .

Similarly, our understanding of the cellular landscape of the adult gut is benefiting from single-cell technologies. We have previously reported regional differences in immune cell activation in the healthy human colon linked to variability in the microbiome composition 6 and expression of ligand-receptor pairs between cell populations has been used to infer cellular communication involving enteric neuronal cells subtypes 7 . Studies comparing inflammatory bowel disease (IBD) samples to healthy controls or non-inflamed tissue have allowed for the identification of a new stromal subtype that expands in disease 89 , clonal expansion of tissue-resident CD8 T cells in disease [10] [11] [12] and correlation between cellular response and clinical treatment 13 . While extensive work has been carried out to profile the intestinal tract at single-cell resolution, a holistic analysis of the gut through space (anatomical location) and time (lifespan) is lacking. expressed in the large intestines. Interestingly, small intestinal BEST4+ cells were marked by high expression of chloride channel and cystic fibrosis gene CFTR ( Figure 1F ), which we also observed at the protein level ( Supplementary Fig. 2F ). Over 60% of cystic fibrosis patients experience intestinal symptoms characterised by inflammation and sticky mucus 15 . This analysis elucidates distinctions between subtypes of BEST4+ epithelial cells found along the intestinal tract, and specifically implicates small intestinal BEST4+ cells in cystic fibrosis pathogenesis.

Together, this represents an in-depth cell census of the gut and a map of cellular diversity across regions and developmental time. This dataset enables us to determine previously unappreciated patterning of gene expression in the BEST4+ enterocyte compartment along the intestine.

Next, we further investigated the composition of the prenatal and postnatal gut epithelium in detail. This lineage is formed of major classes of absorptive cells, secretory and enteroendocrine cells (EECs) (Figure 2A Supplementary Fig. 3B ). NEUROG3, a transcription factor that mediates epithelial cell commitment to EEC lineage 16 , and the absence of mature neuropeptides suggests that the NEUROG3+ EECs may represent early committed EECs.

M cells, which specialise in transferring luminal antigen through the epithelial lining of the gut were very rare and seen only in pediatric and fetal samples within this dataset ( Figure 2B & Supplementary Fig. 3B ). Interestingly, expression of MHCII was restricted to postnatal cells ( Supplementary Fig. 3C ), suggesting that this expression is driven by the presence of the microbiome. All mature epithelial cell types were present from 12 PCW, except paneth cells, which appeared postnatally ( Figure 2B & Supplementary Fig. 3B ).

The gut epithelium represents a point of viral entry, and SARS-CoV-2 has been shown to enter enterocytes via the receptor, ACE2, and protease, TMPRSS2 17 . We used our developing gut epithelial dataset to investigate whether viral entry is possible at all life stages. Indeed, both ACE2 and TMPRSS2 were expressed by enterocytes in early development, childhood and adulthood, with TMPRSS2 also expressed by BEST4+ enterocytes and M cells ( Figure 2C & Supplementary Fig. 3D ). This highlights for the first time the ability for SARS-CoV-2 to infect the fetal gut epithelium. IBD pathogenesis is due to impaired intestinal epithelial cell function (as well as exacerbated immune functions). We sought to investigate how IBD risk may be linked to cell types across the anatomical regions of the gut, by analysing expression of IBD-associated epithelial cell type marker genes. We focused on the pediatric and adult epithelial compartments, since IBD onset is most prevalent at later life stages. Expression of IBD-GWAS genes was largely consistent across regions, with most variability arising from the presence or absence of cells per region ( Figure 2D ). In particular, paneth cells present only in the small intestine during health and well documented as metaplastic in colon of IBD patients 18 , showed specific expression of intelectin, ITLN2, also previously shown to be enriched in the small intestine 19 . Surprisingly, PLCG2, two missense variants of which have been linked to aberrant B cell responses in early onset IBD 20,21 and primary immune deficiency 22 , was specifically expressed by tuft cells.

To more deeply explore the relevance of PLCG2 expression in epithelial cells of IBD patients, we compared epithelial PLCG2 expression across different developmental stages and Crohn's disease. We observed a small increase in expression of PLCG2 in Crohn's disease epithelium ( Supplementary Fig. 3E ). To further confirm the role of epithelial PLCG2 in inflammation, we stimulated primary tissue derived intestinal organoids with inflammatory cytokines TNF-alpha or IFN-gamma and performed scRNAseq. While there was no morphological difference in stimulated organoid lines ( Supplementary Fig. 3F ), we confirm an increase in PLCG2 expression across stimulated organoid epithelial cells ( Supplementary  Fig. 3G, n=3) .

The discovery of expression of PLCG2 by tuft cells at higher levels than B and myeloid cell lineages ( Supplementary Fig. 3H -J) led us to wonder whether tuft cells specifically could be responsive to signals from immune cells. We therefore analysed the expression of genes involved in immune surveillance by these cells. We screened for expression of ITAM-linked receptors that may activate PLCG2 ( Supplementary Fig. 3K ). Among these, only FCGR2A, which is activated in response to IgG and shown to be expressed by epithelial cells in immunized mice 23 , was specifically expressed by a fraction (2.75%) of tuft cells ( Figure 2E ). Using the mouse as a model, we confirmed expression of Fcgr3 (ortholog of FCGR2A) at protein level by approximately 5% of small intestinal tuft cells ( Figure 2F & Supplementary  Fig. 3L ). Receptor tyrosine kinases were also expressed across tuft cells and other epithelial cell types, however since these are known to be mainly linked to PLCG1 activation, whether they are also responsible for PLCG2 activation in tuft cells is difficult to delineate ( Figure  2E ,G & Supplementary Fig. 3J ).

Tuft cells are known to be receptive to chemical signals with essential roles in pathogen responses of the gut [24] [25] [26] [27] . Indeed we confirm expression of the succinate receptor, SUCNR1, and downstream mediators GNAT3, PLCB2, ITPR2, PRKCA and TRPM5 ( Figure  2E ,G), which lead to IL-25-mediated type 2 immune responses. While mouse studies have also implicated TAS1R3 28 , we additionally show expression of TAS1R1, albeit at a lower level. TAS1R1 forms a heterodimer with TAS1R3 for umami sensing to activate the same type 2 immune signaling 26 (Figure 2E ,G).

Our in depth analysis of the gut epithelium adds resolution to its changing composition throughout life, showing regional patterning of IBD-GWAS genes in the adult gut and discovering the likely capacity of tuft cells to sense IgG through PLCG2 activation. This leads us to hypothesize that PLCG2 genetic variants may block the receptiveness of tuft cells and inhibit their downstream signaling to other immune cells, thus contributing to increased risk of IBD.

Similar to the epithelium, the enteric nervous system displayed a shift in development at 10-12 PCW. This provided us the opportunity to investigate the differentiation of human Enteric Neural Crest Cell (ENCC) progenitors. ENCC progenitors are derived from vagal or sacral neural crest cells that colonise the gut during embryonic development. ENCCs must balance proliferation and differentiation into glia and neurons without depleting the precursor pool, as premature differentiation of ENCCs leads to aganglionic colon 29, 30 .

Cycling ENCC progenitor cells, co-expressing RET, SOX10, EDNRB, were most evident in first-trimester development ( Supplementary Fig. 4A, B) . The remaining cycling populations (G2M stages) showed polarised expression of glial and neuronal genes ( Figure 3A -B; Methods), suggestive of distinct proliferating committed progenitors. Both of these committed progenitor populations expressed genes distinct from those observed in differentiated cells of ENCC-lineage. For example, committed glial progenitors expressed high levels of ASPM and SGO2 ( Supplementary Fig. 4C ). Amongst genes expressed by neuronal committed progenitors were GDF10, which acts through TGFβ signalling to activate axon sprouting 31 , multiple Notch ligands (DLL1, DDL3) and Notch enhancer MFNG, altogether implicating NOTCH and TGFβ signalling in human neuroblast differentiation ( Supplementary Fig. 4C ). The developmental transition from ENCC progenitor to glial or neuronal progenitors towards terminal cell types was supported by trajectory analysis ( Figure 3C , Supplementary Fig. 4D -F).

While ENCCs were abundant in the first-trimester, differentiated neuronal and glial cell types were enriched in the second-trimester of development ( Supplementary Fig. 4B ). All terminal glial cells expressed high levels of neural crest precursor genes including FOXD3, MPZ, CDH19, PLP1, SOX10, S100B and ERBB3, but lacked RET ( Supplementary Fig. 4C ). We further defined enteric glia (marked by BCAN, APOE, CALCA, HES5, FRZB), and three subsets of Schwann cell precursor (SCPs) all marked by PLP1, PMP22, CDH19, ITGA4 and CDH2, but varying in expression of TFAP2A, DHH and FABP7 (BFABP) ( Figure 4C ). All three subsets lacked expression of POU3F1 (OCT6), but highly expressed CDH19, a marker restricted to SCPs 32 , supporting their identity as precursors rather than mature Schwann cells. In particular, the first subset, which we labelled 'SCP', was enriched in the firsttrimester and marked by specific expression of TFAP2B and CXCL13. Both transcription factor TFAP2A (AP-2α) and its paralog TFAP2B (AP-2β) have been shown to be expressed in SCPs, but not immature Schwann cells, and regulate timing of SCP differentiation 33 . The second subset, termed 'Myelin Schwann-like', highly expressed myelin associated genes (MGP, MBP, DHH, NTRK2, SBSPON), while the third subset, termed 'Non-myelin Schwannlike', was enriched in TFAP2A, ELN, FABP7 and BMP8B expression ( Figure 3C , Supplementary Fig. 3C ). We visualised BMP8B+ Non-myelin Schwann-like cells in the myenteric plexus ( Supplementary Fig. 4G ), while DHH+ Myelin Schwann-like cells were found primarily in the mesentery ( Figure 3D ). SCPs give rise to 20% of enteric neurons with SCP-derived neurons emerging postnatally in mice 34 . Consistent with this observation, we did not observe differentiation of SCPs into neuronal subsets during embryonic and fetal development in humans ( Figure 3C , Supplementary Fig. 4D -F).

Neuronal differentiation extended from neuroblasts into two main neuronal branches, and one additional sub-branch ( Figure 3C ). Both main branches were found in all intestinal regions ( Supplementary Fig. 4H ). The first branching event divided branch 1, expressing BNC2 and RAMP1, from branch 2, expressing ETV1 and KCNJ5. Branch 2 further divided into sub-branches 2a and 2b, distinct in expression of gastrin releasing peptide GRP, transcription factor NEUROD6 in sub-branch 2a versus calcium binding protein SCGN, neurexophilin-2 NXPH2 and netrin-G1 NTNG1 in sub-branch 2b ( Supplementary Fig. 4C ). We visualise and validate the cells belonging to branch 1, 2a and 2b in the human myenteric plexus using smFISH in Figure 3E .

We sought to relate developing neuron branches to mature enteric neuron subtypes ( Supplementary Fig. 4J ). Branch 1 neurons, in addition to BCN2, expressed TSHZ2, ALK and HTR4 reminiscent of putative excitatory motor neurons 7 . Branch 2a neurons expressed genes for both putative inhibitory motor neurons (NOS1, GFRA1, PTGIR, KCNJ5, DGKB, and NPY) and putative secretomotor/vasodilator neuron markers (VIP, ETV1, SCGN, and UNC5B) 7 , suggesting this branch gives rise to both types of neurons. Branch 2b highly expressed TAC1, PENK, BMPR1B, ADRA2A and NXPH2 reminiscent of putative interneurons 7 . The developing branches expressed low levels of putative sensory neuron marker genes (CALCA, NMU, SST, DLX3), suggesting that this subtype develops later in gestation ( Supplementary Fig. 4J ).

Hirschsprung's disease is characterised by loss or impaired differentiation of the enteric nervous system and can manifest at variable locations of the distal bowel. To identify cell types involved in Hirschsprung's disease, we screened ENCC-lineage cells for the expression of Hirschsprung's disease-associated genes 35, 36 . We further included genes associated with Hirschsprung's disease phenotype in mice 37 , but currently lacking clinical association in humans.

The majority of Hirschsprung's disease-associated genes were expressed across multiple differentiating populations with varying intensity (Supplementary Fig. 4I ). For example, UBE2T was most highly expressed by committed proliferating progenitors, while other Hirschsprung's disease-associated genes were highly expressed in glial (ZEB2, COMT, NTRK3) or neuronal (TLX2, SEMA3D, L1CAM, HOXB5) compartments ( Supplementary Fig.  4I ). Expression levels of disease genes also varied between the three developing neuron clusters (branches 1, 2a and 2b). For example, ASCL1 was highly expressed in branch 1 and 2a, but not in branch 2b. RET was selectively absent from branch 1, and GFRA2 absent from branch 2a.

Interestingly, CDH2, L1CAM and RET were more highly expressed in colonic neuronal lineage cells compared to the equivalent cells in the small intestine, while the opposite was true for ATRX and ASCL1 ( Supplementary Fig. 4I ). We also observe unexpected expression of NTRK3 in pro-myelin SCPs ( Supplementary Fig. 4I ). Mice lacking Ntrk3 or its ligand Nt3 have reduced numbers of intrinsic sensory neurons 38 , suggesting that pro-myelin SCPs may support the differentiation of neurons during development.

To understand the contribution of mesenchymal cell subsets to pathogenesis of Hirschsprung's disease, we investigated the expression of key Hirschsprung's diseaseassociated ligands and their receptors 35, 36 across gut cells of all developmental ages ( Figure  3F ). Among these ligands was GDNF that activates RET and its co-receptor GFRA1, EDN3 and its receptor EDNRB3, and NRG1 and its receptor ERBB2. We observed expression of GDNF and EDN3 primarily by smooth muscle cells (SMC) and interstitial cells of Cajal (ICC) (described among stromal cells below; Figure 3F ). While in adult gut NRG1 is expressed by glial and multiple neuron subtypes 39 , in the fetal gut NRG1 was highly expressed in neuronal branch 2b. Its receptor ERBB2 was expressed across glial and stromal cells including ICC and SMC ( Figure 3F ,G). NRG1 acts as a neurotrophic factor and supports growth, proliferation and differentiation of the enteric nervous system. Our analysis implicates fetal branch 2b neurons, and ICC/SMC interactions, in the pathogenesis of Hirschsprung's disease.

The gut has its own elaborate nervous system as detailed above, and is also a unique organ with respect to having its own immune system. Secondary lymphoid structures (SLO) are key sites of immune surveillance and serve to initiate adaptive immune response to exogenous pathogens. Development of SLOs in the gut, including gut-associated lymphoid tissues (GALT) and mesenteric lymph nodes (mLN), has been reported to begin between 8 to 12 PCW in humans 40 . We observed mLNs starting around 12 PCW, with structures clearly recognisable and disectable at 15 PCW ( Supplementary Fig. 5A ). In order to better understand the fetal gut-resident immune cells and their organisation in humans, we investigate in detail the key cell types involved in lymphoid tissue organisation during embryonic development.

Lymphoid tissue inducer (LTi) cells, categorised as group 3 innate lymphoid cells (ILC3), interact with lymphoid tissue organiser (LTo) cells to initiate SLO formation during development in mice and humans 41 . Sub-clustering of fetal and adult T and innate lymphoid cells revealed three cell clusters matching published characteristics of LTi cells ( Figure  4A ,B). These characteristics include high expression of RORC, KIT, NRP1, TNF, LTA/B, IL7R, and ITGB7 (integrin α 4β7) ( Figure 4A ,C) as well as the absence of productive α β TCR chains ( Supplementary Fig. 5B ). The cluster labelled innate lymphoid cell progenitor (ILCP) was further defined by high expression of chemokine receptors CXCR5 and CCR7, and high expression of cell adhesion molecule SCN1B 42 and serine protease encoded by HPN gene ( Figure 4C ). The second cluster labelled 'LTi-like ILC3' clustered closely with adult ILC3 cells and had the highest expression of TNFRSF11A (RANK) and its ligand TNFSF11 (RANKL), as well as NCR1 (NKp46) and NCR2 (NKp44) supporting an ILC3 phenotype for this cell population ( Figure 4C ). Both ILCP and LTi-like ILC3 cells expressed high levels of NRP1, previously shown to mark human ILC3 cells with LTi cell function 43 .

Lastly, a third population, which we term 'LTi-like NK/T', express IL17A, ITGAE and CCR9, and NK-associated genes including NKG7, PRF1, and GZMA, but lack expression of CXCR5. While we did not capture productive TCRαβ chains in this population, we observed expression of CD3 genes (CD3G, CD3D and CD3E) and TCRδ genes (TRDC, TRDV2), suggesting they are ү δ T cells ( Figure 4C , Supplementary Fig. 5B ). Reassuringly, these three LTi-like cell subtypes were also identifiable in our full-length scRNAseq data from fetal ileum cells ( Supplementary Fig. 5D ). ILCP cells were enriched at 6-10 PCW, while LTi-like ILC3 and LTi-like NK/T cells were present at 6-17 PCW ( Figure 4C ) and expanded in all gut regions throughout development ( Supplementary Fig. 5E ). Only ILCP cells were captured in fetal mLN ( Supplementary Fig.  5E ), suggesting that LTi-like subsets are expanded during GALT, but not mLN development. This observation together with expression of key chemokine signals (CXCR5, CCR7, CCR6) and RNA velocity analysis suggested that ILCP cells are the first LTi-like cells in the developing gut. Hence they likely represent a progenitor state to mature ILC3s ( Figure 4C , Supplementary Fig. 5F ). RNA smFISH staining placed CXCR5+ RORC+ ILCP/LTi-like ILC3 cells in proximity to CXCL13+ LTo cells in fetal proximal gut mucosa, supporting the concept of congregation of these cells beneath the epithelium in the developing gut ( Figure 4D ).

In order to infer whether the developing lymphoid structures were partitioned into B and T cells zones as evident in mature lymphoid tissue, we characterised B cells in fetal and adult tissues. During fetal development, common lymphoid precursor (CLP) cells (CD34, In contrast to adult B cells, we found no evidence of fetal B cell expansion or maturation in response to antigen. This suggests that additional environmental factors, such as microbiota, might be required for the maturation and zonation of B cells in secondary lymphoid organs. Taken together, our data identifies three types of LTi cells orchestrating the earliest development of lymphoid structures, and B cells that are poised for response in these sites in the prenatal intestine.

formation, and can be either mesenchymal (mLTo) or endothelial (eLTo) in origin. LTi and LTo cell interaction leads to activation of canonical and non-canonical NF-κB pathways, and results in LTo secretion of chemokines CCL19, CCL21 and CXCL13 as well as expression of adhesion molecules VCAM-1, ICAM-1 and MAdCAM-1. After visualising RORC+ CXCR5+ LTi-like ILCP/ILC3 cells in proximity to CXCL13+ cells in fetal gut mucosa as mentioned above ( Figure 4D ), we postulate that GALT formation follows similar mechanisms to those previously described in the mLN. To address this, we analyse stromal populations from both gut tissue and mLN.

Within the entire dataset, we identify multiple stromal populations defined by the expression of VIM, COL1A2 and COL6A1 ( Supplementary Fig. 2B ) including myofibroblasts, smooth muscle cells (SMC, SMC1), pericytes, ICCs, mesothelium, as well as populations resembling stromal cells previously described in the colon (labelled as S1-S4) 8 ( Figure 5A ). We further identified fibroblast populations typically defined in mouse lymph nodes 44 , including T reticular cells (expressing CCL21, CCL19 and GREM1) and follicular dendritic cells (FDCs; expressing CXCL13, CR1 and CR2) ( Supplementary Fig. 7A ). T reticular cells and FDCs were found in the terminal ileum ( Supplementary Fig. 7B ). These samples were also enriched for LZ and DZ GC B cells ( Supplementary Fig. 2B , 6F-G), suggesting the presence of cells that organise Peyer's patches.

In prenatal intestine (12-17 PCW) and mLN, we observe a stromal cell population transcriptionally similar to T reticular, FDCs and S4 fibroblasts ( Figure 5A subpanel, 5B & Supp. Figure 7C ). This stromal cell cluster, absent from pediatric and adult samples, was marked by expression of chemokines CCL19, CCL21 and CXCL13 as well as adhesion molecules ( Figure 5C ), and resembles activated mLTo cells described in developing mouse lymph nodes 45 . Other genes expressed in mLTo cells were components of the NF-κB pathway (NFKBIA, NFKB1/2, IL32, REL, RELA/B), S1P transporter spinster homolog 2 SPNS2, involved in T and B cell trafficking 46 and chemokines CXCL2 and CXCL3 ( Supplementary Fig. 7D ).

We further investigate the transcriptional differences of mLTo cells from lymph nodes and intestinal regions to define differences that could lead to differential positioning of SLOs. Interestingly, large intestinal mLTo cells share transcriptional characteristics with mLN mLTo cells, including expression of C7 and SLC22A3 ( Supplementary Fig. 7E ), suggesting similarities between mLN and cryptopatch development. Amongst the genes highly expressed in mLN mLTo cells were angiotensin receptor APLNR, endothelin receptor EDNRA and extracellular matrix protein TNC, all with roles in blood vessel formation. In addition, mLTo cells in mLNs expressed IL33 and FDCSP, also present in lymph node and spleen fibroblastic reticular cells 47 and observed in fibroblasts expanded in ulcerative colitis 8 . Thus, transcriptionally similar mLTo cells are implicated in the formation of human GALT as well as mLN tissue development.

As lymphatic endothelial cells (LECs) connect the mLN to the tissue and can contribute to SLO formation by interacting with LTi and mLTo cells, we investigated endothelial cell heterogeneity. Amongst endothelial cells we observe arterial, venous, capillary and lymphatic endothelial cells ( Supplementary Fig. 2B , 7F-G). LECs separated into six clusters (labelled LEC1-6), all uniformly expressing PROX1 and LYVE1, and TNFRSF11A (RANK), IL7 and CCL21-signals required for initial recruitment of LTi cells ( Figure 5D ; Supplementary  Fig. 7H ,I). LEC1 and LEC5 ressemble subsets described in human lymph nodes 48 . The LEC1 subset is marked by specific expression of C7, TMEM100 and abundant expression of ACKR4 and CAV1, resembling cells lining the subcapsular sinus ceiling in the human lymph nodes 48 , while LEC5 express CLDN11, ANGPT2 and GJA4, shown to be expressed on collecting lymphatic valves (Supp. Figure 7I ) 48 . We further define LEC3 present mostly in the pediatric and adult intestinal regions, and LEC4 specific to developing intestinal regions ( Supplementary Fig. 7J ), which may represent developing and differentiated lymphatic vessels present in the intestinal mucosa.

LEC2 cells emerge during the second trimester (15) (16) (17) and are present in gut regions as well as fetal mLNs (Supp. Figure 7J ). This population expresses TNFRSF9, THY1, CXCL5, and CCL20, similarly to a population of LECs lining the subcapsular sinus floor in human lymph nodes (Supp. Figure 7I) 48 . In addition, they express genes ANO9, PSD3 as well as high levels of adhesion molecules including MADCAM1, VCAM1 and SELE, suggesting their involvement in lymphocyte recruitment ( Figure 5D, Supplementary Fig. 7K ). LEC6 was marked by expression of ADAMTS4, metallothionein genes such as MT1X, MT1E, and is enriched in the adult appendix and mLN ( Supplementary Fig. 7I ,J). Similarly to mLTo cells, LEC2 and LEC6 express targets of the canonical NFkB pathway ( Supplementary Fig. 7L ). As LEC2 is present during development, NFkB pathway activation suggests that the LEC2 population contributes to organisation of GALTs. Finally, histological staining of fetal ileum allows us to identify blood vessels resembling high endothelial venules ( Figure 5E ), which would permit entry of lymphocytes into the tissue.

Mesenchymal and endothelial cells are central to the recruitment and positioning of leukocytes in developing and mature mLN and GALT structures. To determine cell-cell interactions governing early lymphoid structure formation and leukocyte recruiting, we investigate complementary ligand-receptor expression across LECs, mLTo and LTi-like subsets. Amongst ligand-receptors expressed were LTA and LTB (LTα1β2) expressed on all LTi-like cell subsets and LTBR (LTβR) on mLTo and multiple LEC subsets ( Supplementary  Fig. 8A ), characteristic of LTo and LTi cells 41 . Interestingly, we found that LTi-like ILC3 and NK/T cells utilise CD40LG to signal to mLTo and LEC through either CD40 or integrin a5b1 (ITGA5), respectively ( Supplementary Fig. 8A ). The activation of this signalling pathway in endothelial cells has been previously reported to drive angiogenesis 49 , suggesting that LTi-LEC interaction may further enhance local vascularisation. IL2 was expressed by LTi-like ILC3 and NK/T cells and its receptor NFGR was expressed specifically on mLTo cells. This ligand-receptor interaction may facilitate LTi support for mLTo proliferation and survival in GALTs ( Supplementary Fig. 8A ). In addition, cognate receptor expression was evident across immune cell types that would enable their recruitment into the developing SLO ( Supplementary Fig. 8B ). For example, mLTo cells also expressed CCL2 for recruitment of monocytes and dendritic cells via their receptor CCR2 (Supplementary Fig. 8B ).

SLO development is pre-programmed and does not require antigen stimulation. In contrast, ectopic lymphoid structures, organised aggregates of lymphoid cells, develop in response to chronic inflammation 50 . These ectopic structures often involve T and B cell segregation and increased vascularisation. Compared to SLOs, ectopic lymphoid structures are often transient and resolve after antigen clearance. Following on from our observations of lymphoid structure formation during development, we compared similarities and differences between normal organogenesis and ectopic lymphoid structure formation often observed in Crohn's disease (CD) 51 .

We use logistic regression modeling of pediatric CD scRNAseq data 4 to calculate the mean prediction probability of matching fetal and pediatric CD cell types ( Figure 5F ). We observed ILC3 cells in CD that matched fetal LTi-like ILC3 with 60% probability ( Figure 5F ). S4-like, FDCs and T reticular cells in CD matched fetal mLTo ( Figure 5 ), the former two of which were expanded in four out of seven CD donors ( Supplementary Fig. 8C-D) . Curiously, an equivalent stromal population has recently been shown to be expanded in ulcerative colitis 8 . In addition, transcriptionally similar populations are expanded in leukocyte-rich rheumatoid (RA) 52 and cirrhotic liver disease 53 (Supplementary Fig. 8E ).

To determine cell types involved in CD, we calculated the enrichment score for CDassociated GWAS gene expression for each cell type and selected top cell types (FDR 10%; Methods). Adult ILC3, and fetal ILCP and LTi-like ILC3 cells were among the top cell types enriched for expression of CD-associated genes ( Figure 5G ). Finally, we compared cytokine and chemokine expression between fetal and CD counterparts, and showed core genes conserved between these cells ( Figure 5H ). For example, ILC3 and ILCP/LTi-like ILC3 commonly expressed CXCR4, CCL5, CCL20, while CCL4 was specific to adult ILC3. mLTo and CD fibroblasts showed common expression of CCL19, CCL21, CXCL12-14 and CCL2, while CXCL1/2/3/5/6/8 and CCL11 were specific to adult CD cells.

Our identification of equivalent cell types with the same cellular networking in lymphoid aggregate formation in fetal and CD gut suggests a reactivation of the programs for lymphoid formation during disease.

In this study, we present an integrated dataset of over 346,000 single cells of multiple anatomical regions of the human gut throughout development, childhood and adulthood, browsable online at gutcellatlas.org. We discover epithelial and neuronal lineage composition and developmental relationships, and reveal cell-type specific expression of IBD-and Hirschsprung's disease-associated genes, respectively. We identify the three key cell types of mesenchymal, endothelial and innate lymphoid origins that orchestrate lymphoid organ formation in fetal development and show reinitiation of this cellular program during Crohn's disease.

In the epithelial compartment, we are able to subdivide enterocytes in unprecedented detail. BEST4+ enterocytes have previously been shown to function in pH sensing and transport of metals, ions and salts, but have been restricted to the colon 9,14,54 . Here, we show their presence in the small and large intestines, and as early as the first trimester in the small intestines, supporting a consistent requirement for this population along the gut and throughout life. In addition, we screened for expression of the SARS-Cov2 receptor, ACE2, and associated protease, TMPRSS2, in the gut epithelium. We advance on previous observations of their expressed in epithelial cells 55 by showing their expression in enterocytes throughout life, including fetal development. Since amniotic fluid and placental tissue of an infected new mother have tested positive for SARS-Cov2 56 , this could represent a novel route for viral transmission in utero.

Continuing our analysis of the epithelial compartment, we discover that FCGR2A, a receptor activated by the Fc fragment of IgG upstream of PLCG2, and other downstream signaling molecules are expressed by tuft cells. We suggest that this pathway could mediate direct immune sensing by tuft cells. This is surprising since FcG receptors have rarely been observed in non-haematopoietic cells 57 . To our knowledge, FcgrII expression has only been observed in mouse intestinal epithelial cells after immunization with bacterial neurotoxin 23 and by human nasal epithelial cells in response to bacteria-derived ligands 58 . We also note expression of the FCGR2A by our previously described thymic mTEC(II) epithelial population, which expresses low levels of tuft cell transcription factor POU2F3 (developmentcellatlas.ncl.ac.uk). FcG receptors on immune cells mediate response to local IgG via production of inflammatory signals. Although FCGR2A was only expressed by a subset of tuft cells, we expect this would similarly increase during gut inflammation when luminal IgG levels rise. Overall, this data suggests a novel and potentially very impactful immune-sensing role for intestinal tuft cells.

While epithelial cells are responsible for the major absorptive and secretory functions of the gut, the enteric nervous system is necessary for peristalsis. Absence of the ENS is a very serious condition: Hirschsprung's disease is a congenital disorder characterised by the absence of nerve cells. To date, the genes associated with developmental conditions such as Hirschsprung's disease have been mapped in humans using adult cells 59 . We leverage the fetal cells of this dataset to map the expression of Hirschsprung's disease-associated GWAS genes, and show expression of a set of genes including RET and EDNRB in multiple enteric neuronal cell types in early human gut development. We also discover for the first time a signaling circuitry amongst ENS and ICC and muscle cells involving Hirschsprung's disease genes. These observations help explain the difficulty linking HSRC-associated genes with neuronal cell types in mature gut tissue.

Like the ENS, the immune system of the gut is both sophisticated and essential to proper functioning and maintenance of this barrier tissue. Recent single-cell studies have elucidated the diversity and activation of immune cells in the developing human gut 60,61 . Significant species differences exist in location and size of the SLOs 62 , major hubs of immune cell surveillance and activation. However, our understanding of the formation of these structures is mostly derived from animal models. As such, using genetically perturbed mice, a three-cell model was proposed for SLO formation 41 . Here, we identify cell types with characteristics of the three populations: LTi, mLTo, and eLTo (LEC) cells in the developing mucosa and mLN, suggesting a similar process occurs in humans. We provide detailed characterisation of subtypes of the LTi-like cells in both fetal and adult tissues, describing three different subsets of these cells. Based on previous in vitro studies, where RORC-expressing CD56+ CD127+ IL17A+ cells 63 and ILC3 cells 43 have been both shown to activate mesenchymal cells, we propose that all three of the LTi-like subsets defined in this study act in lymphoid tissue initiation.

Previous studies point to immune cell programs that support normal development, however, when overactivated, can lead to inflammation in early life such as necrotising enterocolitis and IBD 61 . Here we describe SLO formation as one such program adopted in chronic inflammation to recruit and retain immune cells at the site of active inflammation. Using direct comparison of LTi subsets, mLTo and eLTo cells in development and CD, we identify conserved chemokine and cytokine signals between these populations, as well as point to CD specific signals expressed by ILC3, stromal and endothelial cells in disease.

Expansion of an inflammatory fibroblast population, termed Stromal 4 cells, has been reported in ulcerative colitis patients 8 and activated fibroblasts have been reported in adult and pediatric CD 4, 13 . Interestingly, we further observe two specialised fibroblast populations in pediatric CD, namely FDCs and T reticular cells. Similar cell types have been described in mouse lymph nodes 44 , but have not been captured in the gut mucosa. Together with the presence of light zone and dark zone B cells in the gut (Supplementary Fig. 2B, 6F-G) , this suggests the presence of stromal populations involved in SLO zonation as observed in Peyer's patches previously [64] [65] [66] [67] .

Direct comparison of these cells in development versus disease allowed us to separate out expression programs driving SLO formation in development or disease. We show that mLTolike fibroblasts in CD have unique expression of cytokines and chemokines (e.g CXCL2-8) associated with increased recruitment of monocytes and neutrophils to the site of inflammation that may be driving ectopic lymphoid formation. In fetal LTi-like cells, we observed expression of IL2 and IL17A. These inflammatory cytokines may also contribute to controlled inflammatory reaction in the developing fetus, aiding establishment of organised lymphoid tissues prior to birth 68 . While others have proposed that small numbers of LTi-like NK cells are present in post-natal tonsils 63 , we do not observe ILCPs or LTi-like NK/T cells in the pediatric or adult gut. Rather, we show transcriptional similarities between ILC3 and fetal LTi-like ILC3 cells. It is possible that ILC3 in CD may act as an initiator cell to activate the mesenchymal cells and recruit other lymphocytes to the site of inflammation.

Overall, our detailed cell annotation of the human intestine across time and space can be mined for cell-type expression of genes associated to any number of intestinal tract diseases. For instance, our identification of cellular programs of prenatal SLO development reactivated during ectopic lymphoid tissue in CD may represent a target for controlling inflammation during this disease, and could also be extended to other inflammatory diseases characterised by ectopic lymphoid structure formation 69 . Finally, this work opens avenues for using principles of prenatal development and differentiation trajectories outlined here for stromal, epithelial, ENS and immune compartments, to engineer human cells and guide cell maturation in vitro, for the purposes of basic research or even screening, and ultimately as therapeutics for regenerative medicine. TF= transcription factor, Prog = progenitor. G) Summary schematic of signalling involved in Hirschsprung's disease as in F. Supplementary Fig. 1 : Data quality control. A) Schematic with tissue processing strategy for second trimester fetal and adult samples. After enzymatic dissociation, either total fraction was loaded onto 10X chip or CD45+/-cell factions were separated using magnetic cell sorting (MACS) and both fractions were loaded on the 10X chip separately. Lymph nodes were processed without enrichment. Second trimester fetal and adult cell samples were processed using 5' v2 10X kits (Methods). B) Pre-processing and quality control of single-cell RNA-seq data generated in this study and described previously 4 . In short, four datasets, namely first trimester fetal, second trimester fetal, pediatric and adult, were preprocessed separately (including quality control and scrublet doublet removal). Firstly, dimension reduction, clustering and annotation by cell lineage was performed on each dataset separately. Each cell lineage was sub-clustered and fine-grained cell type and state annotation was performed based on marker gene expression. The four datasets were then merged together and each lineage was sub-clustered to unify cell type labels where appropriate. UMAP visualisations show combined dataset colored by sample age, enrichment fraction and donor name. Pre-processing of Smart-seq2 sequencing data Cells with greater than 6000 genes and more than 25% mitochondrial reads were excluded, before regression of 'n_counts', 'percent_mito' and "G2M_score". Cells positive for PTPRC expression (logTPM+1 >=0.2) were taken forward for downstream analysis.

Pre-processing of 10x sequencing data 10x Genomics gene expression raw sequencing data were processed using the CellRanger software v3.0.0-v3.0.2 and the 10X human transcriptome GRCh38-3.0.0 as the reference. The 10x Genomics VDJ Ig heavy and light chains were processed using cellranger vdj v.3.1.0 and the reference cellranger-vdj-GRCh38-alts-ensembl-3.1.0 with default settings.

scRNA-seq quality control and processing of 10x sequencing data Pandas (version 0.24.2), NumPy (version 0.25.2), Anndata (version 0.6.19) and ScanPy (version 1.4) python packages were used to pool single-cell counts and for downstream analyses. Single-cell transcript counts for fetal and adult samples were handled separately to control for anticipated differences in cell expression and sample quality. Cells for each dataset were filtered for greater than 500 genes and less than 50% mitochondrial reads, and genes were filtered for expression in greater than 3 cells. Scrublet (version 0.2.1) score cuttoff of 0.25 was applied to assist in doublet exclusion. Additional doublet exclusion was performed throughout downstream processing based on unexpected co-expression of canonical markers such as CD3D (component of the TCR) and EPCAM. Gene expression for each cell was normalised and log transformed. Cell cycle score was calculated using expression of 97 cell cycle genes listed in Tirosh et al. 2015 70 . Cell cycle genes were then removed for initial clustering. Cell cycle score, mitochondrial reads percentage and unique molecular identifiers (UMIs) were regressed prior to scaling the data.

Batch correction of fetal and adult datasets were performed with bbknn (version 1.3.9, neighbours=2-3, metric='euclidean', n_pcs=30-50, batch_key = "fetal_id" or "batch"), and dimensionality reduction and leiden clustering (resolution 0.3-1.5) was carried out and cell lineages were annotated based on marker gene expression for each cluster. Cell lineages were then subclustered and cell cycle genes reintroduced for batch correction and leiden clustering for annotation of cell types and states again based on marker gene expression. The scoring was done using sc.tl.score_genes() function with default parameters to calculate the average expression of selected genes substrated with the average expression of reference genes.

Single-cell count matrices from three organoid growth conditions were combined together using Pandas (version 0.24.2) and NumPy (version 0.25.2) packages. Cells with fewer than 8000 genes and with less than 20% mitochondrial reads were included in the analysis. Genes with expression in fewer than 3 cells were also excluded. For PLCG2 expression comparison we use normalised (sc.pp.normalize_per_cell) and log transformed (sc.pp.log1p) counts. The data was plotted using Seaborn package barplot and swarmplot functions (version 0.11.0).

To identify region-specific subpopulations, we performed compositional analysis between BEST4+ enterocytes from small and large intestine tissue, using an unpublished tool for differential abundance testing on KNN graph neighbourhoods (implemented in the R package miloR https://github.com/MarioniLab/miloR). Briefly, we performed PCA dimensionality reduction and KNN graph embedding on the BEST4+ enterocytes. We define a neighbourhood as the group of cells connected to a sampled cell by an edge in the KNN graph. Cells are sampled for neighbourhood construction using the algorithm proposed previously 71 . For each neighbourhood we then perform hypothesis testing between conditions to identify differentially abundant cell states whilst controlling the FDR across the graph neighbourhoods.

We test for differences in abundance between the cells from small and large intestine tissue in adult samples and fetal samples. To identify markers of SI-specific and LI-specific subpopulations, we performed differential gene expression analysis using a linear model implemented in the package limma 72 , using 1% FDR (implemented in the function findNhoodMarkers of the miloR package).

For neuronal cell trajectory analysis we use scVelo 0.21 package implementation in Scanpy 1.5.1 73 . The data was sub-clustered on fetal neuronal cells and pre-processed using functions for detection of minimum number of counts, filtering and normalisation using scv.pp.filter_and_normalise and followed by scv.pp.moments function. The gene-specific velocities were obtained using scv.tl.velocity(mode='stochastic') and scv.tl.velocity_graph() by fitting a ratio between unspliced and spliced mRNA abundances. The gene specific velocities were visualised using scv.pl.velocity_graph() or scv.pl.velocity_embedding_grid() functions.

In addition, we calculate pseudotime using Monocle 3 0.2.1 (R version 3.6.3) 74-77 In short, we use preprocess_cds(num_dim=100) to normalise and pre-process raw counts, and reduce dimension and cluster cells using functions reduce_dimension(reduction_method = 'UMAP') and cluster_cells(). We then use functions learn_graph() and order_cells( root cell = "GTATTCTGTCATGCCG-1-4918STDY7426908" ) to order cells along the pseudotime. We then overlay the pseudotime onto the UMAP with velocyto arrows using scv.pl.velocity_embedding_grid() with clusters set to monocle3 inferred "pseudotime". To visualise genes that change along the pseudotime we use sc.pl.paga_path() function with pseudotime set to monocle3 pseudotime. This function required calculation of PAGA parameters and dpt_pseudotime with functions as follows: sc.tl.paga(), sc.pl.paga(), sc.tl.draw_graph( init_pos='paga'), sc.tl.dpt() 78 .

Single-cell BCR analyses were performed as described previously 79 . Briefly, poor quality or incomplete VDJ contig sequences were discarded and all IgH sequences for each donor were combined together. IgH sequences were annotated with IgBLAST 80 before isotype reassignment using AssignGenes.py (pRESTO; 81 ). Ambiguous V gene calls were corrected using TIgGER v.03.1 82 before identifying clonally-related sequences with DefineClones.py (ChangeO v.0.4.5; 82 ) using a threshold of 0.2 for nearest neighbour distances. The germline IgH sequence for each clonal family was determined using CreateGermlines.py (ChangeO v.0.4.5) followed by using observedMutations (Shazam v.0.1.11; 82 ) to calculate somatic hypermutation frequencies for individual sequences. Finally, for integration with the singlecell gene expression object, the number of high quality and annotated contigs per Ig chain (IgH, IgK, IgL) was determined for each cell barcode. If multiple unique sequences for a given chain were detected, that cell was annotated as "Multi" and not considered in further analysis. BCR metadata was combined with the scRNA object for downstream analysis and comparison of different B cell populations.

To infer cell-cell communication and screen for ligand-receptors involved we applied CellPhoneDB v2.0 python package 83,84 on the normalised raw counts and fine cell type annotations from the second trimester intestinal samples (12) (13) (14) (15) (16) (17) . We use default parameters and set subsetting to 5000 cells. To identify most relevant interactions, we subset specific interactions based on the ligand/receptor expression in more than 10% of cells within a cluster and where the log2 mean expression of the pair is greater than 0.

Hirschsprung's disease genes were curated from [35] [36] [37] and The Human Phenotype Ontology website (see URL: https://hpo.jax.org/app/, Aganglionic megacolon HP:0002251) and we selected genes with higher than or equal to 0.3 expression in neuronal lineage single cells. We calculated mean expression per cluster and organ and visualised the normalised mean expression (standard_scale="var") using seaborn.clustermap() function (version 0.11.0). Figure 2D Crohn's disease (EFO_0000384) and ulcerative colitis (EFO_0000729) gene names were retrieved from GWAS Catalog of EMBL-EBI (https://www.ebi.ac.uk/gwas/). The list was downloaded on 27 May 2020. We took the intersecting genes from pooled GWAS genes (P value < 1.00E-6) and the top 100 positive and 100 negative differentially expressed marker genes for each epithelial cell type versus all other epithelial cell types. Expression colour scale has a maximum cut off of 1.5 mean normalised counts to increase resolution at lower expression levels. Only cell types present are shown in each respective gut region heatmap.

Cell type enrichment for IBD-GWAS genes IBD GWAS summary statistics of Crohn's disease (CD) and ulcerative colitis (UC) were obtained from the International Inflammatory Bowel Disease Genetics Consortium (IIBDGC) (see URL: https://www.ibdgenetics.org/). The GWAS enrichment analysis of 103 annotated gut cell types for CD and UC was performed using a fGWAS approach 85, 86 , often used for fine-mapping and enrichment analysis of various functional annotations for molecular quantitative trait and GWAS loci. The association statistics (log odds ratios and standard errors) were converted into the approximate Bayes factors using the Wakefield approach 87 . A cis-regulatory region of 1Mb centred at the transcription start site (TSS) was defined for each gene (Ensembl GRCh37 Release 101). The Bayes factors of variants existing in each cis region were weighted and averaged by the prior probability (an exponential function of TSS proximity) estimated from the distance distribution of regulatory interactions 88 . The likelihood of an fGWAS model was given by the averaged Bayes factors across all genomewide genes multiplied by the feature-level prior probability obtained from a linear combination of cell type specific expression and the averaged expression across all cell types as a baseline expression. The enrichment of each cell type was estimated as the maximum likelihood estimator of the effect size for the cell type specific expression. The code of the pairwise hierarchical model (see URL: https://github.com/natsuhiko/PHM) was utilised for the enrichment analysis. The detailed model derivation is demonstrated in the Supplementary Note Section 1.

Cryosectioning, smFISH, and confocal imaging Fetal gut tissue was embedded in optimal cutting temperature medium (OCT) and frozen on an isopentane-dry ice slurry at -60°C, and then cryosectioned onto SuperFrost Plus slides at a thickness of 10 μ m. Prior to staining, tissue sections were post-fixed in 4% paraformaldehyde in PBS for 15 minutes at 4°C, then dehydrated through a series of 50%, 70%, 100%, and 100% ethanol, for 5 minutes each. Staining with the RNAscope Multiplex Fluorescent Reagent Kit v2 Assay (Advanced Cell Diagnostics, Bio-Techne) was automated using a Leica BOND RX, according to the manufacturers' instructions. Following manual pretreatment, automated processing included epitope retrieval by protease digestion with Protease IV for 30 minutes prior to RNAscope probe hybridisation and channel development with Opal 520, Opal 570, and Opal 650 dyes (Akoya Biosciences). Stained sections were imaged with a Perkin Elmer Opera Phenix High-Content Screening System, in confocal mode with 1 μ m z-step size, using a 20× water-immersion objective (NA 0. 16 

Small intestines from C57BL/6 mice were flushed of faecal content with ice-cold PBS, opened longitudinally, cut into 0.5 cm pieces, and washed by vortexing three times with PBS with 10mM HEPES. Tissue was then incubated with an epithelial stripping solution (RPMI-1640 with 2% (v/v) FCS, 10mM HEPES, 1mM DTT, and 5 mM EDTA) at 37°C for two intervals of 20 minutes to remove epithelial cells. The epithelial fraction was subsequently incubated at 37°C for 10 minutes with dispase (0.3 U/mL, Sigma-Aldrich) and passed through a 100µm filter to obtain a single-cell suspension. Cells were blocked for 20 minutes at 4°C with 0.5% (v/v) heat-inactivated mouse serum followed by extracellular staining in PBS at 4°C for 45 minutes with the following antibodies; EpCAM-FITC (1:400, G8. 

An organoid and multi-organ developmental cell atlas reveals multilineage fate specification in the human intestine

High-Resolution mRNA and Secretome Atlas of Human Enteroendocrine Cells

Single-cell messenger RNA sequencing reveals rare intestinal cell types

Single-cell sequencing of developing human gut reveals transcriptional links to childhood Crohn's disease

Publisher Correction: Tracing the temporal-spatial transcriptome landscapes of the human fetal digestive tract using single-cell RNA-sequencing

Distinct microbial and immune niches of the human colon

The enteric nervous system of the human and mouse colon at a single-cell resolution

Structural Remodeling of the Human Colonic Mesenchyme in Inflammatory Bowel Disease

Intra-and Inter-cellular Rewiring of the Human Colon during Ulcerative Colitis

Heterogeneity and clonal relationships of adaptive immune cells in ulcerative colitis revealed by single-cell analyses

Single-cell atlas of colonic CD8 + T cells in ulcerative colitis

Mucosal Profiling of Pediatric-Onset Colitis and IBD Reveals Common Pathogenics and Therapeutic Pathways

Single-Cell Analysis of Crohn's Disease Lesions Identifies a Pathogenic Cellular Module Associated with Resistance to Anti-TNF Therapy

Colonic epithelial cell diversity in health and inflammatory bowel disease

Evidence of intestinal inflammation in patients with cystic fibrosis

Reduced Neurog3 Gene Dosage Shifts Enteroendocrine Progenitor Toward Goblet Cell Lineage in the Mouse Intestine

SARS-CoV-2 productively infects human gut enterocytes

Paneth cell metaplasia in newly diagnosed inflammatory bowel disease in children

The human gastrointestinal tract-specific transcriptome and proteome as defined by RNA sequencing and antibody-based profiling

Genome-wide association study implicates immune activation of multiple integrin genes in inflammatory bowel disease

Monogenic diseases associated with intestinal inflammation: implications for the understanding of inflammatory bowel disease

Severe Autoinflammatory Manifestations and Antibody Deficiency Due to Novel Hypermorphic PLCG2 Mutations

Intraperitoneal Immunization with Cry1Ac Protoxin from Bacillus thuringiensis Provokes Upregulation of Fc-Gamma-II/and Fc-Gamma-III Receptors Associated with IgG in the Intestinal Epithelium of Mice

Distinct ATOH1 and Neurog3 requirements define tuft cells as a new secretory cell type in the intestinal epithelium

Chemical coding and chemosensory properties of cholinergic brush cells in the mouse gastrointestinal and biliary tract

Tuft cells, taste-chemosensory cells, orchestrate parasite type 2 immunity in the gut

The Immune Function of Tuft Cells at Gut Mucosal Surfaces and Beyond

The Taste Receptor TAS1R3 Regulates Small Intestinal Tuft Cell Homeostasis

Altered expression of laminin alpha1 in aganglionic colon of endothelin receptor-B null mouse model of Hirschsprung's disease

Upregulation of the Nr2f1-A830082K12Rik gene pair in murine neural crest cells results in a complex phenotype reminiscent of Waardenburg syndrome type 4

GDF10 is a signal for axonal sprouting and functional recovery after stroke

Identification of a novel type II classical cadherin: Ratcadherin19 is expressed in the cranial ganglia and Schwann cell precursors during development

Developmental regulation and overexpression of the transcription factor AP-2, a potential regulator of the timing of Schwann cell generation

Neuronal Differentiation in Schwann Cell Lineage Underlies Postnatal Neurogenesis in the Enteric Nervous System

Identification of Genes Associated With Hirschsprung Disease, Based on Whole-Genome Sequence Analysis, and Potential Effects on Enteric Nervous System Development

Sporadic Hirschsprung Disease: Mutational Spectrum and Novel Candidate Genes Revealed by Next-generation Sequencing

Development and developmental disorders of the enteric nervous system

Neurotrophin-3 is required for the survival-differentiation of subsets of developing enteric neurons

Expression and function of Neuregulin 1 and its signaling system ERBB2/3 in the enteric nervous system

Development of human lymph nodes and Peyer's patches

Lymph node stromal cells: cartographers of the immune system

Voltage-gated Na+ channel β 1B: a secreted cell adhesion molecule involved in human epilepsy

Neuropilin-1 Is Expressed on Lymphoid Tissue Residing LTilike Group 3 Innate Lymphoid Cells and Associated with Ectopic Lymphoid Aggregates

Single-Cell RNA Sequencing of Lymph Node Stromal Cells Reveals Niche-Associated Heterogeneity

Nestin-Expressing Precursors Give Rise to Both Endothelial as well as Nonendothelial Lymph Node Stromal Cells

The sphingosine-1-phosphate transporter Spns2 expressed on endothelial cells regulates lymphocyte trafficking in mice

Fibroblast derived IL 33 is dispensable for lymph node homeostasis but critical for CD8 T cell responses to acute and chronic viral infection

Single-Cell Survey of Human Lymphatics Unveils Marked Endothelial Cell Heterogeneity and Mechanisms of Homing for Neutrophils

Inhibition of inflammatory lymphangiogenesis by integrin alpha5 blockade

Ectopic Lymphoid Organs and Immune-Mediated Diseases: Molecular Basis for Pharmacological Approaches

Lymphatics, tertiary lymphoid organs and the granulomas of Crohn's disease: an immunohistochemical study

Defining inflammatory cell states in rheumatoid arthritis joint synovial tissues by integrating single-cell transcriptomics and mass cytometry

Resolving the fibrotic niche of human liver cirrhosis at singlecell level

Lineage-specific expression of bestrophin-2 and bestrophin-4 in human intestinal epithelial cells

Relative Abundance of SARS-CoV-2 Entry Genes in the Enterocytes of the Lower Gastrointestinal Tract

Transplacental transmission of SARS-CoV-2 infection

Fc gamma receptor IIb participates in maternal IgG trafficking of human placental endothelial cells

FcγRIII stimulation breaks the tolerance of human nasal epithelial cells to bacteria through cross-talk with TLR4

The Human and Mouse Enteric Nervous System at Single-Cell Resolution

Memory CD4 T cells are generated in the human fetal intestine

Human Fetal TNF-α-Cytokine-Producing CD4 Effector Memory T Cells Promote Intestinal Development and Mediate Inflammation Early in Life

The lymphoid system: a review of species differences

Human fetal lymphoid tissue-inducer cells are interleukin 17-producing precursors to RORC+ CD127+ natural killer-like cells

Germinal Center Centroblasts Transition to a Centrocyte Phenotype According to a Timed Program and Depend on the Dark Zone for Effective Selection

Phenotypic and Morphological Properties of Germinal Center Dark Zone Cxcl12-Expressing Reticular Cells

Mechanosensing by Peyer's patch stroma regulates lymphocyte migration and mucosal antibody responses

Inflammation, a prototype for organogenesis of the lymphopoietic/hematopoietic system

Ectopic lymphoid-like structures in infection, cancer and autoimmunity

Dissecting the multicellular ecosystem of metastatic melanoma by singlecell RNA-seq

Trajectories of cell-cycle progression from fixed cell populations

limma powers differential expression analyses for RNA-sequencing and microarray studies

Generalizing RNA velocity to transient cell states through dynamical modeling

The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells

Reversed graph embedding resolves complex single-cell trajectories

The single-cell transcriptional landscape of mammalian organogenesis

UMAP: Uniform Manifold Approximation and Projection

PAGA: graph abstraction reconciles clustering with trajectory inference through a topology preserving map of single cells

Antibody repertoire and gene expression dynamics of diverse human B cell states during affinity maturation

IgBLAST: an immunoglobulin variable domain sequence analysis tool

pRESTO: a toolkit for processing high-throughput sequencing raw reads of lymphocyte receptor repertoires

Change-O: a toolkit for analyzing large-scale B cell immunoglobulin repertoire sequencing data: Table 1

CellPhoneDB: inferring cell-cell communication from combined expression of multi-subunit ligandreceptor complexes

Single-cell reconstruction of the early maternal-fetal interface in humans

High-resolution mapping of expression-QTLs yields insight into human gene regulation

Joint analysis of functional genomic data and genome-wide association studies of 18 human traits

A Bayesian measure of the probability of false discovery in genetic epidemiology studies

High-resolution genetic mapping of putative causal interactions between regions of open chromatin

We acknowledge the support received from the Wellcome Sanger Cytometry Core Facility, Cellular Genetics Informatics team, Cellular Generation and Phenotyping (CGaP) and Core DNA Pipelines. This . We thank Professor Roger Barker, Xiaoling He, Alexander Ross and Steve Lisgo for access to and handling of fetal tissue; David Fitzpatrick for discussion on developmental intestinal disorders; Tong Li and Ola Tarkowska for image processing and infrastructure support; We acknowledge J. Eliasova for the graphical images. We thank the tissue donors and donor families.

Human fetal gut samples were obtained from the Human Developmental Biology Resource (HDBR, www.hdbr.org). The maternal consent was obtained through Newcastle hospital (the consent REC reference 18/NE/0290-IRAS project ID: 250012). Pediatric patient material used in intestinal organoid culture was obtained with informed consent as part of the ethically approved research study (REC-96/085).Human adult tissue was obtained by the Cambridge Biorepository of Translational Medicine from deceased transplant organ donors after ethical approval (reference 15/EE/0152, East of England-Cambridge South Research Ethics Committee) and informed consent from the donor families. Fresh mucosal intestinal tissue and lymph nodes from the intestinal mesentery were excised within 1 hour of circulatory arrest; intestinal tissue was preserved in University of Wisconsin organ-preservation solution (Belzer UW Cold Storage Solution; Bridge to Life) and mLN were stored in saline at 4 °C until processing. Tissue dissociation was conducted within 2 hours of tissue retrieval.Mouse samples C57BL/6 mice were obtained from Jackson Laboratories (Margate, UK) and maintained in specific pathogen-free conditions at a Home Office-approved facility at the University of Cambridge. Female mice ages 10 to 14 weeks were used. Mice were housed in accordance with the United Kingdom Animals (Scientific Procedures) Act of 1986.

The fetal gut mesentery was cut to lengthen out the tissue and dissected into proximal ileum (PI), middle Ileum (PI), terminal Ileum (TI), colon and appendix. Samples were washed twice with Hanks Buffered Saline Solution (HBSS; Sigma-Aldrich; cat: 55021C) and minced into pieces using a scalpel. The samples were incubated in 2 ml HBSS solution containing 0.21 mg/ml Liberase TL (Roche; cat: 5401020001) or DH (Roche; cat: 5401089001) and 70 U/ml hyaluronidase (Merck; cat: 385931-25KU) for 50 minutes at 37°C shaking every 5 minutes and homogenised every 15 minutes using a pipette. The single cells were passed through a 40-100 μ m sieve and spun down at 400 g at 4°C for 10 minutes. Red cell lysis solution (eBioscience™ 10X RBC Lysis Buffer (Multi-species) was used according to manufacturer's guidelines to remove red blood cells, and the remaining cells were collected in FACS buffer (1% (v/v) FBS in PBS) by centrifugation at 400 g at 4°C for 5 minutes. All gut region samples (except mLN) proceed to MACS enrichment. 

Intestinal organoids from pediatric patients were cultured in Matrigel® (Corning) using media described before. During organoid culture, the media was replaced every 48-72 hours. Organoids were passaged using mechanical disruption with a P1000 pipette and re-seeded in fresh growth-factor reduced Matrigel® (Corning). When comparing culture media, multiple wells were seeded from a single dissociated sample. The organoids were then allowed to grow for 5 days followed by 24 hour treatment with recombinant human protein TNF-α (H8916, Sigma Aldrich) at 40ng/ml or IFN-gamma (PHC4031, Life Technologies Ltd) at 20ng/ml. Processing for single-cell sequencing analysis was performed by removing the organoids from matrigel at passage 3-4 using incubation with Cell Recovery Solution at 4oC for 20 minutes, pelleting the cells, and re-suspending in TrypLE enzyme solution (Thermo Fisher) for incubation at 37 o C for 10 mins. Cells were pelleted again and re-suspended in DMEM/F12. 

Plate-based scRNA-seq was performed with the NEBNext Single Cell/Low Input RNA Library Prep Kit for Illumina (catalog no. E6420L; New England Biolabs). Total cell fractions from dissociated gut sections of donors BRC2033-2034, F67, F72, F78 were snap-frozen in 10% (v/v) DMSO in 90% (v/v) BSA. Cells were thawed rapidly in a 37 °C water bath and diluted slowly with a pre-warmed FACS buffer (2% (v/v) FBS in D-PBS). Cells were pelleted by centrifugation for 5 min at 300 g, washed with 300 µl of D-PBS and pelleted as before. Cells were resuspended in 100 μ l of Zombie Aqua Fixable Viability Kit (1:200 dilution; catalog no. 423101) and incubated at room temperature for 15 minutes. Cells were washed with 2 ml of FACS buffer followed by 300 μ l of FACS buffer and resuspended in a total of 100 μ l of Brilliant Violet 650 mouse anti-human CD45 (Biolegend catalog no. 304043), Alexa Fluor 700 mouse anti-human CD4 (Biolegend, catalog no. 300526), and APC-H7 mouse anti-human CD19 (BD biosciences, catalog no. 560727) and incubated for 20 minutes in the dark at room temperature. Cells were washed twice with 300 µl of FACS buffer. Single, live, CD45+ cells were FACS sorted into wells of a 384-well plate (catalog no. 0030128508; Eppendorf) containing 2 µl of 1× NEB Next Cell Lysis Buffer (New England Biolabs). FACS sorting was performed with a BD Influx sorter (BD Biosciences) with the indexing setting enabled. Plates were sealed and spun at 100 g

Processed single-cell RNA sequencing objects will be available for online visualisation and download at gutcellatlas.org. The code generated during this study will be available at Github https://github.com/teichlab/.