Introduction

The post-embryonic, developing intestine remains a dynamic organ that undergoes important molecular shifts as the crypts develop and the cellular composition matures. Exposures to the microbiome and nutrients cause additional programmatic changes. In the mouse, the postnatal period between day 0 and week 4 approximates the developmental sequence of premature human infants around 28-34 weeks’ gestation, as both species experience additional crypt-villous maturation and the development of Paneth cells during these intervals (15). While many groups have described the developmental sequence of gastrulation (610) and gut tube formation (11, 12), as well as embryonic generation of the foregut (13), midgut (1315), and hindgut (16), less is known about the postnatal developmental sequence.

Studying human development during the second half of pregnancy has been challenging, owing to the difficulty of obtaining tissue from fetal intestine after about 24 weeks. Even when tissue can be obtained from premature infants, it is requisite that a disease is present, causing the need for surgical intervention – whether that be a congenital malformation of the intestine (duodenal atresia, gastroschisis, etc.) or an acute disease process (spontaneous perforation, necrotizing enterocolitis). Mouse models are frequently used to help us understand human intestinal development, however, most studies of intestinal development in mice focus on the in-utero period and important questions remain about postnatal development (17). In general, existing postnatal studies have focused on the role of single genes in modulating postnatal development (1820). While there is a rich body of literature examining the impacts of the microbiome on postnatal intestinal development (18, 21, 22), comprehensive analyses of the postnatal intestinal sequence are lacking, especially compared to the breadth of studies on prenatal development. A study in Balb/c mice showed that cell proliferation is robust in the neonatal mouse and exceeds cell loss until around week three, when more adult-like homeostatic kinetics are attained (23).

Necrotizing enterocolitis (NEC) is a devastating disease of intestinal necrotic injury that often leads to systemic illness and shock physiology, and carries a high mortality rate of about 20 percent (24, 25). NEC most frequently occurs in premature infants, peaking around 30 weeks’ postmenstrual age (26). Postnatal mice from day 0-day 14 are used to model NEC, owing to the similarities in epithelial development with the period of susceptibility in premature humans, which peaks at 30 weeks postmenstrual age (27). After 2 weeks postnatal age, mice become resistant to experimental NEC unless other modifications are performed (2729). We are interested in the natural developmental changes across this time frame, which may confer NEC susceptibility or resistance.

Apart from the possible relationships to NEC, we are also interested in understanding fundamental biological processes that occur during intestinal development in the neonatal period. While technological advances have allowed the advent of single-cell sequencing approaches for understanding mouse transcriptomics, isolation and representative preparation of both intestinal epithelium and mesenchyme simultaneously remains a challenge, which can make comparing between segments convoluted. Furthermore, rare cell types and RNAs present at a lower frequency can be difficult to capture in single-cell sequencing (30). Therefore, we performed a time course analysis of postnatal development using bulk RNA-Seq of the full thickness mouse intestine, which includes epithelium and mesenchymal tissues. We also examined differences in gene expression patterns of the duodenum versus ileum. Findings were confirmed using other modalities and human tissue biopsies were used for comparison in select instances.

Results

The quality of our RNA-Seq dataset was high, and samples clustered together by timepoint

To examine the quality of the RNA-Seq dataset, we examined the total reads per sample; read numbers ranged from 56 to 185 million reads per sample (Supp. Fig. 1A). In general, ileal samples had fewer reads than duodenal ones. We also evaluated the exonic mapping rate, and most of the samples met our target threshold of 75%, although four of the ileal samples had mapping rates slightly below 70% (Supp. Fig. 1B). Principal component analysis (PCA) of all samples showed that duodenal samples and ileal samples intermixed and did not segregate by tissue type in early development but segregated more at later time points (Supp. Fig. 1C), suggesting that these regions were similar in gene expression through the first two weeks of life.

To analyze the clustering within each group, we also performed PCA for samples from each tissue independently. Duodenal samples generally clustered well by time point and followed a temporal sequence in the principal driver of variance, PC1 (Fig. 1A). The samples from the ileum similarly clustered well and followed a temporal pattern in PC1 (Fig. 1B). Hierarchical clustering, which groups the samples by similarity, showed that the 3-week and 4-week samples were distinct from the other groups in the duodenum, while the day zero samples were also dissimilar from other groups (Fig. 1C, 1D). The 2-day through 2-week samples were more closely related to each other. In the ileum, the 3- and 4-week samples again clustered apart from the rest of the group, but the 0d samples clustered similarly to the 2-day through 2-week samples. However, the samples from day 4 did not cluster tightly, particularly in the duodenum, and there was an outlier in the day 6 samples that did not cluster with the other two in duodenum or ileum. Given that the younger timepoints generally clustered tightly, and we had many time points, we opted to drop the day 4 samples and one of the day 6 samples from the rest of our analyses, to focus on the time points that gave us good replication between the samples and had higher reliability.

PCA and hierarchical clustering.

Principal component analysis of all Duodenum samples (A) and all Ileum samples (B). Correlation in gene expression and hierarchical clustering among all Duodenum samples (C) and all Ileum samples (D). . E. & F. Heat map of differentially expressed genes displayed with hierarchical clustering in all Duodenum samples (E) and in all Ileum samples (F). DE genes (time course using LRT) were defined as adjusted P-value <= 0.005.

Analysis of gene expression using significantly differentially expressed (DE) genes (adjusted p-value/q-value <= 0.005) demonstrated that significantly different gene expression in the duodenum was similar in 0-day and 2-day mice, 6-day, 8-day and 2-week mice, and 3-week and 4-week mice, as displayed in heatmap format (Fig. 1E). The similarities between DE gene expression in the ileum were less linear with 0-day and 8-day aligning independently and 2-day, 6-day and 2-week samples appearing more similar in the DE gene profile (Fig. 1F). 3-week and 4-week profiles were also similar, although 4-week samples appeared to have an additional cluster of highly DE genes compared to the 3-week samples (bottom right). See Supplemental File 1 for gene counts and differentially expressed gene lists.

Gene expression trends in the duodenum and ileum highlight key developmental gene expression changes

We next evaluated all differentially expressed genes across time in the duodenum (Supp. Fig. 2; Supp. File 2) and ileum (Supp. Fig. 3; Supp. File 3) and asked whether there were patterns in gene expression trends that were unique in each dataset. Using the top 3000 differentially expressed genes in each tissue (adjusted P-value threshold used <= 0.005), we identified 20 gene sets with correlated gene expression trends. Because we are most interested in the ileum gene changes and how they may be related to NEC susceptibility, we then called out gene set clusters in the ileum that either were low during NEC susceptibility (day 0-day 14) and then increased at week 3-4, or that were elevated in earlier times points and then decreased significantly as the intestine matured. We performed gene ontology (GO) enrichment analyses on these smaller genes sets (Fig. 2) to identify overrepresented GO terms.

Gene ontology enrichment analysis of selected gene expression clusters.

For each gene cluster, the expression pattern, number of genes present, and abundance score are indicated on the left. On the right, the top (top 30, except for cluster 5) overrepresented gene ontology terms in the cluster are indicated.

Notably, genes for immunity, which are already known to be a critical factor in NEC, were present in one of the clusters (cluster 3). Other clusters that increased in the NEC-resistant time points encoded genes for metabolism and intracellular trafficking (cluster 2) and ATP/energy use (cluster 4). Gene groups that were elevated during NEC susceptibility compared to NEC resistant periods included tissue development (angiogenesis, metabolism, catabolism genes) (clusters 1), cell migration genes (cluster 11), DNA repair (cluster 7), and vacuole development (cluster 5). One gene set, cluster 17, was low initially at day 0 and day 2, then increased in the peak NEC time period from day 6 to 2 weeks, then decreased again in weeks 3 and 4, and this was enriched for genes associated with tissue development including extracellular matrix and smooth muscle gene families.

General intestinal markers are stable throughout time points

Given that we used full-thickness biopsies of the intestinal tissue for our analyses, we wanted to ensure that tissue samples were consistent throughout our dataset as an internal validation. We plotted expression of intestine specific markers CDX2 (specification marker) and CDH1 (E-cadherin, epithelial marker) and saw that expression for duodenum and ileum samples was relatively stable over time, while both were increased in duodenum versus ileum (Fig. 3A). The increased villus length in the duodenum likely explains this relative increase. The overall consistency of expression of these ubiquitous markers helps to establish reliable tissue sampling between time points.

Expression of intestinal epithelial and mesenchymal markers.

Gene expression in the duodenum (black) and ileum (red) by normalized gene count across the first month of postnatal life. Highlighted genes include general markers of the intestine (A), intestinal epithelium and secretory populations (B), intestinal stem cells (C), enteroendocrine cells (D), and intestinal mesenchyme (E). Statistical analysis was by 2-way ANOVA compared to day 0. (F) qRT-PCR analysis of select genes from duplicate intestinal samples at the indicated time points. *p< 0.05, *duo – duodenum significant only, *ile- ileum significant only, ***p<0.001.

Epithelial differentiation markers increase over time, while stem cell markers peak around postnatal week 1

One of the genes present in gene cluster 2 was villin, an epithelial marker, suggesting that epithelial expansion is represented by cluster 2. To examine this in more detail, we called out specific markers of epithelial cell sub-types. Expression of epithelial differentiation markers Vil1 (Villin, all epithelium) and Muc2 (Mucin 2, goblet cells) remained relatively low over the first 2 weeks of life, then increased significantly in weeks 3 and 4 (Fig. 3B). Similarly, and as described previously (13, 5), expression of Lyz1 (Lysozyme 1, Paneth cells) also remained low in the first few weeks of life but was statistically increased from week 2 onward. In contrast to the differentiation markers, conventional intestinal stem cell (ISC) markers Lgr5, Olfm4, and Ascl2 peak around 6-8 postnatal days, and are significantly higher in the duodenum than ileum (Fig. 3C).

Meanwhile, expression of markers for enteroendocrine cells, Chga (chromogranin A), Sst (somatostatin), and Gcg (glucagon), did not follow a consistent pattern (Fig. 3D). However, ChgA and Sst were notably increased in the ileum relative to the duodenum.

Mesenchymal markers also demonstrate varied expression by time point

Cluster 11, which had genes that were higher in the NEC susceptibility period, included Pdgfra, a marker for key mesenchymal cell populations. Pdgfra (platelet-derived growth factor A), is expressed in several intestinal mesenchymal cell types as well as vascular cells, and expression peaked at day 0 in the duodenum and then fell off rapidly (Fig. 3E). In the ileum, in contrast, expression was stable until 3 weeks. Expression of Vim (vimentin), a broad marker for the intestinal mesenchyme, was highest at day 0 and then significantly decreased from day 0 by 2 weeks in the ileum and 3 weeks in the duodenum, then stayed lower. Cd81, a marker for stromal telocytes (3133), rose over time in the duodenum but remained stable in the ileum.

Confirmation of key findings by qRT-PCR

To confirm reproducibility of our results across samples and techniques, we collected samples from other mice at select time points and measured RNA of several of the genes of interest using qRT-PCR (Fig. 3F). Indeed, expression of Olfm4, Muc2, and Pdgfra correlated with the gene expression data observed in the bulk RNASeq analysis. Olfm4 expression was highest at day 6, although by qRT-PCR duodenum and ileum levels were similar. Muc2 expression was higher in the ileum throughout and significantly higher at week 3, consistent with RNASeq. Pdgfra expression decreased at the later time points.

Expression of stem cell markers in the developing human intestine

We next sought to compare expression of select epithelial and trophic factors in intestinal epithelium from premature-born infants. With the caveat that samples can only be obtained in the instances of a clinical need for surgery, we used biopsy samples from intact, healthy-appearing tissue from infants who were diagnosed with NEC and had tissue collected from the ileum or jejunum (Table 2). In a 26-week infant, LGR5 expression was already well organized at the crypt base, highlighting the ISCs (Fig. 4A). Expression was still limited to the crypt at 29 weeks but was more widespread at 35 weeks.

Expression of key epithelial genes in human tissues.

A. RNAScope of biopsied human tissue from the indicated gestational ages (Clinical data from the specimens is indicated in Table 2). DAPI is in blue, LGR5 in red. B. RNA expression in human intestinal organoids from 22-week fetus and 37-week infant by qRT-PCR. ** p<0.01.

Average Expression Duodenum and Ileum Defensins by RNASeq

Demographic data from human biopsy samples

We similarly generated human epithelial organoids from human intestinal biopsies from a 22-week fetus and a 37-week newborn (Table 3) and evaluated the expression of intestinal lineage markers (Fig. 4B). Lgr5, the prototypical stem cell marker, was significantly lower in the 37-week organoids, however, Olfm4 expression was higher in the 37-week organoids compared to the 22-week (2nd trimester) ones. It was unexpected that these intestinal stem cell markers did not correlate. Other lineage markers (Lyz1, Agr2, Atoh1) were similar between the groups (not shown), however we did not use differentiation media so it is expected that the culture medium would favor maintenance of stem cells over epithelial differentiation, as this may related to the similar expression of differentiation markers at both time points.

Demographic data from human organoids

Expression of Wnts and other epithelial regulating molecules also vary over time and by region of intestine

Given the patterns of change of epithelial markers, we next investigated regulators of epithelial differentiation. Further, the non-canonical Wnt Wnt5a was listed in cluster 7, which included genes that were increased during the NEC-susceptibility period. Only 5 members of the Wnt family were expressed at appreciable levels in the B6 mouse small intestine: Wnt2b, Wnt3, Wnt4, Wnt5a, and Wnt9b. Our lab has previously shown that Wnt2b expression is lower than Wnt3 in the small intestine compared to the colon, but duodenal and ileal expression were similar (Fig. 5A). Interestingly, Wnt3 expression peaked at the same time point as the ISC markers, at day 6 in the duodenum, while it rose by day 6 in the ileum and then plateaued. Interestingly, while Wnt3A is used in many commercial intestinal organoid preparations, it is not expressed in the mouse small intestine. Wnt5a, the prototypical non-canonical Wnt, was expressed more highly in the ileum, while Wnt9b was predominantly expressed in the duodenum during the first two postnatal weeks.

Expression of intestinal Wnts and differentiation factors.

Gene expression in the duodenum (black) and ileum (red) by normalized gene count across the first month of postnatal life. Highlighted genes include Wnt family members (A) and an array of downstream regulators and trophic factors known to be important for intestinal epithelial regulation (B). Statistical analysis was by 2-way ANOVA compared to day 0. (C) qRT-PCR analysis of Wnt3 from duplicate intestinal samples at the indicated time points. *p<0.05, *duo – duodenum significant only, *ile-ileum significant only. **p<0.01

Expression of beta-catenin (Ctnnb1) was relatively stable over time in both segments, although it was significantly higher in the duodenum at week 2 and there was some fluctuation in the ileum (Fig. 5B). This is not altogether surprising, given that the phosphorylation state of the protein is typically the point of regulation of beta-catenin function, rather than total level of RNA or protein present (3436). Vangl1 and Vangl2 are involved in non-canonical Wnt signaling, and both of these peaked in the duodenum at days 6 and 8 and were higher in the duodenum than ileum (in contrast to the prototypical non-canonical Wnt, Wnt5a, which was expressed at higher levels in the ileum). Vangl2 was not expressed in the ileum. Expression of Notch1 and Notch2, which function in lineage specification, peaked in the duodenum at day 6 but remained stable in the ileum. Dll1 (Delta-like 1), a Notch ligand, was also stable in the ileum but was significantly lower at day 6, day 8, and week 3 in the duodenum.

We again assessed some of these genes using qRT-PCR to assess reproducibility and consistency between assay types (Fig. 5C). Wnt3 expression was indeed highest at day 6, but, similar to what we observed for Olfm4 (Fig. 3F), expression was comparable in the duodenum and ileum.

Immune genes are up-regulated in the mature duodenum and ileum compared to the perinatal intestine

Given that cluster 3 was highly enriched for immune-related genes, we more closely evaluated markers of the immune system. When we compared the genes most differentially expressed in the week 4 samples compared to the day 0 samples, immune genes were enriched among the significantly increased genes in both the duodenum and ileum (Fig. 6A). Several defensins and CD19 were among the increased genes in both tissues. Lyz1 expression was also significantly increased in the ileum. In the duodenum, insulin signaling-related genes (Igf2, Igf2bp1) decreased, while oxygen responsive genes also decreased in both tissues (Hif3a).

Expression of immune genes.

A. Volcano plot of all differentially expressed genes from day 0 to week 4 in the duodenum and ileum. Highlighted genes are enriched for immune regulators. B. Gene expression of all alpha defensins in the duodenum and ileum by normalized gene count across the first month of postnatal life. Expression of adaptive immune cell markers (C) and innate lymphoid cell markers (D) by normalized gene count. (E) qRT-PCR analysis of Rora from duplicate intestinal samples at the indicated time points. (F) Expression of macrophage markers by normalized gene count. We used one-way ANOVA within each group and indicated significant values within the group as compared to day 0. *p<0.05; **p<0.01; ***p<0.001; ****p<0.0001; *duo – duodenum only; *ile - ileum only

RNA expression of antimicrobial defensins and adaptive immune cells increase over time

Cluster 17 and the volcano plots comparing expression across all time points were highly enriched for defensins, so we also looked at these in isolation. We first compared expression of alpha defensins over time, and all defensins expressed in the duodenum and ileum increased over time (Fig. 6B, Table 1). Overall expression of defensins was higher in the duodenum than the ileum, up to 3-fold higher for some of the genes. At week 4, Defa24, Defa3, and Defa23 were the highest expressed in the duodenum, whereas Defa2, Defa23, and Defa 20 were the highest in the ileum. We also evaluated markers for adaptive immune cell populations (Fig. 6C). CD3g (T cells), and CD8 (cytotoxic T cells) all remained low until 4 weeks. CD4 (Th cells, intestinal macrophages) demonstrated low expression throughout the first postnatal month. CD19 (B cells) also increased significantly at 4 weeks. Expression of T cell and B cell markers was similar between duodenum and ileum.

Expression of key gene markers for innate lymphoid cell (ILC) populations also changed over the postnatal time course (Fig. 6D). Tbx21/Tbet remained low until 4 weeks postnatal, suggesting that ILC1 do not populate the small intestine until later in postnatal life. Rora, a marker for ILC2, was higher in the immediate perinatal period and decreased over time. However, while expression fell by postnatal day 6 in the duodenum, levels remained higher in the ileum through week 2. Rorc, a marker for ILC3, remained more stable across development, but there was an increase in both duodenum and small intestine at 2 weeks postnatal. Tbx21 is also expressed in inflammatory T cell populations cells, and Rorc in pro-fibrotic, IL-17 producing T cells, so some of the increase in Rorc expression may also be coming from T cells. Studies on ILC development in the postnatal intestine have been limited. We used qRT-PCR to confirm the expression pattern we saw for Rora, and indeed the ileum expression remained high until 3 weeks postnatal, while the duodenal expression began to decrease by 1 week postnatal (Fig. 6E).

The intestine is a site of active macrophage engagement, both from systemic macrophages and tissue-resident macrophages, which have different ontogenies (37, 38). CD11b, a marker of innate immune cells (39, 40), remained fairly stable across development (Fig. 6F). Expression of macrophage-specific marker Adgr1 (F4/80) was also relatively stable, with the exception of a dip in expression in the duodenum at day 8. Meanwhile Tim4, a marker of tissue-resident macrophages (41, 42), was also stable, suggesting little fluctuation in the size of this population across early development.

Immunofluorescence analyses confirm and contrast our RNA findings

To evaluate whether the RNA changes we saw were consistent with protein changes in the developing mouse, we performed immunofluorescence (IF) of several key markers in the developing ileum. Lyz1 expression, a marker of Paneth cell development, was low until week 3 (Fig. 3B), and IF confirmed that the protein production matched this (Fig. 7, column 1, Sup. Fig. 4A), as others have described previously (43). Interestingly, OLFM4 expression, a marker of intestinal stem cells, (Fig. 7, column 2; Supp. Fig. 4B) also increased steadily over time in the ileum, which was different than the RNA expression pattern, where peak Olfm4 was seen at day 8 (Fig 3C). It was more similar to the human qRT-PCR data, which showed that organoids from older infants had higher Olfm4 expression. We also examined the duodenum by IF and found the same pattern as in the ileum, where protein expression increased over time and was maximal at week 4 (not shown) while the RNA expression was much higher at day 6.

Protein expression of key markers in the ileum using immunofluorescence.

Ileal tissue was stained with E-cadherin (green) and DAPI (blue) as well as Lyz1 (grey, column 1), Olfm4 (red, column 2), MHCII (grey, column 3), and CD19 (magenta, column 4) at the indicated time points. Scale bar = 100 microns.

Recent advances have shown a critical role for intestinal epithelial MHCII expression in several disease processes (4447). Therefore, we also wondered whether MHCII was expressed by the epithelium in early development and whether it could be involved in NEC. Indeed, MHCII positive cells did not emerge until week 4, when MHCII was highly expressed by certain crypt cells, likely to be ISCs based on prior reports and overlap with the OLFM4 distribution (44). Finally, to analyze the production of CD19, which was expressed in early time points but increased significantly at week 4, we stained tissue sections for CD19 (Fig. 5, column 4). Indeed, Peyer’s patches were small in area at day 0 and consistently enlarged over time, as described (48). Interestingly, early postnatal Peyer’s patch (PP) tissue was highly cellular with DAPI+ cells and few CD19+ cells, and the proportion of CD19+ increased over time as a percentage of the total PP lesions, (Supp. Fig. 3C) indicating enlargement of the germinal centers over time. Postnatal development of PPs is not well described (49).

The biological processes associated with genes differentially expressed in duodenum versus ileum changed over time

Given that the ileum has a higher incidence of NEC than the duodenum, we next evaluated changes in gene expression patterns between the duodenum and ileum. Comparing genes differentially expressed in the duodenum versus ileum at all time points showed a predominance of genes related to metabolism in at day 2, day 6, 2 weeks, and 3 weeks (Supp. Fig. 3; Supp. File 4). The exception was at day 6, where there was an increase in genes related to cell proliferation, cell structures (tubules, subunits), and phosphorylation compared to the ileum. At 4 weeks, the differentially expressed genes shifted to show a predominance of immune-related gene families in the duodenum compared to ileum (Supp. Fig. 3).

Components of extracellular adenosine metabolism are altered in the neonatal period

Because our ATP/energy use transcript expression in cluster 4 of our ileum dataset correlated with the NEC susceptibility period from 2 weeks forward, we considered ATP metabolism in NEC. Adenosine deaminase increased over time, and ADA this can contribute to both production of intracellular and extracellular adenosine (50). Others have shown that extracellular adenosine is important for immunomodulation and is a potent anti-inflammatory mediator in the intestine (51, 52). Indeed, a study in mice showed that extracellular adenosine production can modulate the development of NEC (53). Extracellular adenosine is generated by cell surface molecules, where CD39 first catalyzes the production of ADP from ATP, and then CD73 generates adenosine from ADP (Fig. 8A). Receptors for extracellular adenosine are in the Adora family. Indeed, analyses of these genes in the duodenum and ileum showed that Cd73/Nt5e expression was significantly lower in early postnatal life compared to NEC-resistant periods at 2-4 weeks postnatal (Fig. 8B). In contrast, Adora2b expression started higher and decreased significantly over time in the ileum. Adora1 and Cd39/Entpd1 expression was relatively unchanged over time in the ileum but both appeared to be higher in the duodenum compared to ileum and then decrease over time. We used qRT-PCR to confirm these patterns (Fig. 8C) and saw that Cd73/Nt5e did tend to be higher in the ileum and decreased with maturation, whereas Adora2b expression decreased with ileum maturation. Differences between duodenum and ileum were less pronounced by qRT-PCR than by RNA-Seq, as we saw with earlier qPCR comparisons.

Extracellular adenosine metabolism genes are differentially expressed in immature versus mature ileum.

(A) Schematic of extracellular adenosine signaling. (B) Gene expression in the ileum of critical components of the extracellular adenosine pathway by normalized gene count across the first month of postnatal life in mice. (C) qRT-PCR analysis of select genes from duplicate intestinal samples at the indicated time points. (D) 10X magnification images of human epithelial organoids from the indicated age donors (Clinical data from the specimens is indicated in Table 2). (E) Gene expression in the human enteroids by qRT-PCR. Statistical analysis was by 2-way ANOVA compared to day 0. *p<0.05; **p<0.01; ***p<0.001

Again, using the 22-week and 37-week organoids as in Fig. 6A, we analyzed expression of these adenosine metabolism components in human organoids. We also added an additional organoid line obtained from a 27-week infant with spontaneous intestinal perforation (SIP), where there is likely to be some inflammation present in the tissue. The organoids all grew well in culture and expanded normally (Fig. 8D). The expression of CD73 was not different in 22-week versus 37-week organoids, but trended to be lower in the SIP organoids. (Fig. 8E). Expression of Adora receptors was significantly different, with the lowest expression in the near-term organoids and the highest expression in the SIP organoids.

Discussion

In this study we generated an RNA catalogue of the full-thickness duodenum and ileum in the postnatal intestine of mice and analyzed gene expression patterns over time. Our study is unique in its inclusion of many time points throughout the first month of life. Previously published studies on postnatal intestinal development used targeted qPCR and IF, featured fewer time points, and focused on select functions or groups of genes (5456). We also separately analyzed duodenum and ileum, and highlighted that some genes follow a similar pattern in proximal and distal small intestine, while others have distinct RNA expression patterns between the regions.

This study queried whether gene expression patterns from the dataset would inform any new ideas as to why the postnatal ileum is sensitive to NEC in the first two weeks of postnatal life. This approach highlighted genes associated with certain biological processes that are already known to be important in the development of NEC (Fig. 2), including the development of antimicrobial and immune genes. This dataset adds to previous studies by showing that there are developmental changes in the expression of key genes in these biological processes, which may be relevant to the developmental susceptibility to NEC. We also identified additional overrepresented biological processes that remain relatively unexplored (cellular metabolism, DNA repair, vacuole/phagolysosome function).

Our lab is interested in ISCs and epithelial homeostasis, particularly the role of WNT2B in control of epithelial development and homeostasis. Indeed, there were important differences in epithelial marker genes and regulators of epithelial maturation across this time course. Some of these cell lineages (Paneth cells/Lyz1) have been previously implicated in NEC, but the contributions of other lineages and the regulation of intestinal stem cells in the postnatal intestine are not well-studied. We have previously shown that WNT2B is required for normal numbers of ISC (57, 58), and that WNT2B deficiency leads to increased inflammatory cytokine production in the colon (59). We hypothesized that relative ratios of WNT2B versus WNT3 in the small intestine versus colon explain why the small intestine is less affected by WNT2B loss of function, at least in mice. In this dataset, Wnt3 expression was substantially higher than Wnt2b expression in both the duodenum and ileum, with the exception of immediately postnatal at day 0 and day 2 (Fig. 5A). This highlights the need to consider zonation when considering the roles of trophic factors on epithelial maintenance.

Due to a high degree of homology and redundant functions, Wnt3 and Wnt3a are similar in intestinal function, such that most intestinal organoid medias contain Wnt3a. However, our data show that Wnt3a is not expressed in the mouse SI, at least through 1 month of age (Fig. 5). These Wnts are highly homologous and share some functional redundancy (6062), however, their nucleotide sequences are only 55% identical with large regions of dissimilarity and whether they have any functional differences in controlling intestinal homeostasis is yet to be interrogated. Researchers who use organoids should be aware of this nuance. Interestingly, WNT5A, the prototypical non-canonical Wnt, was more highly expressed in the ileum than the duodenum, in contrast to other expressed Wnt family members. The implications of this are not yet clear.

The immune system is rapidly developing in the postnatal period, as the intestine first encounters the microbiome and is exposed to both microbial and food-related antigens. Despite initiation of milk feeding within hours after birth, adaptive immune cell markers did not robustly increase until week 4 (Fig. 6C). Whether this is related to the introduction of standard rodent chow at day 21 in our mice or the resultant decrease in breastmilk consumption, or something else entirely, cannot be determined from our approach. Intestinal macrophage markers remained static across the first postnatal month (Fig. 6E), which may be due to the predominance of tissue-derived macrophages in this organ (41, 63). The changes in expression of markers specific for ILC populations were interesting. Not unexpectedly, the ILC1 markers were expressed on a much lower level than ILC2 or ILC3 markers (Fig. 6D). Expression of ILC2 marker Rora was higher in early time points, and expression was more persistent in the ileum compared to the duodenum. No dedicated gene marker has been identified to date for ILC1 and ILC3, however, so expression of these markers may be confounded by T cell gene expression.

Peyer’s patch size increased over time, as previously described, but interestingly the proportion of B cells within the Peyer’s patches (PPs) increased significantly over time, suggesting rapid growth of the germinal centers in the postnatal period (Fig. 7, Supp. Fig. 4). Others have shown that the postnatal microbiome drives development of PPs, and that early postnatal antibiotic exposure decreases the number of intestinal B cells in flow cytometric analyses, however they only examined samples from P20 onward. The effects of antibiotics, microbial composition, and other postnatal exposures on expansion of GCs requires additional study.

Finally, we considered the changes in gene expression for ATP metabolism (Fig. 8) and how they may relate to extracellular adenosine. We found that expression of CD79 receptor, required for extracellular adenosine production, was low in early development and increased during the NEC-resistant period after 2 weeks. Conversely, expression of adenosine receptors was higher in the NEC-susceptible period. Interestingly, others have shown that adenosine levels are high in neonatal blood (64), however, knowledge about intestinal adenosine is limited. As others have shown that exogenous adenosine, or inhibition of the adenosine receptor, can increase protection against NEC (53, 65), it will be important to characterize baseline ileal adenosine metabolism and whether it is perturbed in prematurity.

The main limitation of this study is that most data were generated by bulk RNA-Seq. We have done some select confirmatory studies using qRT-PCR and immunofluorescence but have not recapitulated the majority of the dataset by other means. This dataset is intended as a starting point for our lab and others when considering novel hypotheses, and secondary approaches should always be used to further interrogate the RNA expression in a given experimental system. Furthermore, total mRNA does not equate to transcriptionally active mRNA or to protein levels, so other assays would be needed to remark on those aspects of translational control.

Another caveat is that we intended to collect data from day 4, but these samples were not of high enough quality, so we omitted that data. Some dynamic changes occur for ISC genes around day 6, so it is not yet clear whether these markers change soon after day 2 or closer to day 6. Finally, this study is done in C57Bl/6 mice, so applicability to other backgrounds cannot be assumed.

It is important to keep in mind that these samples are all under homeostatic conditions, and the absence or abundance of a gene in this dataset may not be reflective of its expression during stress or injury conditions. Furthermore, exposures are changing over the period we queried, so microbial changes and nutrient changes are occurring but should be relatively consistent in all postnatal mice. We were curious about the natural sequence of events in the normal postnatal environment under routine conditions, so we did not manipulate any of these exposures. Of course, the microbiome is well known to play an important role in development of the epithelial transcriptome and in postnatal methylation, as other groups have demonstrated (66).

In sum, the first postnatal month of life is one of the most dynamic times for intestinal growth and crypt expansion. Characterizing the normal developmental sequence can help us to identify developmental aberrations in neonatal intestinal diseases that may drive pathogenesis. We hope this dataset can serve as a reference point for others studying specific genes or biological processes during postnatal health and disease.

Methods

Research Ethics

The study was approved by Boston Children’s Hospital Institutional Review Board for the human subjects research (BCH IRB 10-02-0053 and P00027983) in accordance with the Declaration of Helsinki. All animal experiments were done in accordance with Boston Children’s Hospital Institutional Animal Care and Use Committee (IACUC).

Sex as a biological variable

The role of sex was not considered as a biological variable in this study, and male and female samples were used indiscriminately for all analyses. This was mainly owing to the high cost of RNA-Seq and the difficulty in telling male from female mice at the earliest time points of the study. We do provide data on sex of the human samples, when able.

Reagents

A table of reagents used in this study, including antibodies and qPCR primers, is included in Supplementary Table 1.

Mouse model

C57/Bl6 mice (Jackson labs) mice were used for experiments, and for RNA-Seq experiments, we used large litters and collected mice across postnatal ages to the extent possible (usually 2-3 time points per litter), to limit potential variation in microbial exposures. Mice were housed in standard germ-free conditions and kept with the dam until 21 days postnatal and then transitioned to standard rodent chow. We used both male and female mice for all experiments.

Time course sample collection

Mice were euthanized at the following time points: day 0, day 2, day 4, day 6, day 8, 2 weeks, 3 weeks, and 4 weeks. For the purposes of RNA-Seq, 2cm of proximal intestine immediately after the pylorus was collected for the “duodenum” sample and 2cm of distal intestine immediately before the ileocecal junction was collected for the “ileum” sample. For the histology samples, we isolated the full small intestine and then transected at the mid-way point to create proximal (duodenal-jejunal, called “duodenum” throughout) and distal (“ileum”) intestinal segments. The intestine was cleaned of feculent material using cold PBS and 4% paraformaldehyde (PFA), and then tissue was opened longitudinally, rolled, and placed in a fixation cassette (Swiss roll method). All samples were Swiss-rolled with orientation of the rolls in a proximal to distal fashion, in order to be able to compare consistently between segments.

RNA-Seq

Intestinal tissue was placed into Trizol LS (ThermoFisher) and homogenized with a bead homogenizer. RNA was isolated using the Direct-zol kit (Zymo Research) according to the manufacturer’s instructions. The library was prepared by the Dana Farber Next Generation Sequencing Core Facility using a KAPA library quantification kit (Roche). Sequencing was done on an Illumina NovaSeq 6000 with total RNA.

All samples were processed using a bulk RNA-seq pipeline implemented in the bcbio-nextgen project (bcbio-nextgen,1.2.8-1c563f1). Gene expression was quantified using Salmon (version 0.14.2) (67) using the mm10 transcriptome (Ensembl). Differentially expressed (DE) genes between the duodenum and ileum at specific time points were identified using DESeq2 (version 1.30.1) (68). The DE threshold was set at adjusted P-value <= 0.005. All analyses performed in R utilized R version 4.0.3.

Analysis of Gene Expression Families Over Time

We analyzed gene expression over time on the duodenum and ileum samples separately using the LRT method in DEseq2. After identifying differentially expressed (DE) genes, we the degPatterns() function from the DEGreport package (version 1.34.0, citation - Pantano L (2024). DEGreport: Report of DEG analysis. R package version 1.42.0, http://lpantano.github.io/DEGreport/.) to identify “clusters” of genes with similar expression patterns across the time points in each tissue using the top 3000 DE genes. We visually analyzed the expression patterns to identify gene groups that fluctuated from the time of NEC susceptibility (less than 2week of age) to the time when mice become resistant to NEC (2 weeks of age onward). We performed functional analysis using the R package clusterProfiler (version 3.18.1) (69, 70) on these gene sets to identify GO biological processes or KEGG pathways that are enriched in each cluster of genes. We then expanded the gene list to include other DE genes with the same pattern of expression of the groups we identified using DEGreport (version 1.34.0sak, citation - Pantano L (2024). DEGreport: Report of DEG analysis. R package version 1.42.0, http://lpantano.github.io/DEGreport/.)), and performed enrichment analysis for GO terms (biological process).

Pathway overrepresentation analyses for duodenum vs ileum

Differentially expressed gene lists of duodenal results relative to ileum at various time points were uploaded to WebGestalt 2019 (71). Selected parameters were overrepresentation analysis, KEGG pathway, and Illumina mouse (13).

Biopsy samples from human infants

Biopsied tissue from endoscopic resections or surgeries was fixed in formalin and embedded in paraffin. All samples were collected solely for clinical purposes, and we used them in a retrospective fashion. Unstained slides were prepared by the clinical histology core at Boston Children’s Hospital.

Establishment of Human organoids

Human intestinal organoids were generated using StemCell® protocol that was adapted as follows. Tissue vials containing about 10 biopsy sized pieces of small intestine that were cryopreserved in 10% DMSO in Fetal Bovine Serum (FBS) were quickly thawed in a water bath. The thawed tissue was transferred into a 15 ml conical containing 5 ml of ice-cold Phosphate Buffered Saline (PBS). The conical was centrifuged at 290 rcf for 5 minutes at 4°C. The supernatant was discarded, and the wash step with PBS was repeated. After supernatant was aspirated, tissue with about 1 mL of PBS was transferred into a 1.5 mL Eppendorf under a sterile hood. Using scissors, the pieces of tissue were minced until they were about 1-2 mm in diameter. Tissue suspended in PBS was transferred to a new 15 mL conical and the PBS was aspirated once the tissue settled.

Propagation of Human Organoids

Organoids were resurrected and expanded from frozen prior to use. Briefly, vials were thawed in a 37°C water bath until just starting to liquefy and then rinsed with Dulbecco’s modified eagle’s media/F12. They were plated in 24-well plates within 50μL of Matrigel (Corning), allowed to solidify in a 37°C incubator for 10 mins, and then 500μL of human small intestinal organoid media was added (Supplementary Table 2). Media was changed every 2-3 days and organoids were passaged roughly once a week until used for experiments. All experiments were done at passage 15 or lower.

Quantitative PCR

The Trizol extraction method was used to purify total mRNA from duodenal enteroid samples. Quantity and quality were assessed using Nanodrop 2000 (Thermofisher Scientific, MA), then cDNA was synthesized using the High Capacity cDNA Reverse Transcription Kit (Thermofisher Scientific, MA). Quantitative RT-PCR was performed using TaqMan Universal PCR MasterMix and Taqman Gene Expression Assays (Thermofisher, MA) for the indicated genes (Supplementary Table 1).

In-situ hybridization

RNAscope Multiplex Fluorescent Reagent Kit v 2 (Advanced Cell Diagnostics, ACD) method was used to localize mRNA expression in FFPE intestinal sections according to a modified version of the manufacturer’s protocol. First, the slides were baked at 60° C for 1 hour and then washed twice with Xylene (Sigma Aldrich) and twice with 100% ethanol for five minutes each. Next, the sections were treated with hydrogen peroxide reagent for 10 minutes RT, washed twice with distilled water, and hot-rinsed with boiling distilled water for 15 seconds, followed by incubation with 1X Target Retrieval Buffer at 99°C for 5 minutes in steamer. Afterward, the slides were hot-rinsed in distilled water, dehydrated in 100% ethanol, and dried for 5 minutes RT. The sections were then incubated with Protease Plus Reagent at 40°C for 30 minutes, washed twice with 1X Wash Buffer, and hybridized with LGR5 probe at 40°C for 2 hours. After hybridization according to the manufacturer’s protocols, the sections were washed twice, incubated with DAPI (1:1000) for 15 minutes RT, washed again, and finally sealed with a coverslip. ACD designed the probe used in this study (Supplementary Table 1). Images were taken using a Zeiss LSM 980 with Airyscan microscope and processed using ImageJ Fiji software(72).

Immunofluorescence

Tissues were fixed in 4% PFA overnight at 4°C, washed in 50% ethanol for 1h, and then transitioned into 70% ethanol. Tissues were embedded in paraffin and sectioned (4-6μm). Slides were prepared for staining using Trilogy buffer (Sigma Aldrich) in a pressure cooker for 15 minutes according to the manufacturer’s recommendation. Slides were then rinsed with boiling Trilogy buffer for 5 minutes, washed with deionized water, and then rinsed further in PBS. Blocking buffer was added to each slide (blocking buffer: 0.3% Tween-20, 10% normal donkey serum, 0.05% bovine serum albumin in PBS) and incubated at room temperature for 1 hour, gently shaking. Primary antibody (Supplementary Table 1) was then added in blocking buffer at 4°C overnight. Slides were then washed in PBS and incubated with the secondary antibody for 1 hour at room temperature. DAPI (1:1000) was added for 10 minutes and then the slides were washed again and sealed with a coverslip. Images were taken using a Nikon Eclipse 90i Fluorescence microscope and processed using ImageJ Fiji software (72).

Statistical Analyses

Statistical analysis of RNA-Seq data is described separately above. For experiments comparing two groups (e.g., duodenum vs ileum), the Student’s t test was used in PRISM 9 (GraphPad Software). For experiments assessing expression changes over multiple time points, ANOVA was used and individual points were compared to the day 0 baseline sample. Data were assumed to be normally distributed. For all experiments, p<0.05 was considered significant.

Data availability statement

The synthesized transcript datasets generated and analyzed during the current study are included in this published article and its Supplementary Information files. The raw sequencing data is available in GEO with accession GSE281618.

Acknowledgements

We would like to thank the infants who contributed tissue samples for this study and their families. The authors would also like to thank Nicholas Makoganov for assisting with some of the technical aspects for this manuscript.

Additional information

Author contributions

Authorship was indicated by order of greatest to least contribution, with the exception of the three senior authors, who were listed last in order of least to greatest contribution. All authors read and critically reviewed the manuscript and approved of the final manuscript. In addition:

LFO performed RNAScope experiments and assisted with data collection and manuscript preparation.

SW performed some experimental assays and data analyses.

VSD contributed to experimental assays for immunofluorescence experiments.

RSK performed the bioinformatics analysis of the RNA-seq data and assisted with drafting and editing the manuscript.

SJHS provided mentorship and bioinformatics expertise, aided with data management, and critically edited the manuscript.

AWH assisted with reviewing human infant cases and identifying appropriate pathologic specimens.

OO generated fetal and human neonatal organoids.

JG provided human histology specimens and expertise.

DTB provided mentorship and topical expertise, and critically edited the manuscript.

LK provided human organoids, provided topical expertise, and critically edited the manuscript.

AEO conceptualized the project, performed experiments, and drafted and edited the manuscript.

Abbreviations

  • ANOVA: analysis of variance

  • ATP: adenosine triphosphate

  • ADP: adenosine diphosphate

  • DE: differentially expressed

  • NEC: necrotizing enterocolitis

  • ISC: intestinal stem cells

  • ILC: innate lymphoid cells

  • IF: immunofluorescence

  • PP: Peyer’s patch

  • PFA: paraformaldehyde

  • MCHII: major histocompatibility complex class II

  • RT: room temperature

  • SIP: spontaneous intestinal perforation

  • PCA: principal component analysis

  • FFPE: formalin-fixed paraffin-embedded

Additional files

Supplemental File 1

Supplemental File 2

Supplemental File 3

Supplemental File 4

Supplemental Figures

Supplemental Table 1

Suplemental Table 2