The transcriptomes, connections and development of submucosal neuron classes in the mouse small intestine

Li, Wei; Morarach, Khomgrit; Liu, Ziwei; Banerjee, Sanghita; Chen, Yanan; Harb, Ashley L.; Kosareff, Joel M.; Hall, Charles R.; López-Redondo, Fernando; Jalalvand, Elham; Mohamed, Suad H.; Mikhailova, Anastassia; Linden, David R.; Marklund, Ulrika

doi:10.1038/s41593-025-01962-x

Download PDF

Article
Open access
Published: 29 May 2025

The transcriptomes, connections and development of submucosal neuron classes in the mouse small intestine

Nature Neuroscience volume 28, pages 1146–1159 (2025)Cite this article

7620 Accesses
1 Citations
80 Altmetric
Metrics details

Subjects

Abstract

The enteric submucosal plexus regulates essential digestive functions, yet its neuronal composition remains incompletely understood. We identified two putative secretomotor neuron classes and a previously unrecognized submucosal intrinsic primary afferent neuron class through single-cell RNA sequencing in the mouse small intestine. Using viral-mediated labeling of each class, we uncovered their morphologies and neural projections in the submucosa–mucosa context, finding connections among all classes and an unexpected close association with enterochromaffin cells. Further transcriptome analysis at the postnatal stage and lineage tracing revealed that neuron identities in the submucosal plexus emerge through an initial binary fate split at neurogenesis, followed by phenotypic diversification, akin to the developmental process of the myenteric plexus. We propose a unified developmental framework for neuronal diversification across the gut wall. Our study offers comprehensive molecular, developmental and morphological insights into submucosal neurons, opening new avenues for exploring physiological functions, circuit dynamics and formation of the submucosal plexus.

Diversification of molecularly defined myenteric neuron classes revealed by single-cell RNA sequencing

Article 07 December 2020

Organization of the human intestine at single-cell resolution

Article Open access 19 July 2023

Modeling gut neuro-epithelial connections in a novel microfluidic device

Article Open access 14 November 2023

Main

The enteric nervous system (ENS) comprises the largest component of the peripheral nervous system and stands out for its ability to perform autonomous functions independently of the brain^1,2. In most mammals, the intestinal ENS consists of the following two major concentric layers: the myenteric (Auerbach’s) plexus, located between the circular and longitudinal muscle layers, which drives peristalsis; and the submucosal plexus, near the mucosa, which is thought to mainly regulate luminal secretion and local blood flow^3,4. Although myenteric neurons are better characterized, submucosal neurons were the first component of the ENS to ever be identified (by Georg Meissner in 1857 (ref. ⁵)). In guinea pigs, studies correlating neurotransmitter repertoire with physiological functions have categorized submucosal neuron types. These investigations suggested the presence of cholinergic (ACh⁺) and noncholinergic secretomotor neurons, both potentially contributing to vasodilation^4,6. Based on major molecular markers, equivalents of secreto/vasomotor neurons have also been identified in the submucosa of mice⁷. However, while intrinsic primary afferent neurons (IPANs; sensory neurons) are present in the submucosa of guinea pigs⁸, immunohistochemical analysis in the mouse intestine suggested a lack of submucosal IPANs⁷. Likewise, combined electrophysiological and neurochemical analysis of mouse submucosa neurons failed to identify IPANs⁹. Nevertheless, advillin, a putative IPAN marker, is expressed in small intestine mucosa-projecting submucosal neurons¹⁰, and neurons with IPAN morphologies have been reported in the colonic submucosa¹¹. To advance our understanding of submucosal neuron functions, development and roles in disease, the genetically tractable mouse model system is desirable to implement. Yet, full leverage of the model system requires a comprehensive molecular classification of mouse submucosal neurons.

Both submucosal and myenteric ENS originate from neural crest cells that invade the gut at early developmental time points and then migrate to populate the full extent of the gastrointestinal tract. The progenitor cells initially occupy the outer layer of the gut wall and differentiate into myenteric neurons and glia cells. At midgestational stages, progenitor cells migrate inwards toward the mucosa and differentiate into the submucosal plexus^12,13,14. While myenteric neuron differentiation begins as early as embryonic day (E)10.5, the first neurons of the small intestine submucosa are observed 7 days later. We recently proposed a developmental model delineating a new principle by which myenteric neuron classes diversify, involving two major neurogenic differentiation paths that progress through phenotypic conversions¹⁵. However, the mechanism underlying the specification of neuronal subidentities from submucosal precursor cells remains unknown.

Single-cell RNA sequencing (scRNA-seq) enables unbiased identification and definition of neuronal classes. Myenteric neurons have been extensively characterized by scRNA-seq^15,16 and single-nuclei RNA-seq^17,18, yet isolated submucosal neurons lack specific interrogation. Here we conducted scRNA-seq on submucosal neurons in the small intestine of juvenile mice, identifying three cardinal submucosal (sm) enteric neuron classes (smENCs), including one class of IPANs (smENC1) and two presumed classes of secretomotor/vasodilatory neurons (smENC2 and smENC3). Viral-mediated labeling in transgenic mice specific to each neuronal class revealed their morphological features and intercellular proximities, suggesting extensive neuronal networks in the submucosa–mucosal space and unanticipated closeness to enterochromaffin cells of all three smENCs, but distance to arterioles. Transcriptionally equivalent myenteric neuron classes suggested shared developmental principles between the two ENS plexi. Through scRNA-seq at a postnatal stage combined with lineage tracing from an early neuronal state, we show that, similar to myenteric neurons, submucosal neurons develop through a stepwise process. Our submucosal single-cell transcriptome atlases, cellular connectivity map and delineation of fundamental developmental processes provide a solid foundation and new ideas for further interrogation of submucosal neuron formation, function and involvement in disease.

Results

Three submucosal neuron classes in the small intestine

To achieve an unbiased view of the neuronal diversity within the submucosal plexus in mice, we aimed to conduct scRNA-seq. While we previously have isolated and characterized myenteric neurons by scRNA-seq^15,16, submucosal neurons presented a challenge to isolate due to their small ganglia situated between the muscularis and mucosal layer. Consistent with myenteric neurons, we observed that submucosal neurons (PGP9.5⁺), but not enteric glia (SOX10⁺), showed selective reporter expression in juvenile Baf53b-Cre;R26RtdTom mice (Supplementary Fig. 1a–c). An optimized method for submucosal neuron isolation was developed using flow cytometry capture of fluorescent cells (~0.7% of sorted cells) from mucosa–submucosal small intestine peels of juvenile Baf53b-Cre;R26tdTom mice (Fig. 1a). Initial analysis of the resulting scRNA-seq dataset identified clear neural clusters, as well as smaller non-ENS clusters, including epithelial and mesenchymal cells (Extended Data Fig. 1a,b and Supplementary Table 1). Analysis of Baf53b-Cre;R26RtdTom gut tissue demonstrated TOM⁺/5-HT⁺ enterochromaffin cells (Extended Data Fig. 1c), while the remaining TOM⁺ cells appeared neuronal; thus, the profiles of blood, muscle and lymphatic cells (Extended Data Fig. 1d,e) probably resulted from faulty inclusion during flow cytometry. After excluding non-ENS cells, small clusters with glial or myenteric neuron character still remained (Extended Data Fig. 1f–j), which were removed to focus on 8,341 cells with unambiguous submucosal identities.

**Fig. 1: A molecular definition of neuron classes of the mouse small intestine submucosal plexus.**

We applied different clustering resolutions and determined that the submucosal neurons consist of three cardinal classes, which we denoted smENC1–smENC3 (Fig. 1b and Extended Data Fig. 2a–c). smENC1–smENC3 exhibited markedly different and comprehensive transcriptional profiles (Supplementary Table 2 and Supplementary Fig. 2). Curiously, many of the genes with the most distinct expression have no known role in submucosal neurons, such as the G-protein-coupled receptor Adgrg6 in smENC1, cadherin protein Pcdh7 in smENC2 and acetyl-group-binding transmembrane protein Fibcd1 in smENC3 (Fig. 1c,d). Genes important to conferring neural-specific activity patterns, including membrane trafficking genes, Ca²⁺-regulators and ion channels, displayed unique patterns (Fig. 1e and Supplementary Fig. 2a,b). All types of ion channels were differentially expressed (DE) among smENCs, with smENC3 showing a particularly extensive set of ion channels, including Scn3b, Kcnd2, Cacna2d1 and Kcnj5. All three smENCs expressed Calb2 (calretinin), none expressed Calb1 (calbindin), while only smENC3 expressed Scgn (Supplementary Fig. 2a). The most DE adhesion molecules included Gpr85, Adgrb2 and Slit2 for smENC1, smENC2 and smENC3, respectively (Fig. 1e and Supplementary Fig. 2c). Selectively expressed transcription factors included Tbx2 in smENC1, Bnc2 in smENC2 and Etv1 in smENC3 (Fig. 1e and Supplementary Fig. 2d). A large diversity of ligands and receptors were selectively or combinatorially expressed (Fig. 1f,g and Supplementary Fig. 2e,f). Of note, smENC1 and smENC2 were identified as cholinergic (Chat and Slc18a3) and CGRP⁺ (Calcb), with smENC1 additionally expressing Nmu and Edn1 and smENC2 expressing Sst and Tac1. In contrast, smENC3 was noncholinergic and characterized by Vip, Npy and Th/Dbh. A wide range of receptors were identified, indicating responsiveness of smENC1 to glutamate (Grin1 and Grin3a) and BMP/Activins (Bmpr2 and Acvr1/Acvr2), smENC2 to vasoactive intestinal peptide (Vipr2) and endothelin (Ednrb) and smENC3 to somatostatin (Sstr1) and glutamate (Gria2, Gria3 and Gria4). All smENCs appear responsive to ACh (for example, Chrna3, Chrnb4, Chrm1 and Chrm2) and are capable of purinergic signaling (for example, Adora1, P2rx2 and/or P2ry6). smENC1 and smENC3 may be responsive to noradrenaline (Adra2a). For a more comprehensive overview of gene categories likely relevant to the function of each smENC, please view Supplementary Fig. 2.

Two putative secretomotor/vasodilatory classes and one IPAN

Prior immunohistochemical characterization of mouse submucosal neurons⁷ alongside comparative investigations of guinea pigs recognized that submucosal ganglia mainly consist of cholinergic (ChAT, CGRP, SST and CALR) and noncholinergic (VIP, NPY and CALR) secretomotor neurons. These neurons likely correspond to smENC2 and smENC3, respectively (Fig. 1h). Additionally, cholinergic neurons lacking secretomotor or IPAN characteristics were identified⁷, constituting 10% of the neurons. We surmise that smENC1 corresponds to these nonidentifiable neurons, but in contrast to the previous analysis, we propose that they are submucosal IPANs, as their transcriptional profile showed significant similarities to myenteric IPANs (myENC6; ref. ¹⁵), including Nmu, Adgrg6 and Syt15 (Supplementary Table 2, Fig. 1h and Supplementary Fig. 3).

The distinction between noncholinergic secretomotor and secretomotor/vasodilatory neurons has been proposed to correlate with tyrosine hydroxylase (Th) expression⁷. To investigate whether Th expression marks a distinct subcluster within smENC3, we performed finer resolution clustering, revealing four subclusters each of smENC2 and smENC3 (denoted a–d), while smENC1 remained as one cluster (Fig. 1i, Extended Data Fig. 3a and Supplementary Table 3). Th correlated with other catecholaminergic genes and was most evidently expressed in smENC3b and smENC3d, although not exclusively (Fig. 1j,k). Thus, the catecholaminergic profile did not align with distinct subcluster(s) of smENC3. Further analysis revealed few uniquely expressed subcluster genes, yet several genes exhibited graded expression across the subclusters (for example, Cd9 and Zbtb20; Extended Data Fig. 3a). Paip2b was more clearly enriched in smENC3d (Extended Data Fig. 3a,b), and we confirmed PAIP2B by immunohistochemistry in submucosal neurons primarily located in anterior small intestine (Extended Data Fig. 3c).

smENCs are transcriptionally similar to myenteric ENCs

smENC1 showed clear transcriptional resemblance to myENC6, but gene expression parallels were also evident between smENC2, smENC3 and certain myenteric neuronal classes (Supplementary Fig. 3a–c). Marker genes such as Sst, Vipr2 and Avpr1a were expressed in smENC2 and myENC5, while Vip, Scgn and Prokr1 were expressed in smENC3 and myENC11. To further validate and explore these inferred similarities, the submucosal scRNA-seq dataset was integrated with our published myenteric scRNA-seq dataset¹⁵ (Fig. 2a). Indeed, smENC1 and myENC6 were found to segregate together, myENC5 aligned with smENC2 and myENC11 juxtaposed smENC3. Additionally, label transfer analysis confirmed the similarities between the myenteric–submucosa pairs (Fig. 2b).

**Fig. 2: Submucosal neuron classes show transcriptional equivalence to three myenteric neuron classes.**

While presumed IPANs in the myenteric (myENC6) and submucosal (smENC1) plexi show transcriptional similarities, they reside in distinct gut compartments and may differ in their specific functions, which should be reflected in their gene expression profiles. Screening for differences, we found that Galr1, Npas1, Cdh8 and Nefm predominantly expressed in myENC6 (Fig. 2c,d), while Clic5, Cyp26b1 and Pkp1 associated with smENC1 (Fig. 2c,e). To assess the validity of Nmu as a marker for smENC1, we performed RNAscope with the transcriptionally more robust smENC1/myENC6 marker Adgrg6 (Fig. 2c). We found that nearly all Adgrg6⁺ neurons co-expressed Nmu in both plexi (Fig. 2f,g). Therefore, we concluded that Nmu-Cre mice, previously used to examine myENC6, also serve as a suitable transgenic tool for smENC1. Using Nmu-Cre;R26RtdTom mice, we confirmed the exclusive expression of Galr1 in myENC6 (Fig. 2h,i). Furthermore, neurofilament medium chain (NF-M; Nfem), previously detected in 63% of myENC6 (ref. ¹⁵), was not found in smENC1 (Fig. 2j). The absence of NF-M expression in smENC1 could explain unsuccessful attempts to identify IPANs in the mouse submucosa, as it was considered an IPAN marker⁷. We furthermore confirmed that most smENC1 but only a small subset of myENC6 expressed Pkp1 (Fig. 2k,l). Thus, we identify several markers that could be used to distinguish myenteric and submucosal IPANs, while Adgrg6 and Nmu are reliable markers for small intestine IPANs.

Class morphologies and proportions across the small intestine

We have unveiled three cardinal submucosal clusters, signified by Nmu (smENC1), Sst (smENC2) and Vip (smENC3). To validate the three clusters in vivo and determine their distribution across the small intestine, we performed immunohistochemistry for SOM, VIP and pan-neuronal PGP9.5 in adult Nmu-Cre;R26RtdTom peels. NMU-TOM, SOM and VIP rarely colocalized within the same neuron (0,09% neurons co-expressed >1 marker); thus, they are specific smENC1–smENC3 markers (Fig. 3a–c and Supplementary Fig. 4a–c). Notably, cells expressing either NMU-TOM, VIP or SOM constituted nearly all PGP9.5⁺ neurons (96–98%; Fig. 3c), suggesting the essentially complete classification of submucosal neurons by smENC1–smENC3. Although all regions contained the three neuron classes, the proportion of smENC1 was higher in the duodenum and jejunum (∼40%) than in the ileum (∼10%; Fig. 3b,c). Conversely, a complementary pattern was observed among smENC2 neurons (duodenum/jejunum, ∼18% and ileum, ∼38%) and smENC3 neurons (duodenum/jejunum, ∼38% and ileum, ∼51%).

**Fig. 3: Submucosal neuron class morphologies and distribution across the small intestine.**

We next acquired Sst-Cre and Vip-Cre mouse lines to target smENC2 and smENC3 selectively (see validation in Supplementary Fig. 4d–g), while Nmu-Cre mice were chosen for smENC1 (Fig. 2f,g). To elucidate the cellular morphology of the three smENCs and compare it with their transcriptionally equivalent myENCs, the class-specific Cre mice were administered adeno-associated virus (AAV.PHP.S) carrying Cre-inducible yellow fluorescent protein (DIO-EYFP; Fig. 3d), and the soma size was assessed in the three Cre lines crossed with R26RtdTom. smENC1 neurons displayed a Dogiel type II morphology (smooth cell body), consistent with myENC6 neurons (Fig. 3e,f and Supplementary Fig. 4h)¹⁵. smENC1s were slightly smaller than myENC6s (Fig. 3g), and 40% displayed more than one long process (Fig. 3e,m). smENC2 exhibited a Dogiel type I morphology (lamellar dendrites; Fig. 3h,m) and were larger than myENC5 (Fig. 3i,j and Supplementary Fig. 4i). Similarly, smENC3 displayed Dogiel type I morphology, in agreement with myENC11 (Fig. 3k–m and Supplementary Fig. 4j). A proportion of both smENC2 and smENC3 had two long processes, and some displayed irregularities at the presumed axonal start (Fig. 3h,k and Supplementary Fig. 4i,j), consistent with previous studies of secretomotor neurons in guinea pigs⁶. Among the submucosal neuron classes, the soma of smENC3 (both TH⁺ and TH⁻) was the largest considering both length and width (Fig. 3n).

Taken together, while smENC1–smENC3 are found in the entire small intestine, they distribute unevenly. The transcriptionally similar neuron pairs in the myenteric and submucosal plexus displayed similar morphologies, but differences in their soma size were observed.

Response patterns of smENC1 suggest an IPAN identity

To further confirm that the smENC1 morphology was aligned with an IPAN identity, low-titer AAV-DIO-tdTOM was administered to Nmu-Cre mice, achieving very sparse single-neuron labeling. From a single soma, we observed an extensive axonal network in the submucosal plane and collateral projections within the mucosal region, penetrating multiple villi (Extended Data Fig. 4a and Supplementary Video 1). To assess the ability of smENC1 to respond to mechanical deformation of the epithelia, full-thickness tissue from Nmu-Cre;LSL-GCaMP6f was imaged at a controlled stage set-up for sequential motorized mucosal poking (Extended Data Fig. 4b). Recordings of myenteric (myENC6; n = 42) or submucosal (smENC1; n = 34) plexus within 24 epochs (measurement associated with poking at a defined spot) identified the following three response categories: nonresponding, quiescent responding and spontaneous (Extended Data Fig. 4c). In both plexi, more than half of the neurons were quiescent responders, which were analyzed further (Extended Data Fig. 4d–j) and categorized based on response latency (Extended Data Fig. 4d). An immediate response was measured in more than half of the epochs per cell in both plexi (Extended Data Fig. 4e). Mean latencies, distribution of response latencies and number of responding epochs per cell were comparable between plexi (Extended Data Fig. 4f–i). Although only one spike was recorded in most epochs, considering the full distribution, larger number of spikes were observed per epoch in smENC1s than in myENC6s (P = 0.0387; Extended Data Fig. 4j). Taken together, smENC1 extends processes throughout multiple villi and is, similar to myENC6, responsive to mucosal poking, with an overall similar response pattern, suggesting that IPANs exist in both ENS plexi of the mouse small intestine.

smENCs communicate through homotypic and heterotypic connections

Neuron-class-specific fluorescent labeling in Cre mice allowed intricate morphology analysis. However, as fluorescence is distributed across the neuronal projections, we reasoned that the methodology could be extended to investigate connections between smENCs by assessing axon-to-soma relationships. Nmu-Cre, Sst-Cre and Vip-Cre mice were injected with AAV-DIO-EYFP, followed by staining for complementary smENC1–smENC3 markers. All three neuron types exhibited homotypic connections, that is, nerve processes between neurons belonging to the same class (Fig. 4a–d). Furthermore, a variety of heterotypical connections (Fig. 4d) were indicated by the proximity of bouton-like structures (puncta) from one class near the soma of another class. For instance, puncta from smENC1 were detected close to neurons expressing VIP (TH⁺ and TH⁻; Fig. 4e,f) or SOM (Fig. 4g). Likewise, puncta from smENC2 juxtaposed Nmu⁺ (Fig. 4h) and VIP⁺ neurons (TH^+/−; Fig. 4i,j), while puncta on axons emanating from smENC3 (TH^+/−) appeared to contact Nmu⁺ (Fig. 4k) and SOM⁺ neurons (Fig. 4l,m). Although we found correlation between boutons and the presynaptic protein synapsin 1 (SYN1), not all puncta colocalized clearly and may not represent release sites (Extended Data Fig. 5a). To label presynaptic sites with higher fidelity, we used AAV to deliver Cre-inducible EGFP-fused synaptophysin (Extended Data Fig. 5b), which has previously been used to mark synapses in the central nervous system (CNS)¹⁹. A high correlation was confirmed between Syp-EGFP puncta and SYN1 in the three class-specific Cre lines injected with AAV-DIO-SYP-EGFP (Supplementary Fig. 5). To investigate homotypic connections, we co-injected AAV-DIO-SYP-EGFP and AAV-DIO-tdTOM into Sst-Cre and Vip-Cre mice or AAV-DIO-SYP-EGFP into Nmu-Cre; R26RTom mice. We identified EGFP⁺ puncta in close proximity to TOM⁺ soma in all three cases, again suggesting homotypic neuronal connections (Extended Data Fig. 5c,f,i). Further single injections of AAV-DIO-SYP-EGFP complemented with smENC markers confirmed all heterotypic connections (Extended Data Fig. 5d,e,g,h,j,k). Finally, a nanobody against PSD95 was applied to visualize postsynaptic sites. Colocalization of presynaptic and postsynaptic sites corroborated submucosal neuron–neuron synaptic interactions (Extended Data Fig. 5l–q).

**Fig. 4: Neuronal projections of submucosal classes in the submucosal plane.**

We reasoned that ligand–receptor analysis within the scRNA-seq profiles could help predict potential signaling pathways used in both homotypic and heterotypic connections. Through manual interrogation (Supplementary Fig. 2) and CellChat (v2), a tool for quantitative analysis of cell–cell communication (Extended Data Fig. 6), we identified, for instance, Calcb–Calclr in homotypic smENC1 connections and Sst–Sstr1 between smENC2 and smENC3, while smENC3 to smENC2 connections may use Vip–Vipr2 signaling (Fig. 4n). Additionally, nonsynaptic signaling interactions included Ptger2 (prostaglandin E), Bmp4, Tgfb1, Fgf1, Wnt9a and Gdf10 signaling (Fig. 4n and Extended Data Fig. 6). The cell–cell contacts may form and exhibit unique interaction patterns through synaptic adhesion molecules. CellChat predicted that each neuron pair has the possibility to form unique connections by using combinatorial patterns of neurexins/neuroligins, contactins, Ncam, semaphorins, Eph–Ephrin and Slitrk/Ptpr (Fig. 4o).

In summary, our connectome analysis implied comprehensive interconnectivity in the submucosal plexus, and we propose potential adhesion molecules within these connections along with signaling molecules that could mediate cell–cell communication.

Connectivity of smENCs within the submucosa–mucosal region

IPANs are believed to receive information from epithelial cells, in particular enterochromaffin cells, while secretomotor neurons are thought to give efferent input, mainly to secretory enterocytes of the crypts. To assess the proximity between the smENCs and epithelial cells (Fig. 5a), we visualized the innervation pattern in class-specific Cre mice crossed with R26RtdTom mice. We observed an extensive network of axons within the villi (Fig. 5b,c and Extended Data Fig. 4a) not only from smENC1 (presumed IPANs) but also from smENC2 and smENC3. Similarly, axonal tracts from all three smENCs surrounded the crypts, harboring Paneth cells (marked by wheat germ agglutinin (WGA); Fig. 5b). Within villi, WGA instead labeled goblet cells, which also appeared close to projections from all three smENCs (Fig. 5c). We next focused on the axonal proximity to enterochromaffin cells due to their acclaimed interaction with IPANs. We immunostained enterochromaffin cells for 5-HT and compared them with DCLK1⁺ tuft cells in the class-specific R26RtdTom mice (Fig. 5d,e). Surprisingly, TOM⁺ axons from all three submucosal classes were found in close and comparable proximity to enterochromaffin cells, which often displayed neuropod structures²⁰, suggesting communication with all three smENCs, not just smENC1 (Fig. 5d,f,g and Supplementary Fig. 6a–c,g–l). In contrast, axons from all three smENCs showed larger distances to tuft cells (Fig. 5e,h,i and Supplementary Fig. 6d–f).

**Fig. 5: Neuronal projections of submucosal classes within the mucosa.**

To predict plausible signaling pathways used by epithelial cells and smENCs, we incorporated a public small intestine epithelial scRNA-seq atlas²¹ (Fig. 5j and Extended Data Fig. 7a) and applied CellChat (v2) analysis, complemented by manual interrogation. Supporting interaction between the submucosal neurons and 5-HT⁺ enteroendocrine cells, we detected serotonin receptors in all three smENCs, not just in presumed IPANs (Fig. 5k, Supplementary Fig. 2 and Extended Data Fig. 7b). Moreover, enteroendocrine cells expressed Tac1, indicating tachykinin production, which could potentially be perceived by Tacr1-expressing smENC2. The overall analysis identified Neurturin–Gfra2, Edn–Ednrb, Wnt, BMP and Prostaglandin E (Ptgs1/Ptges2—Ptger3/Ptger4) as potential signaling pathways emanating from enterocytes, tuft, goblet, stem/TA cells and Paneth cells (Fig. 5k, Extended Data Fig. 7b and Supplementary Fig. 7a–e). Reciprocal signaling from smENCs to epithelial subpopulations included Fgf1–Fgfr2/Fgfr3, Vip–Vipr1, Nrg1/Nrg2–Erbb2/Erbb3, Tnr–Sdc1/Sdc4 and Mif–Cd44/Cd74 (Fig. 5l and Extended Data Fig. 7c).

Neural control of vascular tone is thought to necessitate proximity to blood vessels, in particular arterioles. We investigated the class-specific reporter mice at the level of submucosal arterioles, mucosal capillaries and villus capillaries in an unbiased manner. Axons from all three smENCs occasionally aligned with capillaries in the mucosal space and extended throughout the villus edges (Fig. 5n,o and Supplementary Videos 2–4); however, submucosal arterioles were not obviously innervated (Fig. 5m). In contrast, arterioles are densely innervated by TH⁺ sympathetic fibers, which establish a vasoconstrictive tone^7,22,23. Taken together, our proximity analysis revealed extensive communication possibilities of smENC1–smENC3 with various epithelial cell types, all showing closer connections to enterochromaffin cells than anticipated. However, there was a lack of specific anatomical closeness to suggest substantial direct vasodilatory activities for any of the smENCs.

Postnatal scRNA-seq captures submucosal neuron development

We have previously investigated myenteric neuron development by connecting transcriptomes of juvenile and embryonic ENS. Twelve myenteric neuron class identities (myENC1–myENC12) were found to be gradually acquired within two major neuronal trajectories, named branches A and B. These findings were based on scRNA-seq data dominated by myenteric neuron differentiation¹⁵. To instead investigate submucosa neuron differentiation, we performed scRNA-seq of ENS isolated from the small intestines of Wnt1-Cre;R26RtdTom mice at postnatal day (P) 7 (Fig. 6a), when submucosal neuron differentiation predominate¹². We obtained 5,740 healthy ENS cell profiles and performed cluster analysis (Extended Data Fig. 8a,b and Supplementary Table 4). Clusters were annotated to generic differentiation states (glia/progenitors, neuroblasts or neurons) based on their relative expression of Sox10, Ascl1 and Elavl4 (Fig. 6b,c and Extended Data Fig. 8b,c). Glia/progenitor clusters were signified by the previously identified¹⁵ transcription and signaling factors Foxd3, Mef2c, Pdgfa and Metrn (Extended Data Fig. 9a). Three glia/progenitor clusters corresponded to the cell cycle phases (Fig. 6d and Extended Data Fig. 8c), while one cluster was identified as Schwann cell precursor (SCP) based on previously defined SCP markers in ENS¹⁵, including Col14a1, Gfra3 and Dhh (Fig. 6e and Extended Data Fig. 8c). The remaining glia/progenitor clusters were categorized as presumed glia and expressed generic (Plp1) and subtype glial markers (Gfap; Fig. 6f). Neuroblasts expressed previously identified genes, including Mfap4 and Bcl11b¹⁵, while transcriptional regulators Actl6b, Phox2a and Zcchc12 marked neurons (Extended Data Fig. 9b,c).

**Fig. 6: scRNA sequencing at P7 identifies submucosal development within branches A and B.**

The overall manifold projection mirrored that of earlier stages¹⁵, with no additional major neurogenic root. Neuroblasts appeared to undertake a binary decision, activating transcriptional programs equivalent to the previously described branch A (Etv1) or branch B (Casz1; Fig. 6g,h and Extended Data Fig. 9d). These branch markers were also detected at the juvenile stage (Fig. 6g,h), indicating that smENC3 arises within branch A and smENC1 and smENC2 within branch B. To identify the emergence of smENC1–smENC3, label transfer analysis was conducted by mapping the combined cardinal submucosal and myenteric transcriptome cluster profiles on the P7 cells (Fig. 6i). Cells with predicted identities of smENC1–smENC3 expressed their respective cardinal marker Adgrg6, Sst or Vip (Fig. 6j). Myenteric class markers¹⁵ were also confirmed in cells with the corresponding predicted identities (Extended Data Fig. 9e). Our analysis revealed overall conserved differentiation pathways for myenteric and submucosal neurons.

Identity acquisition via diversification in neuron branches

During myenteric neuron differentiation, we previously discovered that maturing cells transition through ‘prototypic’ continuous states at the root of each neurogenic branch. These states resemble the transcriptional profiles of myENC8/9 (branch A; pENC8) and myENC1/myENC2 (branch B; pENC1)¹⁵ and subsequently transform into the various neuron class identities. We next investigated whether submucosal ENS differentiation follows a similar stepwise process. In branch A, the pENC8 marker NOS1 was associated with the smENC3 differentiation pathway (Fig. 7a), while in branch B, the marker NDUFA4L2 (pENC1) was linked to a subbranch giving rise to smENC1 and smENC2 (Fig. 7b). To validate submucosal prototypic states, we analyzed NOS1 and NDUFA4L2 protein expression in developing submucosal tissue (Fig. 7c–h). At P2–3, ~15% of submucosal neurons expressed NOS1, but this decreased sharply to <1% at juvenile (P23–25) and adult stages (W8–12; Fig. 7c–e). Similarly, ~20% of the submucosal neurons expressed NDUFA4L2 at P2–3, while expression was essentially lost at juvenile and adult stages (Fig. 7f–h). These results show that submucosal neurons transiently exhibit prototypic characteristics, which are lost as they differentiate into smENC1–smENC3.

**Fig. 7: smENC identities are formed through neuronal phenotypic transitions.**

To more firmly determine if smENC3 arises from pENC8 transient neurons, we lineage-traced cells expressing Nos1 in Nos1-Cre;R26EGFP mice (Fig. 7i). At P2–3, 94% of EGFP⁺ submucosal neurons expressed NOS1, whereas this dropped to only 6% in adult mice, confirming a dramatic reduction in NOS1 expression during submucosal neuron development (Fig. 7j–l). Notably, neurons that no longer expressed NOS1 corresponded to smENC3, as all EGFP⁺ NOS1⁻ cells expressed VIP and over one-third of VIP⁺ neurons expressed EGFP (Fig. 7m–o). These findings suggest that smENC3 originates from a transient branch A prototypic identity that at least partially includes NOS1 expression (Extended Data Fig. 10i). Furthermore, no neurons were found to co-express EGFP and the smENC2 marker SOM (Fig. 7p), confirming that transient NOS1 expression is exclusive to branch A neurons.

In parallel, Nos1-Cre;R26EGFP mice were also used to confirm neuronal conversions during myenteric plexus development. Our previous study in this layer indicated that the generation of myENC12 depends on differentiation of NOS1^+/− pENC8 (ref. ¹⁵). Examining the myenteric plexus, we found that 13% of EGFP⁺ myenteric neurons had lost NOS1 expression by weeks 8–12 (Extended Data Fig. 10a–c) consistent with a previous study²⁴. To investigate whether transient NOS1⁺ cells differentiate further into myENC12, we performed costaining with the class markers CALB/NTNG1 and the subclass marker 5-HT. In total, 45% of NTNG1/CALB⁺ neurons and 45% of 5-HT⁺ neurons were found to express EGFP (Extended Data Fig. 10d–g). As expected, no EGFP⁺ cells co-expressed CALR, a marker for branch B-derived myenteric neurons (that is, myENC1, myENC2, myENC5 and myENC6; Extended Data Fig. 10h). These results demonstrate the emergence of myENC12 from NOS1-expressing prototypic branch A neurons in the myenteric plexus, representing a nitrergic to a cholinergic phenotype switch (Extended Data Fig. 10i).

Collectively, our findings strongly suggest that all neurons in the small intestine arise within two primary developmental trajectories (branches A and B; Fig. 7q). The fundamental principle appears to be that differentiating progenitor cells adopt one of two prototypic states and subsequently diversify into distinct neuron classes. However, compared to the myenteric plexus, the final cellular composition of the submucosal plexus exhibits less neuronal diversity, likely due to its later developmental timeline (Fig. 7q).

Discussion

This study offers an in-depth examination of the submucosal neuron classes in the mouse small intestine with regard to their molecular profiles, morphological attributes and connectivity. Furthermore, we reveal a stepwise differentiation process through which submucosal neuron identities arise, suggesting its applicability to the entire developing ENS across the gut wall.

Our scRNA-seq analysis revealed three cardinal submucosal neuron classes characterized by Nmu/ChAT/Calcb (smENC1), Sst/ChAT/Calcb (smENC2) and Vip/Npy/Th^+/⁻ (smENC3). We provide morphological, transcriptional and functional evidence to support the classification of smENC1 as submucosal IPANs, concluding that the lack of NF-M in submucosal IPANs may have hindered proper identification in previous studies. Interestingly, smENC3 expresses genes for noradrenaline production and uptake, but the lack of Dopa decarboxylase (Ddc) suggests they cannot produce noradrenaline under steady-state conditions. Whether certain stimuli enable full biosynthesis or whether they reflect an evolutionary vestige remains unclear. Nevertheless, Th has been proposed to mark functionally distinct VIP neurons⁷, yet subclustering did not associate Th with a molecular subtype. Overall, subclustering revealed no smENC1 subgroup, while smENC2 and smENC3 split into subclusters that lacked unique marker genes, implying that subcluster heterogeneity stems from regional-, maturation- or activity-level variation.

Capitalizing on class-specific markers, we combined Cre mouse lines with transgenic reporters and viral-mediated labeling to uncover connectivity in the submucosa–mucosal space. A surprising finding was the presence of serotonin receptors in all three smENCs and the equal proximity of their projections to enterochromaffin cells, together implying responsiveness to luminal stimuli mediated by serotonin-secreting enterochromaffin cells. Secretomotor neurons, previously thought to rely on IPANs for sensory-evoked activation, may thus possess intrinsic sensory capabilities. IPANs share major neurochemical coding with cholinergic smENC2 (Calcb and Ach), indicating that the ‘division of labor’ may not be as rigid as previously thought and that submucosal neuron classes could respond to slightly different sensory modalities or with different thresholds to the same stimuli. Conversely, smENC1 may execute cholinergic secretomotor functions through axon reflex mechanisms²⁵, supported by our observed axonal closeness between all smENCs and crypt enterocytes. VIP⁺ smENC3 projections in the crypt area, combined with the receptor Vipr1 in enterocytes, fit the presumed major role of VIP⁺ secretomotor neurons in regulating secretion^4,26,27. Tachykinins and acetylcholine in smENC2 corroborated a secretomotor role also of smENC2 (refs. ^4,26,28). Conversely, none of the smENCs exhibited obvious connections with arterioles, aligning with several studies of the mouse intestine^22,23. Hence, vasodilatory actions^4,29 of mouse smENCs may be indirect, possibly mediated through connective tissue or capillaries.

Our systematic connectivity analysis indicated homotypic communication, potentially amplifying local activity, along with extensive heterotypic connections, likely ensuring coordinated and balanced overall responses to mechanical, chemical and osmotic stimuli. Between the two secretomotor types, Sst–Sstr1 and Vip–Vipr2 may mediate communication, although experimental validation is needed. Indeed, the ENS field would benefit from clarifying communication modes (synaptic versus volume transmission) and contact points (axon–soma versus axon–axon or axon–dendrite). Our combined visualization of presynaptic and postsynaptic proteins supports that submucosal neurons communicate at least in part through axon–soma synaptic contacts. We hope that the refined molecular profiles and connections of submucosal neurons presented here may help unravel the cellular and molecular interplay maintaining balance between absorption and secretion at homeostasis versus pathological conditions. However, for any translational purposes, we recognize the importance of careful comparison between the mouse and human submucosal neurochemical coding.

Clear transcriptional similarities were discovered between pairs of myenteric and submucosal neuron classes (smENC1–myENC6, smENC2–myENC5 and smENC3–myENC11). We propose to denote these pairs ‘transcriptional siblings’—a term that highlights molecular equivalence and similarity between cell types that in the adult body are anatomically and/or functionally distinct. One may wonder how similar neurons may emerge in both plexuses. Seminal studies from the Gershon Lab, later expanded by the Young Lab, established that neurons with different phenotypes are generated in birth waves that peak at different developmental time periods^24,30. Although the neural birth-dating would need reinterpretation in the light of recent discrete molecular identity codes, some insights can be drawn from the original studies. For instance, 5-HT⁺ neurons, which we identified as a myENC12 subtype¹⁵, are the earliest generated neuron type (birth peak at E10–E12)^24,30, while scRNA-seq identified myENC5 and myENC11 as late-arising neurons¹⁵. Therefore, it may be the overlapping neurogenesis window at late embryonic-to-early postnatal stages that allows transcriptionally similar cells to be formed in both plexi. Conversely, early born neuron types, such as myENC12 neurons¹⁵, are limited to the myenteric plexus, as they emerge before submucosal development has begun. We suggest a shared developmental program in the ENS, where neuron identities are shaped by the timing of generation. Deciphering the exact mechanisms and molecular programs involved in temporal cell fate specification will be a large endeavor for future studies of ENS development.

Regardless of developmental stage, ENS stem cells undergoing neurogenesis appear to make a binary choice between branch A or B, proceeding through branch-specific prototypic states followed by final identity acquisition. As proof-of-principle for the further diversification of immature but committed neurons, we demonstrate that both smENC3 and myENC12 arise from transient NOS1⁺ branch A prototypic neurons. The transient NOS1 expression could not be attributed to all smENC3 and myENC12 neurons, which either means that not all prototypic cells express NOS1 or that some neurons differentiate in a direct manner. Given that nitrergic neuron deficiencies are hallmarks of several enteric neuropathies, including achalasia, diabetic gastroparesis and Chagas disease³¹, there is considerable clinical interest in producing stem cell-derived enteric NOS1⁺ neurons in vitro^32,33. However, the transient NOS1 expression in both myenteric and submucosal plexus development should raise awareness of the likely occurrence of nonstable NOS1 expression also in cell culture conditions, and measures may be needed to ensure the longevity of phenotypic profiles in manufactured enteric neurons.

Submucosal neurons are linked to a wide range of conditions, including infectious (for example, cholera), inflammatory (for example, inflammatory bowel disease) and functional (for example, irritable bowel syndrome) disorders^34,35, all of which lack satisfactory treatments. We anticipate that the molecular profiles and the differentiation pathways uncovered in this study could pave the way for submucosal neuron derivation from pluripotent stem cells. Possibly by addressing the relative temporal age of progenitors, the generation of late-born neurons could be facilitated. Combined with human gut organoids, submucosal neurons could advance research on neural mechanisms in fluid control at homeostasis or in models of neurocristopathies and inflammatory/infectious conditions.

Methods

Mice

The generation of the Baf53b-Cre (JAX, 027826), Nmu-Cre (MMRRC, 036643-UCD), Vip-IRES-Cre (JAX, 010908), Sst-IRES-Cre (JAX, 013044), Cck-IRES-Cre (JAX, 12706), Wnt1-Cre³⁶ and Nos1-Cre³⁷ (JAX 017526) mouse strains have been previously described. Strains were crossed with Ai14 (R26RtdTom; JAX, 007908), R26EGFP (JAX, 004077) or Ai148D (LSL-GCaMP6f; JAX, 030328). Wild-type mice used in the study were C57/B6J. Animals were group housed, with food and water ad libitum, under 12-h light/12-h dark cycle conditions, 22 °C ambient temperature and 50% humidity. Animal experiments were approved by the local ethics committee in northern Stockholm (Stockholm Norra Djurförsöksetiska Nämnd, Jordbruksverket; N5264/18, N6626-2019 and N5237-2023). Procedures described in Extended Data Fig. 4 were reviewed and approved by the Mayo Clinic Institutional Animal Care and Use Committee. The Mayo Clinic is American Association for Accreditation of Laboratory Animal Care (AAALAC)-accredited (000717) and has Assurance with the Office of Laboratory Animal Welfare, OLAW (A3291-01). All experiments included both male and female mice.

Tissue preparation for scRNA-seq

Juvenile tissue preparation

A total of eight Baf53b-Cre x R26RtdTom mice (both sexes), postnatal day (P)24, were used for scRNA-seq. The tissue was kept in dissection solution—Tris–HEPES recovery solution containing 76 mM Tris–HCl, 19.5 mM Tris base, 2.5 mM KCl, 1.2 mM NaH₂PO₄, 20 mM HEPES, 30 mM NaHCO₃, 25 mM glucose, 5 mM sodium ascorbate, 2 mM thiourea, 3 mM sodium pyruvate, 0.5 mM CaCl₂ and 10 mM MgSO₄ (pH 7.3–7.4). The dissection solution was equilibrated in 95% O₂, 5% CO₂ for 30 min before the start of the experiments and held on ice at all steps. Three fire-polished Pasteur pipettes with decreasing opening size (from ~70% to 10% of the original opening size) were coated in a 1% BSA solution for 1 h at room temperature. The small intestines were cut into 5 cm pieces, and each segment was flushed inside with ice-cold dissection solution. The mesentery was removed, and the intestinal pieces were opened along the mesenteric border, then transferred to prewarmed epithelial strip buffer, consisting of 5 mM EDTA, 1 mM DTT, 5% fetal calf serum (FCS) in Hanks’ balanced salt solution (HBSS) and incubated at 37 °C for 30 min with regular shaking. The pieces were transferred to cold dissection solution, and epithelial cells were brushed off and transferred to a new cold dissection solution. The tissue pieces were pinned, mucosa side down, onto a Sylgard (Dow Corning) coated dish. The smooth muscle layers, containing the myenteric plexus, were peeled off from the submucosa using watchmaker forceps, Peyer’s patches were removed and the submucosa/mucosa layers were stored in dissection solution on ice until all pieces were peeled. Peels were cut into 2–5 mm² pieces and placed in prewarmed digestion solution (1 mg ml⁻¹ Liberase TH Research Grade (Roche), 0.2 mg ml⁻¹ DNaseI (Worthington Biochem, LK003172) and 12 U ml⁻¹ dispase (Corning, 354235) in dissection solution at 37 °C. The tissue was distributed into tubes with 5 ml of digestion solution in each. After 15 min, the tissue was dissociated with the wide-tip Pasteur pipette by pipetting slowly ten times. After an additional 15 min, the tissue was dissociated with a medium-tipped Pasteur pipette ten times. The digestion was stopped by adding 10% BSA and 5 mM EDTA (final concentration) and doubling the volume of dissection solution. A narrow-tipped Pasteur pipette was used for final slow trituration 20 times on ice. The cell suspension was applied to a MACS SmartStrainer (70 m; Miltenyi Biotec) followed by 2 ml fresh cold dissection solution. Single cells were collected by centrifugation at 150g for 10 min at 4 °C. The pellet was resuspended in cold DPBS with 1% BSA. Tomato⁺ cells were sorted by flow cytometry on a BD Biosciences InFlux (v7) cell sorter and collected in ice-cold PBS containing 0.04% BSA (see Supplementary Fig. 8 for representative plots). We captured ~0.7% tdTOM⁺ cells out of the total singlets. Cell gating was performed on BD FACS Software (v.1.0.0.650).

Postnatal tissue preparation

Four Wnt1-Cre;R26RtdTom P7 mice (both sexes) were used for scRNA-seq experiments. During all stages of the dissociation protocol, the tissue was kept in ice-cold DMEM/F12 medium (Invitrogen), except if mentioned otherwise. Small intestines were dissected and flushed with a blunt needle. The mesentery and Peyer’s patches were removed. The intestines were opened longitudinally, washed with PBS, cut into 2–3 cm pieces and put into predigestion solution (5 mM EDTA and 2% BSA in PBS) for 10 min at 37 °C. Epithelial cells were shaken off by vortex, and tissues were washed in DMEM/F12 twice. The tissues were then merged in dissociation solution (0.75 mg ml⁻¹ Liberase TH Research Grade, 3 U ml⁻¹ dispase and 0.2 mg ml⁻¹ DNAse I in DMEM/F12) at 37 °C for 20 min with shaking every 5 min. The enzyme mixture was replaced with DMEM/F12 medium containing 2% BSA and 5 mM EDTA. The cells were manually triturated using three fire-polished Pasteur pipettes with decreasing opening size (from ~70% to 10% of the original opening size) that were previously coated with 1% BSA solution for 1 h at room temperature. The single-cell suspension was filtered through a DMEM/F12-equilibrated 40 μm cell strainer (Miltenyi Biotec) and centrifuged at 150g for 5 min at 4 °C. Cells were resuspended in DPBS containing 1% BSA and 1% DRAQ7 (Biostatus). tdTom⁺ cells were sorted on a BD Influx (v7) Cell Sorter and collected in ice-cold PBS containing 0.04% BSA. We captured 7.2% tdTom cells out of the total Draq7 singlets. The gating of cells was performed on BD FACS Software (v.1.0.0.650).

scRNA-seq

The sampling was carried out with the 10x Genomics Chromium Single-Cell Kit (10x Genomics) version 2 (P7) and version 3 (P24) at the Eukaryotic Single-Cell Genomics (ESCG). Cell suspensions were adjusted to 500–1,000 cells per μl and added to 10x Chromium reverse transcription mix to achieve loading target numbers between 2,500 and 5,000 cells. The manufacturer’s instructions were followed for downstream cDNA synthesis, library preparation and sequencing (NovaSeq S2 reagent kit for P24 and NovaSeq S1 for P7) using the NovaSeq 6000 platform (Illumina). Cellranger (v.3.0.1) was used.

Virally mediated neuron labeling

ssAAV-PHP.S/2-shortCAG-dlox-EYFP(rev)-dlox-WPRE-hGHp(A) (AAV-DIO-EYFP), ssAAV-PHP.S/2-shortCAG-dlox-tdTomato(rev)-dlox-WPRE-hGHp(A) (AAV-DIO-tdTomato) and ssAAV-PHP.S/2-shortCAG-dlox-rSyp1_EGFP-dlox-WPRE-bGHp(A) (AAV-DIO-SYP_EGFP) were constructed at the Viral Vector Facility of the Neuroscience Center Zurich (v309-PHP.S and v167-PHP.S). Plasmid pAAV-pCAG-DIO-rSyp1_EGFP was a kind gift from S. Arber. Virus was injected intravenously into the tail vein of 8–12 week-old mice. Virus titer and injected dosage are listed in Supplementary Table 5c. Animals were killed 2.5–4 weeks after injection. The procedure was different for obtaining the single-neuron morphology (‘Whole-neuron imaging’).

Tissue preparation for histochemical analysis

Juvenile and adult small intestine sections

The small intestine was dissected from juvenile and adult mice, the mesentery was removed and the duodenum, jejunum and ileum were separated. The tissues were prepared as tubes or Swiss rolls. For tubes, the intestinal pieces were flushed clean with ice-cold PBS and fixed overnight in 4% paraformaldehyde (PFA) in PBS at 4 °C. Tissues were then washed three times with PBS and incubated in 30% sucrose in PBS overnight at 4 °C and embedded in Optimal Cutting Temperature (OCT) Cryomount (Histolab), snap frozen and stored at −80 °C. The tissues were cut at 20 μm or 100 μm. For Swiss rolls, briefly, the tissues were cut open along the longitudinal axis and placed flat, mucosa side down, on the bottom of a Petri dish. A toothpick was placed at one end, and the tissue was rolled up to form a roll. The rolls were pinned to a Sylgard-coated dissection dish and fixed in 4% PFA at 4 °C overnight. Tissues were then washed three times with PBS and incubated in 30% sucrose in PBS overnight at 4 °C and embedded in OCT Cryomount (Histolab), snap frozen and stored at −80 °C. The tissue was cut at 14 μm or 20 μm.

Submucosal plexus peels

Duodenum, jejunum and ileum intestinal segments from juvenile and adult mice were cleaned of mesentery and opened lengthwise along the mesenteric border. The intestinal contents were rinsed out with ice-cold PBS. The tissue was stretched and pinned with the mucosa side down on a Sylgard-coated dish. Longitudinal muscle-myenteric plexus-circular muscle layer was peeled off, and the remaining tissue containing submucosal plexus was fixed in 4% PFA at 4 °C overnight. Tissue was then washed three times with PBS and cut into 1 cm² segments and either used directly for histochemical analysis or kept in PBS at 4 °C.

Postnatal small intestine sections

Small intestines were dissected out from P2–3 pups and fixed in 4% PFA in PBS at 4 °C for 2 h. Samples were subsequently washed in PBS and cryoprotected in 30% sucrose in PBS at 4 °C overnight. The tissue was then dissected into distinct anatomical regions and embedded in OCT Cryomount and kept frozen at −80 °C. Samples were sectioned at 10–12 μm.

Immunofluorescence analysis

Frozen tissue sections were air-dried at room temperature for 1 h, rinsed with PBS and incubated in blocking solution containing donkey anti-mouse Fab (Jackson Laboratories, 715-007-003) diluted 1:50 in PBS for 1 h at room temperature. Sections were then blocked with 2% normal donkey serum (NDS; Jackson Laboratories) and 0.1% Triton X-100 (Sigma) in PBS for 1 h and incubated with primary antibodies diluted in the same solution overnight at 4 °C. They were washed three times with PBS (10 min each) and incubated with secondary antibodies at room temperature for 1 h. Tissue was washed three times with PBS (10 min each) and mounted in DAKO mounting medium (Agilent) containing DAPI (4′,6-diamidino-2-phenylindole).

Submucosal peels and 100 μm floating transverse sections were incubated in blocking solution containing donkey anti-mouse Fab diluted 1:50 in PBS and 0.3% Triton X-100 overnight at 4 °C. The blocking solution was replaced by incubation buffer containing 2% NDS, 1% BSA and 0.3% Triton X-100 in PBS for about 8 h at 4 °C. The tissue was incubated in primary antibodies (and nanobody) diluted in incubation buffer for 48 h at 4 °C, washed three times with PBS (30 min each) and then placed in the secondary antibodies diluted in incubation buffer for 2 h at room temperature. The tissue was washed three times with PBS (15 min each) and mounted on glass slides using DAKO mounting medium containing DAPI. WGA staining was performed after antibody staining using fluorescein-conjugated WGA (5 μg ml⁻¹; Thermo Fisher Scientific, W834) for 30 min at room temperature, followed by PBS washing. For a complete list of primary and secondary antibodies, see Supplementary Table 5a,b.

RNAscope

RNAscope on 20 μm sections was performed according to the manufacturer’s instructions (RNAscope Multiplex Fluorescent Reagent Kit v2 Assay; Bio-Techne). Protease III was used for tissue digestion. The following probes were used for the analysis: Mm-Adgrg6 (467351-C2), Mm-Nmu (446831-C3), Mm-Galr1 (448821-C2) and Mm-Pkp1 (1247231-C1). Immunohistochemistry (as described above) was applied after completion of RNAscope. RNAscope was also conducted on submucosal peel tissue based on protocols from ref. ³⁸. In brief, tissue was prepared as for immunohistochemistry, but villi were scraped off to allow better probe penetration. Experiments were carried out in 48-well plates. Protease IV was used for tissue digestion. Probes were incubated overnight for better hybridization.

Confocal imaging

Imaging of z-stacks with intervals of 1 μm was performed using Zeiss LSM800 or LSM980 confocal microscopes using ×20, ×40 or ×63 objectives with Zeiss ZEN 3.3 software. The images were then processed in Fiji image analysis software (v.2.0.0-rc-69/1.52i; National Institutes of Health (NIH)). Color was, in some cases, changed, and each channel was individually adjusted. Distance analysis (Fig. 5f–i) was performed in Imaris (v10). Briefly, three-dimensional (3D) surfaces were created for enterochromaffin cells, tuft cells and nerve processes individually, and the shortest distance between enterochromaffin cells or tuft cells and nerve processes was calculated.

Whole-neuron imaging

Nmu-Cre mice aged 14 weeks (males and females) were injected intraperitoneally with 100 µl of the final dose 1.5–3 × 10¹¹ vector genomes PHP.S-DIO-tdTomato. Viral doses were determined in preliminary experiments to label <0.5% of enteric neurons and used to achieve sparse labeling. After 18–23 weeks, mice were killed and the ileum was collected, opened along the mesenteric border, pinned flat and fixed in 4% PFA at 4 °C for 24 h. The 5 cm segments of ileum were cleared for at least 10 days in Ce3D as described³⁹ with the addition of 2.5% DABCO-33LV (Sigma) as antifading agent. Tissues were mounted on glass slides and imaged using a Nikon AXR microscope at a resolution of 0.43 × 0.43 × 0.63 µm. Nikon Elements (v.5.41.01) was used to stitch tiles, and the image was imported to the Volume Render subroutine of ANALYZE software (v.12; Mayo Foundation). Using the Trace Tool, starting with the neuron cell soma, a series of clipped volumes were set to subjective binarization thresholds until the annotator considered the nerve fibers as fully isolated objects, and these voxels were selected and manually assigned to a new object. This was repeated for every branched neurite until further fluorescence was not observed. The resulting single binarized object contained the 3D reconstructed neuron. Neurites were subjectively assigned to tissue regions, identified by autofluorescence. One reconstructed submucosal ileal Nmu⁺ neuron was chosen to represent several observations, as its fibers did not overlap with fibers from other labeled Nmu⁺ neurons. The 3D color-coded reconstructed neuron was rendered in different orientations to create the three-plane visualizations. In addition, sequences of renderings were generated using ANALYZE (v.14) at different orientations and zoom levels to generate video sequences.

Epithelial deformation functional experiment

Experiments were performed on eight Nmu-Cre;LSL-GCaMP6f mice (six males and two females) at 18–22 weeks. Segments of the ileum were excised into normal Krebs solution with 10 μM nifedipine and 100 nM atropine to minimize contractions. Tissues were opened along the mesenteric border and pinned flat on a thin silicon elastomer in a glass bottom dish that was mounted on an inverted IX70 microscope (Olympus) equipped with an ES107 motorized stage and PS3H122R focus motor with Optiscan III controller (Prior), X-Cite X-LED1 light source (Excelitas Technologies) with QuadCube and GFP filters (Semrock), ×10 UPlanFL and ×20 and ×40 LUCPlanFLN objectives and ORCA-Flash4.0 CMOS camera (Hamamatsu) with two parallel Camera Link outputs to a Firebird Frame Grabber (Active Silicon). The 76 mm platinum iridium concentric bipolar microelectrode with a 2–3 μm tip (World Precision Instruments) was used to deliver focal mechanical stimuli from an S88 stimulator (Grass) and was positioned in X and Y with two stacked CONEX-SAG-LS16P closed-loop linear piezo stages with 16 mm range at 25 nm resolution (Newport), as well as in Z with a P-625.1CD closed-loop linear piezo stage with 500 μm range at <1 nm resolution (Physik Instrumente). Automated acquisitions (50 fps) controlled these stages to maneuver the electrode to deform 24 points in a 4 mm² area over 19 min. Each epoch started with the lowering of the electrode 25 μm into the mucosa and a further lowering of 25 μm 10 s later. The electrode was raised after an additional 15 s. The baseline period of the subsequent epoch, during which the electrode was raised, lasted 20 s and overlapped with the poststimulus period. Peripherals were controlled through Clampex (v11.2) and MetaMorph (v7.10) with custom protocols. Image series were motion corrected via a template matching plugin for ImageJ Fiji and mean intensity values were extracted from manually curated ~150 μm² regions of interest that were transferred to a custom JupyterLab (v3.5.3; Python 3 kernel) notebook for trace feature extraction, especially ‘spike’ detection, analysis and visualization. Cells in the submucosal or myenteric plexus were subjectively categorized into nonresponding, quiescent responding and spontaneous cells based on ongoing activity and responses during any of the 24 mucosal deformations. Note that all spontaneous cells had at least one epoch with increased firing frequency during mechanical deformation, so they could also be regarded as responsive but were not included in subsequent analysis. Epochs in quiescent responding cells were further analyzed. Each stimulus epoch was analyzed for the number of spikes and the latency to the first spike. Epoch latencies were categorized into immediate (<0.5 s from the start of the mechanical stimulus), delayed (>0.5 s after the mechanical stimulus but during the 25-s period of mechanical deformation) and poststimulus (response within 20 s after the deformation was released). Comparisons were made between submucosal (total 34 neurons) and myenteric (total 42 neurons) Nmu⁺ neurons.

Cell counting, statistical analysis and reproducibility

Counting of fluorescent cells was performed on confocal images or directly under a Zeiss fluorescent microscope. GraphPad Prism (versions 9 and 10) was used to generate bar, dot and box-and-whisker plots. In all graphs, error bars indicate mean ± s.d. of the mean. Cell sizes were measured on z-stack images. The density plots of the shortest distance were generated using RStudio (Fig. 5f–i). No statistical methods were used to predetermine sample sizes, but sample sizes are similar to those reported in previous publications. Data collection and analysis were not performed blind to the conditions of the experiments. Samples with poor immunohistochemical staining quality were excluded from further analysis. Data distribution was assumed to be normal, but this was not formally tested.

Analysis of scRNA-seq data

Clustering analysis

For the juvenile submucosal (P24) dataset, we started with 11,428 cell-containing droplets from CellRanger (3.0.1). Cells with less than 200 unique molecular identifiers (UMIs) and genes that were expressed in less than three cells were excluded. From the frequency distribution, we applied a global thresholding—(1) excluded cells with a fraction of mitochondrial genes more than 0.2, (2) retained cells with a number of genes within the range of 500–6,500 and (3) excluded cells with more than 40,000 UMIs. Given that our data also contained non-ENS cell types, a more fine-tuned filtering was required. We performed clustering across a range of resolutions (0.2, 0.6, 1.0 and 1.8), resulting in plausible clusters, from each of which we removed cells with the total number of UMIs below the 1st percentile and above the 99th percentile. This step helped further eliminate unhealthy cells and homotypic doublets or triplets and was followed by the usage of the DoubletFinder (v2.0.3) algorithm. After these procedures, we retained 9,515 nonambiguous high-quality cells for further analysis. Initially, we normalized the expression data using the regularized negative binomial regression method implemented in SCTransform(). For principal component analysis (PCA), we used 3,000 variable genes and excluded sex-specific genes (Xist, Gm13305, Tsix, Gm8730, Eif253y, Ddx3y, Uty and Kdm5d) from the variable genes to prevent artifacts in the dimensionality reduction. However, these genes were retained for subsequent analyses to ensure their inclusion in the dataset. In the PCA space, we constructed a K-nearest neighbor graph based on Euclidean distance and identified distinct cell communities using the Louvain algorithm (res. = 1.8). From here, we performed iterative clustering and excluded non-ENS cells at each round. After each round, the raw data containing more purified enteric neurons were analyzed again. At level 1, we removed non-ENS clusters (Extended Data Fig. 1a–e). In further iterations, we identified and eliminated clusters with glial identity (Sox10) and suspected contaminating myenteric cells based on their similarity in DE genes to major myenteric plexus clusters, their small fraction compared to other clusters in the dataset and validating histochemical experiments (Extended Data Fig. 1f–j). Finally, only submucosal neurons (8,341) underwent third-level 30 principal component (PC) clustering, where res. = 0.04 resulted in three cardinal clusters and res. = 0.4 resulted in nine clusters. Visualization of the clusters was performed on two dimensions using the RunUMAP() function (min.dist = 0.5, n_neighbors = 30L, umap.method = ‘umap-learn’, metric = ‘correlation’).

Likewise, for the postnatal day (P) 7 datasets, we analyzed 8,179 starting cells. After removing cells with more than 0.1 fraction of mitochondrial genes, with more than 6,500 and less than 500 genes and with more than 40,000 UMI counts, we performed first-level clustering of the remaining 6,062 cells in the same manner as we did for the juvenile dataset. Here we excluded sex-specific genes (Xist, Gm13305, Tsix, Gm8730, Eif253y, Ddx3y, Uty and Kdm5d) and a set of immediate early genes (Fos, Jun, Junb and Egr1) from variable genes before computing PCAs. We then constructed a K-nearest neighbor graph and identified communities of cells using the Louvain algorithm. At this step, we refined the cell population by removing five non-ENS clusters. One cluster corresponded to fibroblast-like cells described previously in embryonic scRNA-seq in ref. ¹⁵. One cluster likely constituted myocytes (Myh11, Cnn1 and Ckm); three other clusters corresponded to hematopoietic (Lcp1), endothelial (Lyve1 and Cldn5) and epithelial (Muc13) cells. The remaining ENS cells (5,740) were then subjected to SCTransform() and second-level clustering (50 PCs, res. = 1.8), yielding 24 clusters. We then performed differential expression analysis based on DE genes, categorizing each cluster as one of the following major cell types: neuron (Elavl4^high), neuroblast (Ascl1^high), dividing progenitors (for example, Cenpa and Top2a), enteric glia (S100b and Plp1) and SCP cells (Dhh, Col14a1 and Gfra3). Clusters considered to share a generic differentiation state were merged (Fig. 6b). The clusters were visualized in two dimensions using the RunUMAP() function (min.dist = 0.4, n_neighbors = 60L, umap.method = ‘umap-learn’, metric = ‘correlation’).

Sex inference, dataset integration and identity mapping

Sex inference (predicting the sex of cells in the scRNA-seq datasets) was performed using classifySex() in the Oshlack/speckle (v0.0.3) package⁴⁰. The classifier had been trained on mouse data. Cells with zero counts on Xist and the sum of the Y chromosome genes were annotated as N/As. To integrate different scRNA-seq datasets we used the raw counts from myenteric neurons¹⁵ and from our submucosal neurons (myENC1–myENC12 and smENC1, smENC2a–smENC2d; smENC3a–smENC3d) for integration analysis. Each dataset was normalized separately using SCTransform() to obtain 3,000 variable genes. Subsequently, we merged the two lists of variable genes from each dataset to obtain a shared set of variable genes, which were then used to compute PCAs (50 PCs) in the Harmony algorithm. Cells in this integrated reduced space were visualized using uniform manifold approximation and projection (UMAP; min.dist = 0.4, n_neighbors = 45L, umap.method = ‘uwot’, metric = ‘cosine’; Fig. 2a). Label transfer was performed using the FindTransferAnchors() function followed by the TransferData() function, with myENC1–myENC12 as the reference dataset. Finally, the mapping identities were presented as Alluvial plot using the R package ggalluvial (Fig. 2b). For label transfer to the P7 dataset, we first built a reference dataset by combining our previously published myenteric neuron dataset¹⁵ with the submucosal neuron dataset, as well as enteric glia¹⁶, representing the most comprehensive mature reference dataset so far. After merging the three datasets from the mature ENS, we performed SCTransform (method = ‘glmGamPoi’). Using 3,000 variable genes, we computed PCAs and built a reference reduced space (100 PCs). We retained annotations for published myenteric neuron classes and submucosal classes (as defined in this study) as labels for transfer to the P7 dataset. For this, we used the FindTransferAnchors() function and then the TransferData() function (with smENC1–smENC3, myENC1–myENC5, myENC10, myENC12 and enteric glia as labels). To avoid competition between transcriptionally similar neuron classes, myENC5, myENC6 and myENC11 were masked from the reference dataset when performing label transfer. P7 cells were assigned according to their maximum predicted score (Fig. 6i).

Ligand–receptor analysis

To infer possible signaling communication between submucosal neuron classes and between these neurons and intestinal epithelial cell types, we performed ligand–receptor analysis with CellChat (v2)^41,42. Epithelial cell clusters were retrieved from the scRNA-seq dataset GSE92332 (ref. ²¹). CellChat identifies overexpressed ligands or receptors in one cell group and then identifies overexpressed ligand–receptor interactions if either ligand or receptor is overexpressed. Based on the extent of overexpression combined with prior knowledge of the interactions between signaling molecules (ligands, receptors and their cofactors), CellChat determines and weights biologically significant cell–cell communication. We then calculated the overexpressed genes and the corresponding interactions using the CellChat functions identifyOverExpressedGenes() and identifyOverExpressedInteractions(), respectively. To identify the significant interactions between submucosal neurons (smENC–smENC), we used 20% truncated means, while 10% truncated means were used to identify smENC–epithelial cell communications. A ten percent threshold was used to detect ligands expressed in epithelial cell clusters that were small and heterogeneous (for example, enterochromaffin cells). Circle plot was generated using netVisual_circle() to show the total number of interactions among all cell groups. The communication probability of interactions was calculated using computeCommunProb(). Significant interactions (L–R pairs) from senders to the receiver cell groups were plotted using netVisual_bubble().

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Data availability

Raw sequence and processed data for submucosal small intestine Baf53b-Cre;R26RtdTom (P24) and small intestine Wnt1-Cre;R26RtdTom (P7) are available on the Gene Expression Omnibus database under the identifier GSE263422. We have additionally used mouse small intestine epithelial cell scRNA-seq from GSE92334 (ref. ²¹), myenteric mouse small intestine neuron scRNA-seq from GSE149524 (ref. ¹⁵) and glial cells from mouse small intestine myenteric scRNA-seq from SRP135960 (ref. ¹⁶). Source data are provided with this paper.

Code availability

The data analysis code for P24 and P7 is available at https://github.com/UMarklundLab/SubmucosaENS. Code used for calcium imaging analysis is available at https://github.com/Linden-Laboratory/Mucosal-Deformation.

References

Sharkey, K. A. & Mawe, G. M. The enteric nervous system. Physiol. Rev. 103, 1487–1564 (2023).
Article CAS PubMed Google Scholar
Furness, J. B. The enteric nervous system and neurogastroenterology. Nat. Rev. Gastroenterol. Hepatol. 9, 286–294 (2012).
Article CAS PubMed Google Scholar
Furness, J. B. The Enteric Nervous System, Second Edition (Blackwell Publishing, 2006).
Vanner, S. & Macnaughton, W. K. Submucosal secretomotor and vasodilator reflexes. Neurogastroenterol. Motil. 16, 39–43 (2004).
Article PubMed Google Scholar
Neckel, P. H. Annotated translation of Georg Meissner’s first description of the submucosal plexus. Neurogastroenterol. Motil. 35, e14480 (2023).
Article PubMed Google Scholar
Furness, J. B. et al. Morphologies and projections of defined classes of neurons in the submucosa of the guinea-pig small intestine. Anat. Rec. A Discov. Mol. Cell Evol. Biol. 272, 475–483 (2003).
Article PubMed Google Scholar
Mongardi Fantaguzzi, C. et al. Identification of neuron types in the submucosal ganglia of the mouse ileum. Cell Tissue Res. 336, 179–189 (2009).
Article CAS PubMed Google Scholar
Pan, H. & Gershon, M. D. Activation of intrinsic afferent pathways in submucosal ganglia of the guinea pig small intestine. J. Neurosci. 20, 3295–3309 (2000).
Article CAS PubMed PubMed Central Google Scholar
Foong, J. P. et al. Properties of cholinergic and non-cholinergic submucosal neurons along the mouse colon. J. Physiol. 592, 777–793 (2014).
Article CAS PubMed PubMed Central Google Scholar
Melo, C. G. S. et al. Identification of intrinsic primary afferent neurons in mouse jejunum. Neurogastroenterol. Motil. 32, e13989 (2020).
Article CAS PubMed Google Scholar
Wong, V., Blennerhassett, M. & Vanner, S. Electrophysiological and morphological properties of submucosal neurons in the mouse distal colon. Neurogastroenterol. Motil. 20, 725–734 (2008).
Article CAS PubMed Google Scholar
McKeown, S. J., Chow, C. W. & Young, H. M. Development of the submucous plexus in the large intestine of the mouse. Cell Tissue Res. 303, 301–305 (2001).
Article CAS PubMed Google Scholar
Uesaka, T., Nagashimada, M. & Enomoto, H. GDNF signaling levels control migration and neuronal differentiation of enteric ganglion precursors. J. Neurosci. 33, 16372–16382 (2013).
Article CAS PubMed PubMed Central Google Scholar
Jiang, Y., Liu, M. T. & Gershon, M. D. Netrins and DCC in the guidance of migrating neural crest-derived cells in the developing bowel and pancreas. Dev. Biol. 258, 364–384 (2003).
Article CAS PubMed Google Scholar
Morarach, K. et al. Diversification of molecularly defined myenteric neuron classes revealed by single-cell RNA sequencing. Nat. Neurosci. 24, 34–46 (2021).
Article CAS PubMed Google Scholar
Zeisel, A. et al. Molecular architecture of the mouse nervous system. Cell 174, 999–1014 (2018).
Article CAS PubMed PubMed Central Google Scholar
May-Zhang, A. A. et al. Combinatorial transcriptional profiling of mouse and human enteric neurons identifies shared and disparate subtypes in situ. Gastroenterology 160, 755–770 (2021).
Article CAS PubMed Google Scholar
Drokhlyansky, E. et al. The human and mouse enteric nervous system at single-cell resolution. Cell 182, 1606–1622 (2020).
Article CAS PubMed PubMed Central Google Scholar
Li, L. et al. Visualizing the distribution of synapses from individual neurons in the mouse brain. PLoS ONE 5, e11503 (2010).
Article PubMed PubMed Central Google Scholar
Kaelberer, M. M. et al. A gut–brain neural circuit for nutrient sensory transduction. Science 361, eaat5236 (2018).
Article PubMed PubMed Central Google Scholar
Haber, A. L. et al. A single-cell survey of the small intestinal epithelium. Nature 551, 333–339 (2017).
Article CAS PubMed PubMed Central Google Scholar
Wang, L., Yuan, P. Q. & Tache, Y. Vasculature in the mouse colon and spatial relationships with the enteric nervous system, glia, and immune cells. Front. Neuroanat. 17, 1130169 (2023).
Article CAS PubMed PubMed Central Google Scholar
Fu, Y. Y. et al. 3-D imaging and illustration of mouse intestinal neurovascular complex. Am. J. Physiol. Gastrointest. Liver Physiol. 304, G1–G11 (2013).
Article CAS PubMed Google Scholar
Bergner, A. J. et al. Birthdating of myenteric neuron subtypes in the small intestine of the mouse. J. Comp. Neurol. 522, 514–527 (2014).
Article CAS PubMed Google Scholar
Bertrand, P. P. & Thomas, E. A. Multiple levels of sensory integration in the intrinsic sensory neurons of the enteric nervous system. Clin. Exp. Pharm. Physiol. 31, 745–755 (2004).
Article CAS Google Scholar
Cooke, H. J., Sidhu, M. & Wang, Y. Z. Activation of 5-HT1P receptors on submucosal afferents subsequently triggers VIP neurons and chloride secretion in the guinea-pig colon. J. Auton. Nerv. Syst. 66, 105–110 (1997).
Article CAS PubMed Google Scholar
Cox, H. M. et al. Multiple Y receptors mediate pancreatic polypeptide responses in mouse colon mucosa. Peptides 22, 445–452 (2001).
Article CAS PubMed Google Scholar
Cooke, H. J. et al. Substance P as a mediator of colonic secretory reflexes. Am. J. Physiol. 272, G238–G245 (1997).
CAS PubMed Google Scholar
Vanner, S. & Surprenant, A. Cholinergic and noncholinergic submucosal neurons dilate arterioles in guinea pig colon. Am. J. Physiol. 261, G136–G144 (1991).
CAS PubMed Google Scholar
Pham, T. D., Gershon, M. D. & Rothman, T. P. Time of origin of neurons in the murine enteric nervous system: sequence in relation to phenotype. J. Comp. Neurol. 314, 789–798 (1991).
Article CAS PubMed Google Scholar
Rivera, L. R. et al. The involvement of nitric oxide synthase neurons in enteric neuropathies. Neurogastroenterol. Motil. 23, 980–988 (2011).
Article CAS PubMed Google Scholar
Burns, A. J. et al. White paper on guidelines concerning enteric nervous system stem cell therapy for enteric neuropathies. Dev. Biol. 417, 229–251 (2016).
Article CAS PubMed PubMed Central Google Scholar
Fattahi, F. et al. Deriving human ENS lineages for cell therapy and drug discovery in Hirschsprung disease. Nature 531, 105–109 (2016).
Article CAS PubMed PubMed Central Google Scholar
Fung, C. et al. Cholinergic submucosal neurons display increased excitability following in vivo cholera toxin exposure in mouse ileum. Front. Physiol. 9, 260 (2018).
Article PubMed PubMed Central Google Scholar
Brierley, S. M. & Linden, D. R. Neuroplasticity and dysfunction after gastrointestinal inflammation. Nat. Rev. Gastroenterol. Hepatol. 11, 611–627 (2014).
Article PubMed Google Scholar
Danielian, P. S. et al. Modification of gene activity in mouse embryos in utero by a tamoxifen-inducible form of Cre recombinase. Curr. Biol. 8, 1323–1326 (1998).
Article CAS PubMed Google Scholar
Leshan, R. L. et al. Leptin action through hypothalamic nitric oxide synthase-1-expressing neurons controls energy balance. Nat. Med. 18, 820–823 (2012).
Article CAS PubMed PubMed Central Google Scholar
Herath, M. et al. Quantitative spatial analysis of neuroligin-3 mRNA expression in the enteric nervous system reveals a potential role in neuronal–glial synapses and reduced expression in Nlgn3^R451C mice. Biomolecules 13, 1063 (2023).
Article CAS PubMed PubMed Central Google Scholar
Bossolani, G. D. P. et al. Comparative analysis reveals Ce3D as optimal clearing method for in toto imaging of the mouse intestine. Neurogastroenterol. Motil. 31, e13560 (2019).
Article CAS PubMed Google Scholar
Jin, X., Oshlack. A. & Phipson, B. phipsonlab/cellXY: cellXY (v1.1.0). Zenodo https://doi.org/10.5281/zenodo.14065058 (2024).
Jin, S. et al. Inference and analysis of cell–cell communication using CellChat. Nat. Commun. 12, 1088 (2021).
Article CAS PubMed PubMed Central Google Scholar
Jin, S., Plikus, M. V. & Nie, Q. CellChat for systematic analysis of cell–cell communication from single-cell and spatially resolved transcriptomics. Nat. Protoc. 20, 180–219 (2025).
Article CAS PubMed Google Scholar

Download references

Acknowledgements

Cell sampling was performed at the ESCG core facility at Science for Life Laboratory, Sweden, funded by the Swedish Research Council. The authors acknowledge support from Science for Life Laboratory, the Knut and Alice Wallenberg Foundation, the National Genomics Infrastructure funded by the Swedish Research Council and Uppsala Multidisciplinary Center for Advanced Computational Science for assistance with massively parallel sequencing and access to the UPPMAX computational infrastructure. U.M. was supported by the Swedish Research Council (Vetenskapsrådet; 2020-01129 and 2022-01570), the European Research Council (divENSify; 101045026), the Knut and Alice Wallenberg Foundation (KAW; 2020.0109), the Strategic Research Area in Neuroscience, the Strategic Research Area in Stem Cell Research and Regenerative Medicine and the Brain Foundation (Hjärnfonden; FO2023-0130 and FO2021-0105). D.R.L. and U.M. were supported by NIH (5R01DK129315; Subaward—KAR-287434). Y.C. and A.M. were supported by the Wenner-Gren Foundation. U.M. thanks the Stellenbosch Institute for Advanced Studies (South Africa) for the creative space and fellowship during her sabbatical in 2023. The imaging was performed in the Biomedicum Imaging Core with support from the Karolinska Institutet. We thank the Viral Vector Facility of the Neuroscience Center Zurich for excellent service and production of AAV vectors. We thank animal staff from the Department of Comparative Medicine at Karolinska Institutet (Komparativ Medicin Biomedicum) for help with animals. We thank R. Hogenboom for optimization and M. Herath for advice on the RNAscope protocol. We thank V. Knoflach, F. Memic and A. Joly for laboratory assistance. We thank P. Ernfors (KI, Sweden) for collaboration on the use of Nmu-Cre and Sst-Cre transgenic mice and F. Lallemend (KI, Sweden) for collaboration on the use of Vip-Cre mice. We acknowledge J. Frisén and S. Giatrellis for access and assistance with flow cytometry at Karolinska Institutet. We acknowledge the staff of the Mayo Clinic Biomedical Imaging Resource for assistance in creating visualizations of the 3D structure of the single neuron and D. Szent-Györgyi for help in engineering the mucosal stimulation assay.

Funding

Open access funding provided by Karolinska Institute.

Author information

These authors contributed equally: Wei Li, Khomgrit Morarach.

Authors and Affiliations

Division of Molecular Neurobiology, Department of Medical Biochemistry and Biophysics, Karolinska Institutet, Stockholm, Sweden
Wei Li, Khomgrit Morarach, Ziwei Liu, Sanghita Banerjee, Yanan Chen, Fernando López-Redondo, Elham Jalalvand, Anastassia Mikhailova & Ulrika Marklund
Enteric Neuroscience Program (ENSP), Mayo Clinic, Rochester, MN, USA
Ashley L. Harb, Joel M. Kosareff, Charles R. Hall, Suad H. Mohamed & David R. Linden
Department of Physiology and Biomedical Engineering, Mayo Clinic, Rochester, MN, USA
Ashley L. Harb, Joel M. Kosareff, Charles R. Hall, Suad H. Mohamed & David R. Linden
RIKEN Center for Integrative Medical Sciences, Yokohama, Japan
Fernando López-Redondo

Authors

Wei Li
View author publications
Search author on:PubMed Google Scholar
Khomgrit Morarach
View author publications
Search author on:PubMed Google Scholar
Ziwei Liu
View author publications
Search author on:PubMed Google Scholar
Sanghita Banerjee
View author publications
Search author on:PubMed Google Scholar
Yanan Chen
View author publications
Search author on:PubMed Google Scholar
Ashley L. Harb
View author publications
Search author on:PubMed Google Scholar
Joel M. Kosareff
View author publications
Search author on:PubMed Google Scholar
Charles R. Hall
View author publications
Search author on:PubMed Google Scholar
Fernando López-Redondo
View author publications
Search author on:PubMed Google Scholar
Elham Jalalvand
View author publications
Search author on:PubMed Google Scholar
Suad H. Mohamed
View author publications
Search author on:PubMed Google Scholar
Anastassia Mikhailova
View author publications
Search author on:PubMed Google Scholar
David R. Linden
View author publications
Search author on:PubMed Google Scholar
Ulrika Marklund
View author publications
Search author on:PubMed Google Scholar

Contributions

K.M., W.L., D.R.L. and U.M. conceptualized and designed the project. W.L., K.M., Z.L., Y.C., F.L.-R., S.B., A.L.H., E.J., J.M.K., C.R.H., S.H.M., A.M., D.R.L. and U.M. conducted data acquisition. K.M., W.L., Z.L., S.B., A.L.H., E.J., J.M.K., C.R.H., S.H.M., A.M., D.R.L. and U.M. analyzed and interpreted the data. W.L., Z.L., S.B., K.M., J.M.K., C.R.H., D.R.L. and U.M. drafted the manuscript. D.R.L. and U.M. secured funding.

Corresponding author

Correspondence to Ulrika Marklund.

Ethics declarations

Competing interests

The authors declare no competing interests.

Peer review

Peer review information

Nature Neuroscience thanks Michael Gershon and the other, anonymous, reviewer(s) for their contribution to the peer review of this work.

Additional information

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Extended data

Extended Data Fig. 1 Identification of non-ENS cells, glial and plausible myenteric cells in the P24 submucosal scRNA-seq dataset.

a, Uniform manifold approximation and projection (UMAP) of level 1 clustering of P24 submucosal cells. See Supplementary Table 1 for differentially expressed genes of the clusters. Jagged boxes are further analyzed in b–f. b, Feature plots indicating epithelial cells based on Cdh1 and various subclass markers and immune-related cells (H2-Aa). c, Analysis of Baf53b-Cre x R26RtdTom mice revealed TOM⁺/5-HT⁺ double-positive cells, showing transgene expression in serotonergic enterochromaffin cells. Image represents observations from 10 fields of view from 2 mice. Feature plots displaying marker genes of clusters representing mesenchymal-derived lineages, including smooth muscle cells (d) and lymphatic epithelial cells (LECs; e). f, Feature plot of Sox10 indicating an enteric glia cell cluster. g, UMAP of cluster level 2, displaying putative contaminating glial and myenteric classes as indicated in feature plots. h,i, Representative submucosa peels showing lack or little ENK (h) and NDUFA4L2 (i) expression in PGP9.5⁺ neurons. 0/2583 neurons expressed ENK; n = 3 mice. 2/7615 neurons expressed NDUFA4L2; n = 5 mice. j, CCK expression was investigated in Cck-IRES-Cre mice injected with AAV-DIO-EYFP or AAV-DIO-tdTOM. Only 1 of 2353 submucosal neurons expressed a fluorescent reporter; n = 4 mice (age, 8–12 weeks). Scale bars, 20 µm.

Extended Data Fig. 2 Quality control of final submucosal neuron scRNA-seq dataset.

a, Box plots showing number of UMIs, detected genes and percentage of mitochondrial genes per cell for each of the cardinal smENCs. Box-and-whisker plots indicate max–min (whiskers), 25–75 percentile (boxes) with median as a center line. smENC1, 503 cells; smENC2, 3770 cells; smENC3, 4068 cells. b, Feature plot depicting cells from the three libraries distributed evenly over the clusters. Note, the same single-cell suspension was used for all libraries. c, Pie chart depicting the number of cells of each sex, and a bar graph verifying presence of female and male cells in all clusters, and their similar ratios between clusters. UMI, unique molecular identifier; smENC, submucosal enteric neuron class.

Extended Data Fig. 3 Differential gene expression in submucosal neuron subclusters.

Supportive data related to Fig. 1i-k. a, Dot plot displaying examples of most differentially expressed genes in submucosal neuron subclusters. b, Feature plot showing preferential expression of PAIP2B in smENC3d. c, Validation of PAIP2B protein expression in a submucosal PGP9.5⁺ neuron (indicated with a star) in duodenum at P28. Image represents observations from 9 fields of view from 5 mice. Scale bar, 20 µm.

Extended Data Fig. 4 Morphological and functional interrogation of smENC1 suggests an IPAN identity.

a, 3D reconstructed image showing one smENC1 neuron from a male mouse (32 weeks), with color-coded fibers within the submucosal plane (endings and fibers), the pericryptal and cryptal spaces and villi. Video representation is available in Supplementary Video 1. b, Schematic drawing of mucosal deformation experimental set-up enabling optical recording of neuronal Ca²⁺ responses in full-thickness tissue collected from ileum of Nmu-Cre;LSL-GCaMP6f adult mice (18–22 weeks). Part of panel b is created with BioRender.com. c, Bar graph indicating the response type of myenteric (n = 42 neurons) and submucosal (n = 34 neurons) NMU⁺ cells (p = 0.2838; chi-square). d, Representative traces of individual epochs from quiescent responding cells demonstrating immediate, delayed and poststimulus responses. e, Bar graph indicating the number of epochs in quiescent responding cells with either type of response pattern in each plexus (p = 0.2503; chi-square). Immediate response myENC6:47/93 epochs; smENC1:40/66 epochs. f, Scatter plot showing the response latency for each epoch and mean (±s.e.m.) in quiescent responding cells of each plexus (p = 0.1235, two-tailed Student’s unpaired t test). g,h, Heat map visualizations of the percentage of quiescent responding cells responding during deformation at each of the 24 stimulated spots in the submucosal (g; n = 23 cells) or myenteric plexus (h; n = 23 cells), which can be used as an estimate of the average mucosal receptive field. i, Distribution histogram indicating the percentage of cells responding to different numbers of epochs in each plexus (p = 0.1984; chi-square). j, Distribution histogram indicating percentage of epochs with different numbers of spikes, showing more spikes in responding smENC1 neurons than myENC6 neurons per epoch (p = 0.0387; chi-square).

Source data

Extended Data Fig. 5 Close proximity of presynaptic and postsynaptic sites indicates synaptic homotypic and heterotypic connections of smENC1–3.

a, Representative picture showing correlation between bouton-like structure (green puncta) and pre-synaptic marker SYN1 (white arrowhead: SYN1⁺, yellow arrowhead: SYN1⁻). b, Schematic drawing showing the viral-mediated strategy to label pre-synaptic sites in submucosal neuron classes. See Supplementary Fig. 5 for validation of strategy. Panel b is created with BioRender.com. c,f,i, Representative pictures addressing potential homotypic communication showing Syp-EGFP puncta surrounding TOM⁺ neuron somas in animals co-injected with AAV-DIO-Syp_EGFP and AAV-DIO-tdTOM (f,i) or AAV-DIO-Syp_EGFP (c). d,e,g,h,j,k, Representative pictures addressing potential heterotypic communication showing Syp-EGFP puncta in close proximity to somas of different neuron classes. l–q, Representative pictures showing proximity of presynaptic (Syp-EGFP) and postsynaptic (PSD95) sites indicating homotypic (n,p) or heterotypic communication (l,m,o,q). White arrowheads—colocalization between Syp-EGFP puncta and PSD95. Asterisk indicates identified potential target cell. Representative images were based on observations in 5 fields of view from each mouse (n = 3 in a,d,e,g,h,j and k; n = 2 in c,f,i,l–q). Scale bars, 20 µm, except for in the enlarged pictures, 5 µm. Mice age, 8–12 weeks.

Extended Data Fig. 6 CellChat analysis addressing ligand–receptor communication between smENC1–3.

Dotplots of CellChat analysis displaying possible ligand–receptor interactions between pairs of smENC1–3. The dot color shows the relative score of communication probability value of each interaction pair. The dot size indicates the level of significance of the interaction. The Benjamini–Hochberg (BH) correction is applied to adjust p-values to control the false discovery rate (FDR).

Extended Data Fig. 7 CellChat analysis addressing ligand–receptor communication between epithelial cell types and smENC1–3.

Related to Fig. 5. a, Reconstructed UMAP of scRNA-seq of epithelial cell types²¹. b, Dot plots of enteroendocrine cells as senders to smENC1–3. c, Dot plots of smENC1–3 as senders to the epithelial cell types. The dot’s color and size indicate the probability score of each ligand–receptor pair interaction and the corresponding p-values, respectively. The Benjamini–Hochberg (BH) correction is applied to adjust p-values to control the false discovery rate (FDR). Note, dot plots for remaining epithelial cells as senders can be viewed in Supplementary Fig 7.

Extended Data Fig. 8 ENS clusters in scRNA-seq dataset at P7 and their differentiation state.

Supportive data related to Fig. 6. a, Inferred male and female cell distribution shown on a pie chart and feature plot. b, UMAP of postnatal day 7 ENS clusters. Legend indicates proposed differentiation state of each cluster as determined based on differentially expressed genes (Supplementary Table 4 and analysis in c). c, Box plots indicating defining genes used to determine the differentiation states of clusters. Box-and-whisker plots indicate max–min (whiskers), 25–75 percentile (boxes) with median as a center line.

Extended Data Fig. 9 Expression of marker genes of different maturity stages or cell types at P7.

Supportive data related to Fig. 6. Feature plots displaying marker genes of different maturity stages or cell types in the developing ENS at P7—a, progenitors; b, neuroblasts; c, neurons; d, early branch markers; e, myenteric neuron class markers (myENCs).

Extended Data Fig. 10 myENC12 is formed through phenotypic transition of NOS1⁺ branch A neurons.

Supportive data related to Fig. 7. Representative images (a) and quantification (b,c) of EGFP⁺ myenteric neurons not maintaining NOS1 expression in Nos1-Cre x R26EGFP mice at W8–12. 12.84 ± 2.50% (n = 3 mice). No difference was detected between intestinal regions. Each animal was assigned a unique symbol in c. Representative images (d,e) and quantification (f,g) of EGFP⁺ myenteric neurons expressing NTNG1/CALB (44.84 ± 5.62%; n = 5 mice) or 5-HT (45%; 20 cells in 5 mice) in Nos1-Cre x R26EGFP mice at W8–12. No difference was detected between intestinal regions. Each animal was assigned a unique symbol in g. h, Representative image showing no overlap between EGFP and the branch B neuron marker CALR. 0 EGFP⁺ of 1903 CALR⁺ cells (n = 3 mice). i, Schematic drawing indicating observed and suggested EGFP expression, suggesting transient NOS1 expression within branch A in both myenteric and submucosal development. Two-tailed Student’s t test was performed to determine statistical significance. Each dot (n) represents one animal. Data are presented as mean ± s.d. NS, P > 0.05. Scale bars, 10 μm. jej, jejunum; ile, ileum.

Source data

Supplementary information

Supplementary Information

Supplementary Figs. 1–8.

Reporting Summary

Supplementary Data

Supporting data for Supplementary Figs. 1, 4 and 5.

Supplementary Video 1

Video depicting single-cell 3D reconstruction of a smENC1 neuron imaged in the ileum of an Nmu-Cre mouse, injected with AAV-DIO-tdTom. Color-coding of fibers indicates their estimated ___location within the submucosa–mucosa space (Extended Data Fig. 4a).

Supplementary Video 2

Video depicting 3D stack of ileal villi imaged in a Nmu-Cre x R26RtdTom mouse, immunostained for capillaries (CD31⁺LYVE1⁻) and the central lacteal (CD31⁺LYVE1⁺; Fig. 5o).

Supplementary Video 3

Video depicting 3D stack of ileal villi imaged in a Sst-Cre x R26RtdTom mouse, immunostained for capillaries (CD31⁺LYVE1⁻) and the central lacteal (CD31⁺LYVE1⁺; Fig. 5o).

Supplementary Video 4

Video depicting 3D stack of ileal villi imaged in a Vip-Cre x R26RtdTom mouse, immunostained for capillaries (CD31⁺LYVE1⁻) and the central lacteal (CD31⁺LYVE1⁺; Fig. 5o).

Supplementary Table 1

DE genes in level 1 clusters of P24 submucosa.

Supplementary Table 2

DE genes in smENC1–smENC3.

Supplementary Table 3

DE genes in subclusters of smENC1–smENC3.

Supplementary Table 4

DE genes in P7 submucosal ENS clusters.

Supplementary Table 5

Primary antibodies, secondary antibodies and AAV virus used in the study.

Source data

Source Data Figs. 2, 3, 5 and 7

Statistical source data.

Source Data Extended Data Figs. 4 and 10

Statistical source data.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Reprints and permissions

About this article

Cite this article

Li, W., Morarach, K., Liu, Z. et al. The transcriptomes, connections and development of submucosal neuron classes in the mouse small intestine. Nat Neurosci 28, 1146–1159 (2025). https://doi.org/10.1038/s41593-025-01962-x

Download citation

Received: 19 May 2024
Accepted: 27 March 2025
Published: 29 May 2025
Issue Date: June 2025
DOI: https://doi.org/10.1038/s41593-025-01962-x

Subjects

Abstract

Similar content being viewed by others

Main

Results

Three submucosal neuron classes in the small intestine

Two putative secretomotor/vasodilatory classes and one IPAN

smENCs are transcriptionally similar to myenteric ENCs

Class morphologies and proportions across the small intestine

Response patterns of smENC1 suggest an IPAN identity

smENCs communicate through homotypic and heterotypic connections

Connectivity of smENCs within the submucosa–mucosal region

Postnatal scRNA-seq captures submucosal neuron development

Identity acquisition via diversification in neuron branches

Discussion

Methods

Mice

Tissue preparation for scRNA-seq

Juvenile tissue preparation

Postnatal tissue preparation

scRNA-seq

Virally mediated neuron labeling

Tissue preparation for histochemical analysis

Juvenile and adult small intestine sections

Submucosal plexus peels

Postnatal small intestine sections

Immunofluorescence analysis

RNAscope

Confocal imaging

Whole-neuron imaging

Epithelial deformation functional experiment

Cell counting, statistical analysis and reproducibility

Analysis of scRNA-seq data

Clustering analysis

Sex inference, dataset integration and identity mapping

Ligand–receptor analysis

Reporting summary

Data availability

Code availability

References

Acknowledgements

Funding

Author information

Authors and Affiliations

Contributions

Corresponding author

Ethics declarations

Competing interests

Peer review

Peer review information

Additional information

Extended data

Supplementary information

Source data

Rights and permissions

About this article

Cite this article

Share this article

Search

Quick links