Single-Cell Analysis Reveals Regional Reprogramming During Adaptation to Massive Small Bowel Resection in Mice.

BACKGROUND & AIMS: The small intestine (SI) displays regionality in nutrient and immunological function. Following SI tissue loss (as occurs in short gut syndrome, or SGS), remaining SI must compensate, or "adapt"; the capacity of SI epithelium to reprogram its regional identity has not been described. Here, we apply single-cell resolution analyses to characterize molecular changes underpinning adaptation to SGS.
METHODS: Single-cell RNA sequencing was performed on epithelial cells isolated from distal SI of mice following 50% proximal small bowel resection (SBR) vs sham surgery. Single-cell profiles were clustered based on transcriptional similarity, reconstructing differentiation events from intestinal stem cells (ISCs) through to mature enterocytes. An unsupervised computational approach to score cell identity was used to quantify changes in regional (proximal vs distal) SI identity, validated using immunofluorescence, immunohistochemistry, qPCR, western blotting, and RNA-FISH.
RESULTS: Uniform Manifold Approximation and Projection-based clustering and visualization revealed differentiation trajectories from ISCs to mature enterocytes in sham and SBR. Cell identity scoring demonstrated segregation of enterocytes by regional SI identity: SBR enterocytes assumed more mature proximal identities. This was associated with significant upregulation of lipid metabolism and oxidative stress gene expression, which was validated via orthogonal analyses. Observed upstream transcriptional changes suggest retinoid metabolism and proximal transcription factor Creb3l3 drive proximalization of cell identity in response to SBR.
CONCLUSIONS: Adaptation to proximal SBR involves regional reprogramming of ileal enterocytes toward a proximal identity. Interventions bolstering the endogenous reprogramming capacity of SI enterocytes-conceivably by engaging the retinoid metabolism pathway-merit further investigation, as they may increase enteral feeding tolerance, and obviate intestinal failure, in SGS.

Here, single-cell RNA sequencing reveals interactions between the retinoid metabolism pathway and ‘regional reprogramming' of distal small intestinal epithelium to a proximal identity following proximal small bowel resection. This provides novel insight into physiological adaptation to short gut syndrome.

The small intestine (SI) absorbs nutrients necessary to sustain life, and displays regional specialization for absorption of specific nutrients along its cephalocaudal axis from duodenum to jejunum to ileum. The majority of nutrient absorption occurs in the duodenum and jejunum, while the ileum is primarily responsible for absorbing bile, vitamin B12, and fat-soluble vitamins. The ileum is also more prone to inflammatory disorders relative to proximal intestine, due in part to higher bacterial load.

A variety of diseases require surgical resection of significant lengths of SI. These range from congenital anomalies and necrotizing enterocolitis in children to trauma, embolism, and malignancy in adults. The resulting loss of SI can cause short gut syndrome (SGS), or the inability of the SI to completely support the metabolic demands of a patient. Management options for SGS are limited, comprising parenteral nutrition (PN), intestinal lengthening procedures, and ultimately small bowel transplant, all incurring significant morbidity and mortality.1, 2

Our murine model of SGS is based on small bowel resection (SBR), in which 50% of the proximal SI is surgically removed. Sham surgery consists of transection and anastomosis, without removal of SI, and acts as a control for exposure to anesthesia, laparotomy, and intestinal transection. This model elicits villus lengthening in the remnant ileum of SBR but not sham mice, yielding increased mucosal absorptive surface area to compensate for lost tissue. Importantly, the degree of this “structural adaptation” response correlates with “functional adaptation,” as evidenced by increased oral tolerance and weight gain in mice. This model is relevant to clinical SGS, as structural adaptation correlates with oral tolerance and weaning from PN observed in human patients. At the same time, structural adaptation does not intrinsically predict functional adaptation, as perturbed weight gain and steatorrhea affect mice deficient in CXCL5 after SBR, despite normal structural adaptation. This leads us to conclude that additional factors beyond simple tissue hyperplasia are at play, which are likely underscored by molecular changes at the single-cell level.

While structural and, to a lesser extent, functional adaptation following SBR has been characterized, relatively little is understood about the molecular changes that accompany the adaptation process. In this respect, messenger RNA (mRNA)– and protein-level expression analyses offer crucial insight, as clinically appreciable adaptation may require that cells assume molecular identities mimicking those of the resected region. Studies of adaptation at the structural and functional level lack the resolution to explore this possibility. Here, in the case of proximal SBR, we hypothesize that remnant ileum (distal SI) upregulates gene and protein expression patterns characteristic of the jejunum (proximal SI) at the single-cell level, a process we term “regional reprogramming.”

Clinical therapies to induce regional reprogramming could enhance an SGS patient’s ability to tolerate oral intake and wean from PN by augmenting the innate functionality of epithelial cells. It is possible this approach may actually be more effective than the previous “holy grail” of SGS research, which has primarily focused on inducing structural adaptation. Enhanced structural adaptation is intrinsically more metabolically demanding (tissue growth) and may or may not affect the key absorptive, metabolic, and immunological pathways specifically deficient in a SGS patient.

To test our hypothesis and address gaps in our understanding of adaptation to SGS, we employed high-throughput single-cell RNA sequencing (scRNA-seq) to characterize gene expression changes of distal SI epithelium during adaptation following massive proximal SBR. This allowed us to dissect population heterogeneity within the epithelium and characterize the regionalization pathways critical in the adaptive response. Here, we show that following SBR, the SI epithelium regionally reprograms toward mature proximal enterocyte identity, accompanied by increased proximal SI nutrient processing gene expression. These changes are punctuated by the increased expression of the proximal SI transcription factor, Creb3l3, a key candidate for reprogramming distal SI gene regulatory networks to a more proximal identity. Analysis of upstream pathways suggests a role for retinoic acid (RA) signaling in driving the adaptation response. Together, our single-cell analyses have enabled high-resolution characterization of the molecular changes and regional reprogramming that underlie adaptation.

Results

Cellular Heterogeneity of the Small Intestinal Epithelium Is Captured by scRNA-seq

Structural adaptation—or villus growth— reliably occurs by day 7 after SBR, making it a commonly utilized experimental end point in SBR literature. As such, we chose day 7 after sham or SBR surgery as our experimental endpoint. We confirmed typical SBR structural adaptation, with villi height increasing by 86.19 ± 19.14 μm (P < .01) relative to sham (Figure 1A). Epithelial cells from animals demonstrating structural adaptation were harvested from SI in an area equidistant from the anastomosis, dissociated into single cells, and processed via high-throughput droplet-based scRNA-seq using the 10x Genomics platform. In total, we sequenced 19,245 cells from 9 independent biological replicates (sham: 8209 cells, n = 5 replicates; SBR: 11,036 cells, n = 4 replicates). A mean of 1767 and 1763 genes per cell in sham and SBR, respectively, and 6754 and 6111 transcripts per cell in sham and SBR, respectively, were detected (Figure 1B).

Figure 1

Experimental design, quality control, and single-cell analysis. (A) A 50% proximal SBR and sham operation were performed on mice. Seven days after surgery, the intestine distal to the anastomosis (ileum) was harvested and equal amounts of tissue equidistant from the anastomosis were used to generate single-cell epithelial suspensions. An area immediately adjacent to this was prepared for histological examination. Typical structural adaptation of SBR mice (lengthened villi) was confirmed (P = .003), with a representative hematoxylin and eosin image of SI tissue from a sham vs SBR mouse shown (20× image acquired using Nikon Eclipse 80i). Scale bar = 100 μm. Epithelium from mice demonstrating structural adaptation was prepared for scRNA-seq analysis. (B) A mean of 1767 and 1763 genes per cell in sham and SBR, respectively, and 6754 and 6111 transcripts per cell in sham and SBR, respectively, were detected. (C) UMAP of integrated biological replicates identified 16 unique cell clusters. (D) Cell cycle states projected onto the UMAP. (E) Representative plot of SBR experimental replicates demonstrated similar gene expression profiles. Correlation coefficient (R) of average gene expression are as shown between these biological replicates. Total biological replicates were 5 sham and 4 SBR (n = 3 “sham1,” n = 1 “sham2,” n = 1 “sham3,” n = 3 “SBR1,” n = 1 “SBR2”). (F) The same 16 clusters were identified in both sham and SBR, in all replicates, as described in panel E. Distribution of all cells across clusters 0–15 (from left to right), by replicate, is shown.

To cluster and visualize cells based on their transcriptional similarity, we used the R package, Seurat.7, 8 and Uniform Manifold Approximation and Projection (UMAP). UMAP analysis and plotting revealed 16 clusters of transcriptionally distinct cell types/states (Figure 1C), where clustering was not driven by numbers of detected genes or transcripts (not shown). Scoring and projection of cell cycle state onto the UMAP plot revealed clustering of cells in S and G1 phases, corresponding to stem cells and transit-amplifying (TA) cells (Figure 1D). Gene expression between equivalent biological replicates was highly correlated, demonstrating a high degree of consistency between the independent biological replicates (Figure 1E). Furthermore, cells from every cluster were represented in each biological replicate, demonstrating consistency of cell capture (Figure 1F).

To assign cell identity to each cluster in an unsupervised manner, we used a computational method based on quadratic programming (QP) to score individual cell identity against an existing single-cell atlas of well-annotated SI cell types. This reference atlas contains stem cells, TA cells, early and late enterocyte progenitors, immature proximal and distal enterocytes, mature proximal and distal enterocytes, goblet cells, Paneth cells, enteroendocrine cells, and tuft cells, annotated based on an extended list of previously identified markers. Scoring cell identity using QP is beneficial because it can capture cells in transitional states, rather than assigning cell identities as binary values. This is ideal for assessing developmental and disease processes where cell identity can be considered as continuous rather than discrete. This flexibility is not found in most classification systems that definitively assign discrete identities to cells, without considering their transitional states. We have previously demonstrated the efficacy of QP in placing cells into an identity continuum during lineage reprogramming. Considering enterocytes progressively differentiate while migrating along the villus axis, and our hypothesis that cell identity is reprogrammed during adaptation, QP represents an appropriate method for quantifying any changes in cell identity that accompany SBR.

Cell identity scores generated by QP were projected onto the UMAP plot, enabling cell cluster identity to be annotated (Figure 2A). This confirmed that all major cell types of the SI epithelium, as previous, were captured by our single-cell analysis, including low numbers of tuft cells and enteroendocrine cells (0.3% and 0.03%, respectively, as a proportion of all captured cells). Projection of all identity scores onto the UMAP plot demonstrated this clustering and visualization method does indeed retain both local and global information, capturing the differentiation trajectory from stem cells to mature enterocytes (Figure 2B). Further examination of these clusters revealed that expression of the proliferation marker Mki67 is enriched in stem, TA, and progenitor cells, and is downregulated as cells begin to differentiate. Conversely, expression of mature enterocyte marker alkaline phosphatase (Alpi) increases in concert with differentiation or maturation, with highest expression colocalizing in areas identified by QP as mature enterocytes (Figure 2C). Furthermore, based on this differentiation trajectory, we posited the most terminal cluster on the UMAP plot would be enriched for markers specific to villus tip cells, such as adenosine deaminase (Ada), which we confirmed (Figure 2C). Finally, one population (cluster 15) did not receive a definitive QP score for any of the reported SI epithelial lineages. This cluster is enriched for intraepithelial lymphocyte marker expression, including Cd45 (Figure 2C), suggesting these cells represent intraepithelial lymphocytes captured alongside the SI epithelium.

Figure 2

Annotation of cell identities using QP. (A) QP-based identity scores of intestinal epithelial populations, projected onto UMAP. Cell populations (from left to right) include stem, early TA, TA G1, TA G2, early enterocyte progenitors, late enterocyte progenitors, immature proximal enterocytes, immature distal enterocytes, mature proximal enterocytes, mature distal enterocytes, goblet cells, and Paneth cells. Low percentages of tuft and enteroendocrine cells (0.3% and 0.003%, respectively) were identified and are not shown. (B) Aggregated QP scores provide a summary of cell identities within the UMAP, demonstrating a maturation trajectory from stem cells to mature enterocytes. (C) Projection of transcript enrichment for selected QP validation markers, clockwise from top left: proliferation marker antigen KI-67 (Mki67) is appropriately enriched in the stem, TA, and progenitor regions of UMAP; alkaline phosphatase (Alpi) expression increases as enterocyte maturation occurs; Cd45 identifies a population enriched for intraepithelial lymphocytes (IELs); adenosine deaminase (Ada), a villus tip marker, localizes at the termination point of the developmental trajectory. Color scale bar indicates relative intensity of cell identity scoring (A) or gene expression (C) across the UMAP.

In summary, our single-cell analyses identified all major SI epithelial cell populations, again demonstrating the strength of our unsupervised QP cell classification approach in concert with UMAP visualization to identify cells differentiating along a defined trajectory. Together, this provides a comprehensive picture of SI epithelial heterogeneity, under both sham and SBR conditions.

Quantification of SI Epithelial Cell Composition Changes Following SBR Reveals Regional Reprogramming Toward Mature Proximal Enterocyte Identity

Changes in cell type composition of the SI epithelium accompany adaptation. However, reports on the precise nature of these changes have been conflicting, with studies describing a relative expansion of either enterocytes or secretory lineages following SBR.14, 15, 16, 17, 18 Discrepancies between these studies may arise due to differences in resection location, amount, the use of differing animal models, and experiment conditions or durations. Furthermore, these previous studies were limited by a lack of resolution to assess cell identity in an unbiased manner. Here, we use our unbiased single-cell resolution classification of cell identity to precisely quantify changes in epithelial composition following SBR. We assessed the distribution of sham- vs SBR-derived cells by projecting the densities onto the UMAP plot (Figure 3A). First, looking at the villus enterocyte population (comprising immature and mature enterocytes), we found a significant increase following SBR, as a total of all cells surveyed (68.9% ± 3.1% of sham events sampled, vs 76.8% ± 0.1% of SBR (P < .05) (Figure 3B). This increase in villus enterocytes, at single-cell resolution, is in agreement with immunohistochemistry analysis of our tissue samples. Here, we identified villus enterocytes as nongoblet (mucin-2, or MUC2 expressing) nucleated villus cells divided by all nucleated cells of the crypt villus axis, revealing that enterocytes comprise 71.6% of sham epithelium vs 74.9% of SBR epithelium (Figure 3C). From our single-cell analysis, we did not observe significant changes (P = .06) in the proportion of non–enterocyte-differentiated cells (Paneth, goblet, enteroendocrine, and tuft cells) comprising sham (22.3% ± 5.8%) and SBR (13% ± 0.3%) SI epithelium (Figure 3B).

Figure 3

Relative expansion of mature proximal enterocytes occurs in SBR mice. (A) Density of sampled cells from sham and SBR epithelium, projected onto the UMAP, demonstrates a population shift toward mature proximal enterocytes in SBR. Color scale bar indicates relative density. (B) Graphical representation of how epithelial lineages (as identified by QP), contribute to the total composition of events sampled by scRNA-seq in sham vs SBR. (C) Immunohistochemical analysis was performed on tissue samples (n = 3 sham and n = 3 SBR) to confirm the relative expansion in villus enterocytes (ie, immature and mature enterocytes) as a percent of total epithelium in SBR, as predicted by scRNA-seq. Representative images of immunohistochemistry for mucin 2 (Muc2, goblet cell marker) in sham vs SBR are provided (20× image acquired using Nikon Eclipse 80i). Scale bar = 100 μm.

Continuing to focus on the enterocyte lineage, as this absorptive cell type is central to adaptation, we next quantified enterocyte differentiation between sham and SBR populations. We found that sham samples contained 9.0% ± 5.5% immature and 59.9% ± 7.2% mature enterocytes, compared with 4.5% ± 2.6% immature and 72.3 ± 2.7% mature in SBR (Figure 3B). Though these differences in immature and mature enterocyte composition between sham and SBR did not reach significance (P = .2 and 0.06, respectively), this observation is supported by a previous report of enhanced metabolically mature enterocyte migration after SBR.

Considering our hypothesis that cells adopt a different regional identity to aid adaptation, we next quantified the balance between proximal and distal enterocyte identities in sham vs SBR samples. As expected considering the tissue was harvested from ileum, a large percentage of cells (58.5% ± 7.5%) from sham samples received high mature distal enterocyte scores with very few cells (1.4% ± 0.6%) scoring as mature proximal enterocytes (Figure 3B). In contrast, in SBR samples we found a significant increase in the percentage of cells receiving high mature proximal enterocyte scores (11.7% ± 4.1%, P < .05) (Figure 3B). This shift toward mature proximal enterocyte identity in SBR suggests a transcriptional “proximalization,” or regional reprogramming, of ileum in response to proximal SI resection.

Regional Reprogramming of Distal SI After SBR Is Accompanied by Increased Proximal SI Nutrient Processing Gene Expression

From our previous analyses, changes in SI epithelial composition following SBR center primarily on a shift toward mature proximal enterocyte identities, suggesting that these changes are a key driver of the adaptive response. Thus, we next focused on characterizing the transcriptional changes underlying the shift toward mature proximal enterocyte identity in SBR.

To identify significant transcriptional changes after SBR, we performed differential gene expression analysis, identifying 174 differentially expressed genes between all sham and SBR cells. The 10 most significantly expressed genes in SBR, relative to sham, are shown in Table 1. This list of SBR-associated transcripts is enriched for signature genes of proximal SI nutrient processing function, including apolipoprotein A-IV (Apoa4), fatty acid binding protein 1 (Fabp1), apolipoprotein C-III (Apoc3), lactase (Lct), and epoxide hydrolase 2 (Ephx2)11, 20 (Figure 4). In contrast, distal SI transcript fatty acid binding protein 6 (Fabp6) was significantly depleted in SBR (0.62 average log fold change depleted, P < .001) (Figure 4).

Table 1

Top 10 Genes Upregulated in SBR Relative to Sham Epithelium, With Average logFC and Adjusted P Values

Gene	logFC	P
Apoa4	2.05	<.001
Fabp1	1.76	<.001
Apoc3	1.46	<.001
Lct	1.11	<.001
Rbp2	1.10	<.001
Fabp2	0.96	<.001
Sepp1	0.92	<.001
Leap2	0.79	<.001
Ephx2	0.71	<.001
Apob	0.69	<.001

logFC, log fold change; SBR, small bowel resection.

Figure 4

Identification of signature proximal small intestine transcripts that increase after SBR. Violin plots (left) show relative expression of transcripts in sham vs SBR populations and UMAP plots (right) show relative transcript expression levels within cell populations.

Increased Apoa4 and Rbp2 expression following SBR has been previously described.21, 22, 23 Contrary to our findings, one of these studies reported no significant changes in ileal Fabp1 mRNA expression. However, this study examined whole tissue preparations, rather than epithelium at single-cell resolution. Thus, to validate the transcriptional changes revealed by our scRNA-seq analysis, we surveyed Fabp1 expression via RNA fluorescence in situ hybridization (RNA-FISH) on histological sections of sham and SBR animals at 7 and 70 days postsurgery, with the latter analysis designed to investigate whether the observed changes are stable, a property rarely investigated in the context of SBR. At day 7, Fabp1 showed a 1.5 average log2-fold AU increased expression per nucleated villus cell in SBR, relative to sham (P < .001), and this response was maintained through day 70 postsurgery (1.4 average log2-fold AU increased expression, P < .01) (Figure 5A).

Figure 5

Validation of proximal small intestine markers that increase after SBR. (A) RNA FISH for Fabp1 shows significant upregulation of transcripts (fluorescein signal) per nucleated cell (DAPI) at days 7 and 70 after surgery. Day 7: n = 15 sham images, n = 15 SBR images (3 biological replicates). Day 70: n = 12 sham images, n =15 SBR images (3 biological replicates) (images acquired using Olympus FV1200 Confocal Microscope). (B) Immunohistochemistry staining for FABP1, FABP6, and SEPP1 in sham and SBR mice at postoperative day 7 shows qualitative increases in FABP1 and SEPP1, and decrease in FABP6, in SBR mice (images acquired using Nikon Eclipse 80i with Ds-Ri2 camera). (C) Immunofluorescence staining images for FABP1 (n = 12 sham images, n = 15 SBR images, 3 biological replicates), FABP6 (n = 18 sham images (4 biological replicates), n =11 SBR images (3 biological replicates)), and SEPP1 (n = 22 sham images (5 biological replicates), n = 18 SBR images (4 biological replicates)) were computationally analyzed to confirm significant protein-level changes corresponding to mRNA changes. White arrows indicate areas of intense SEPP1 expression in SBR (images acquired using Nikon Eclipse 80i with Ds-Ri2 camera). (D) Top to bottom: Western blot analysis of APOC3 in sham and SBR epithelial lysates, normalized to GAPDH and quantified (n = 3 sham and n = 4 SBR mice); qPCR validation of upregulated SBR genes Rbp2 and Ephx2 from SI tissue (n = 3 sham and n = 4 SBR mice). RNA-FISH images are at magnification 60×, scale bar = 30 μm. Immunohistochemistry stains are at 20×, scale bar = 100 μm. IF are at 40×, scale bar = 100 μm. All graphs are presented as mean ± SD. *P < .05, **P < .01, ****P < .0001.

To further validate our findings, we performed quantitative immunofluorescence, western blotting, and quantitative polymerase chain reaction (qPCR) for selected proteins and transcripts, including FABP1, FABP6, APOC3, Rbp2, and Ephx2. Representative immunohistochemistry images of FABP1 and FABP6 show qualitative changes in these proteins consistent with our scRNA-seq analysis (Figure 5B). Quantification of this immunostaining panel using immunofluorescence demonstrated significant increases in FABP1 (1.2 average log2-fold AU, P < .001) in SBR relative to sham (Figure 5C). In contrast, FABP6 was significantly decreased in SBR compared with sham (–1.2-fold relative fluorescent intensity, P < .05) (Figure 5C). Western blotting and quantification of APOC3 showed a 1.92-fold increase in SBR compared with sham (P < .05), and qPCR confirmed upregulation of Rbp2 and Ephx2 (Figure 5D).

Sepp1 upregulation during adaptation to SBR is a novel finding, warranting further investigation. SEPP1 is a secreted glycoprotein with important immunomodulatory and antioxidant effects in the intestine.24, 25, 26 Representative immunohistochemistry images of SEPP1 show qualitative changes consistent with our scRNA-seq analysis (Figure 5B). Quantification of this immunostaining via immunofluorescence demonstrated significant increases in SEPP1 (1.4-fold relative fluorescent intensity, P < .05) in SBR relative to sham (Figure 5C). Furthermore, SEPP1 has been shown to suppress inflammation-associated tumorigenesis, in part through its effects on macrophage polarization, more specifically, by suppressing M2 associated Ym1 expression.25, 27 Interestingly, this is consistent with RNA-sequencing analysis we performed on subepithelial tissue from sham and SBR mice, showing that Ym1 is the second most depleted transcript in SBR mice (–6.19-fold depleted, P < .05, unpublished). This suggests a role for SEPP1 in mitigating oxidative-stress induced injury during adaptation, which is known to occur after SBR, consistent with the fact that FABP1 directs fatty acids toward oxidative metabolism, possibly via immunoregulatory effects on macrophages, and merits further study.

Together, these orthogonal validations confirm our scRNA-seq results, demonstrating that gene expression programs to support proximal SI nutrient processing, and to counteract the associated oxidative stress, are engaged following SBR. These mRNA- and protein-level changes underlie the regional reprogramming to mature proximal identity we observe in SBR.

Proximal SI Transcription Factor Creb3l3 Shows Stable Upregulation Accompanied by Expanded Villus Axis Zonation in SBR Mice

We next aimed to identify candidate transcription factors responsible for driving the observed regional reprogramming and shift toward a mature proximal nutrient processing profile following SBR. Of 44 transcription factors previously shown to be differentially expressed between proximal vs distal enterocytes, only one associated with proximal identity, cAMP responsive element binding protein 3 like 3 (Creb3l3), was upregulated in SBR (0.46 average log-fold increase in SBR, P < .0001). Creb3l3 is a master regulator of lipid metabolism, fitting with increased lipid metabolism after SBR. An additional transcription factor, Kruppel-like factor 4 (Klf4) was also increased in SBR (0.26 average log-fold increase in SBR, P < .001). Though Klf4 was identified as a distal enterocyte transcription factor by Haber et al, another group reported highest Klf4 expression in duodenum and jejunum, and described that Klf4 plays a role in the maturation of intestinal stem cells, as well as in the differentiation of absorptive lineages, such as enterocytes. It is worth noting the possibility that a wider variety of proximal transcription factors were upregulated immediately after surgery and stabilized to baseline by day 7, when structural adaptation was complete. However, since these were the only previously identified regional enterocyte transcription factors significantly differentially expressed in SBR at day 7, we focused our analysis on these factors as putative drivers of stable regional reprogramming.

First, we sought to determine whether proximal SI transcription factor Creb3l3 expression was transiently upregulated in SBR, or whether it was critical to maintaining a long-term adaptive response. As Creb3l3 regulates lipid metabolism, and the increased demand for ileal lipid absorption should persist indefinitely after SBR, we expected its expression to remain elevated. Indeed, RNA-FISH for Creb3l3 at days 7 and 70 demonstrated its significant and long-term upregulation following SBR (day 7: 1.3 average log2-fold AU increase in SBR, P < .001; day 70: 1.2 average log2-fold AU increase in SBR, P < .01) (Figure 6A).

Figure 6

Dissecting genetic underpinnings of epithelial proximalization following SBR. (A) RNA-FISH for Creb3l3 shows significant upregulation of transcripts (fluorescein signal) per nucleated cell (DAPI) at days 7 and 70 after surgery. Images are at 60×, scale bar = 30μm, acquired using Olympus FV1200 Confocal Microscope. Day 7: n = 15 sham images, n = 14 SBR images (3 biological replicates). Day 70: n = 9 sham images, n =10 SBR images (2 biological replicates). (B) Projection of Creb3l3 and Klf4 transcript enrichment onto the UMAP plot shows increased expression with enterocyte maturation. Color scale bar indicates relative intensity of gene expression. (C) Violin plots showing differential expression of Creb3l3 and Klf4 between sham and SBR. (D) The length down from the villus tip was measured for appreciably higher Creb3l3 (left, 3 sham and 3 SBR biological replicates) and ADA (right, 3 sham and 4 SBR biological replicates) expression, and represented as a percent of total villus length. Representative immunohistochemistry staining for ADA in sham and SBR mice at postoperative day 7 is shown. Images are at 20×, scale bar = 100 μm. (E) Relative Creb3l3 expression in SI from day 3 postoperative sham and SBR mice (n = 3 sham and n = 5 SBR mice) was measured using qPCR. (F) Interactome of genes upregulated in SBR epithelium. Genes in red are involved in RA signaling. All graphs are presented as mean ± SD. ns, not significant. *P < .05, **P < .01, ***P < .001, ****P < .0001.

Visualization of Creb3l3 expression via RNA-FISH also provided valuable information on the localization of its expression within the SI epithelium. According to a recent study, 83% of enterocyte genes demonstrate spatial zonation during homeostasis, where Creb3l3 and Klf4 transcripts were found to localize to the upper villus. This is consistent with our UMAP analysis, which showed increased Creb3l3 and Klf4 expression as enterocytes mature and migrate toward the villus tip (Figures 2A and C and 6B and C). The earlier study also reported that upper villi are collectively enriched in transcripts associated with lipoprotein biosynthetic processes, including SBR enriched transcripts Apoa4 and Apoc3. Given the relative expansion of mature proximal enterocytes as a percent total epithelium in SBR, the associated enrichment of lipid processing transcripts, and the enrichment in mature proximal enterocyte transcription factor Creb3l3 within SBR cells, we hypothesized that RNA-FISH would reveal a relative expansion of Creb3l3 expressing cells along the villus axis after SBR. Indeed, sham mice exhibited Creb3l3 expression along 50.7% ± 3.9% of their villi, vs 94.1 ± 1% in SBR (p < .01, Figure 6D), suggesting preferential expansion of enterocytes with transcriptional profiles typically found in villus tip cells, following SBR. To explore this further, we performed immunohistochemistry for an additional villus tip marker, adenosine deaminase (ADA). This confirmed increased relative concentration and distribution of ADA in the upper villi of SBR mice (22.4% ± 3.5% of upper villi with relatively intense expression in sham vs 34.7% ± 3% in SBR, P < .05) (Figure 6D). Together, these results led us to conclude that stereotypical villus axis zonation patterns are at least partially abrogated during adaptive challenge, likely underscoring the expanded mature enterocyte population in SBR.

Interactome Analysis Indicates Regional Reprogramming Is Driven by RA Signaling

Although Creb3l3 was elevated at day 7 and sustained through day 70 after SBR—suggesting a continued role in maintaining adaptation—no significant differences in Creb3l3 expression were observed at day 3 postsurgery (Figure 6E). This suggested that inductive upstream signaling at earlier stages of adaptation may be critical to driving proximalization, which is subsequently mediated and sustained by Creb3l3. To investigate this, we generated in silico interactomes from single-cell gene expression profiles of all analyzed sham and SBR cells to determine which differentially expressed genes were most strongly coexpressed, thereby inferring gene-gene relationships and pathways (Figure 6F). This approach identified a network of interacting genes induced by SBR, including Sepp1, Apoa1, Fabp2, and Rbp2, the upregulated expression, several of which we confirmed previously (Figures 4 and 5). We used a total of 59 genes from the SBR interactome to perform gene list functional enrichment analysis (5 genes from the interactome were excluded from analysis as they were absent from the database). This analysis generated a list of pathways, including “lipid digestion, mobilization, and transport” (P = 2.290 × 10–11) and “digestion of dietary carbohydrate” (P = 2.857 × 10–8), in addition to “PPAR signaling pathway” (P = 3.933 × 10–7) and “retinoid metabolism and transport” (P = 1.935 × 10–8).

In the context of this study, retinoid metabolism was of particular interest because it has been shown to play a key role in structural adaptation.21, 22, 23, 32, 33 Retinoids are derived from vitamin A, which must be obtained from the diet, and are mostly absorbed in proximal SI. RA is the intracellularly bioactive hydrolysis derivative of vitamin A, and it interacts with retinoid X receptor (RXR) and RA receptor (RAR) heterodimers, which bind to RA response elements (RAREs) within the nucleus, to drive effects. Notably, mice deficient in dietary vitamin A do not adapt after SBR, and it was observed that RA drives adaptation in part via regulation of enterocyte proliferation, migration, and apoptosis.32, 33 However, relatively little detail on the molecular changes induced by RA has been revealed so far.

To investigate how RA induces the transcriptional changes accompanying SBR, we mined our dataset for genes differentially expressed between sham and SBR treatments which were either putative targets of RA signaling based on the literature, or contain a predicted RARE as determined by FIMO Motif Search (Table 2). As a result, we found 45 genes differentially expressed between sham and SBR that likely respond to RA signaling (Table 2). Several of these genes, including Plb1 (phospholipase B1), Rbp2 (retinol binding protein 2), Apoa1 (apolipoprotein A1), Apoa4 (apolipoprotein A4), Apob (apolipoprotein B), and Apoc3 (apolipoprotein C3) (highlighted in red in Fig 6F), are active in retinoid metabolism and transport, 4 of which were are also found in our list of the top 10 most differentially expressed genes in SBR (Table 1). Importantly, Rbp2 (1.1 log fold change enriched, P < .001) is preferentially induced by RA in differentiated cells, such as the expanding enterocyte population we observe after SBR. Furthermore, motif analysis indicated that Creb3l3, our main proximal transcription factor, contains a RARE and can be activated by RXRα based on ENCODE transcription factor targets (Table 3). Klf4, an additional transcription factor identified in our dataset, is also influenced by RA signaling through RARα (Table 2). Together, these findings confirm that RA signaling is induced in cells responding to SBR, placing these signals upstream of the key transcriptional changes we observe, supporting a crucial role for RA signaling in adaptation.

Table 2

Genes Increased in Small Bowel Resection vs Sham That Are Putative Responders to RA Signaling

Apoa4a	Apoba	Gdaa	Mgama	Creb3l3a
Slc43a2a	Treha	Gpx4a	Gka	Cdhr5a
Fabp165b	Ace266b	Slc2a267b	Mme68b	Pdzk169b
Nudt470b	Plb171b	Dio167b	App72b	H2-Q267b
Klf439b	Lap370b	Ace72b	Clca4b73b	Slc5a170b
Mical170b	Fos74b	Chka75b	Rnf12876b	Fbln177b
Prap170b	H2-K167b	Vnn178b	Gls79b	Malat180b
Aqp181b	Clec2e82b	Slc6a1970b	Apol7a70b	Rbp283b
Neat184b	Ogdh85b	Egr186b	Apoa187b	Apoc388b

Predicted to contain retinoic acid response elements from Find Individual Motif Occurrences.

Putative direct and indirect retinoic acid signaling targets from the literature.

Table 3

Genes Identified in Small Bowel Resection vs Sham That Are Target Genes of Retinoid X Receptor Alpha From ENCODE Transcription Factor Targets Dataset (Human)

Fabp1	Apoc3	Sepp1	Leap2	Ephx2	Apob	Apoa1	Ace2	Slc2a2	Dnase1
Pdzk1	Creb3l3	Slc43a2	Acsl5	Pls1	Chka	Ano6	Rnf128	Fbln1	Prap1
Rfk	Khk	Malat1	Pepd	Gsdmd	Neat1	Ogdh	Dhrs1	Gpx4	Egr1
Gk	Acox1	Cdhr5	Nudt4

Source: http://amp.pharm.mssm.edu/Harmonizome/gene_set/RXRA/ENCODE+Transcription+Factor+Targets.

The results presented here so far, together with previous findings, support a model in which RA signaling induces expression of proximal SI transcriptional programs, mediated by Creb3l3 and possibly Klf4, to drive regional reprogramming as an adaptive response. As RA is a hydrolysis derivative of vitamin A, we next investigated potential RA signaling sources, and did not find differential expression of RA biogenesis genes between sham and SBR populations. This would suggest SBR epithelium does not upregulate RA biogenesis genes to increase RA signaling but, rather, may be responding to increased vitamin A exposure after SBR. This is consistent with a prior report showing epithelial RA production is dependent on substrate concentration, implicating an “active saturable enzymatic conversion process.” It is also consistent with a previous study which demonstrated decreasing concentrations of vitamin A derivatives from proximal to distal SI tissue, suggesting that a gradient of bioavailability may impact the regional effects of these compounds within the SI. Based on these observations, we propose a simple model, where after SBR, the ileum becomes exposed to dietary stimuli, including vitamin A, which previously would have been processed in the resected, more proximal, SI. We thus hypothesize that SBR disrupts the endogenous RA gradient of the small intestine, thereby promoting RARE-mediated transcription in ileal enterocytes due to novel exposure to RA.

Discussion

Here, we report epithelial single-cell analysis of ileum following proximal SBR, showing expansion of cells identified as mature proximal enterocytes in SBR vs sham mice. SBR enterocytes differentiated significantly on the basis of solute and nutrient transporters typically associated with proximal SI, especially with regard to lipid metabolism. Because the duodenum and jejunum absorb the majority of nutrients under normal conditions, we propose that this “regional reprogramming,” driven by transcriptional proximalization, is a critical component of the adaptation response to SBR. This is a novel principle, as the typically studied structural adaptation, while increasing absorptive surface area, does not necessarily facilitate functional adaptation. Rather than simple tissue hyperplasia, as others have suggested, we demonstrate that enterocyte-level alterations in transcriptional profiles occurs after SBR, in order to mimic proximal SI function. These findings highlight the significant contribution that single-cell analysis makes toward our understanding of organ pathophysiology. Of note, SBR epithelium also retained a mature distal enterocyte population (Figure 2B and 3B), and differential gene expression results suggest largely preserved ileal bile acid metabolism and cobalamin absorptive function. For example, bile acid metabolism genes Slc10a2, Nr5a2, Slc51b, Abcc3, and Nr1h4, and also cobalamin metabolism genes Lrp2 and Tcn2, showed no significant differential expression after SBR. Enterocyte basolateral bile acid transporter Slc51a was depleted in SBR, (0.29 average log fold change, P < .0001), while cobalamin receptor Cubn was enriched (0.54 average log fold change, P < .0001). These results suggest a hybrid proximal-distal identity of ileal epithelium following SBR.

While investigating causative signaling mechanisms driving these proximalization changes, we identified upregulation of 2 transcription factors associated with proximal SI: Creb3l3 and Klf4. Furthermore, upstream analysis reiterated the importance of retinoid metabolism to the adaptation response, and the depth of our analysis allowed identification of previously undescribed transcriptional changes likely mediated by RA after SBR, including Creb3l3 and Klf4. These transcription factors localized to villus enterocytes rather than crypt-based cells (Figure 6A and B). We did not identify strong evidence for perturbed intestinal stem cell regional identity based on previously identified markers, consistent with previous observations that retinoic acid exerts its effects on differentiated cells. We therefore conclude that a critical component of the adaptation response to massive proximal SBR is transcriptomic “proximalization” of distal SI enterocytes, and that this is driven at least in part by RA signaling which is upstream of “proximalization” transcription factors and signaling cascades.

RA is derived entirely from the diet. Following proximal SBR, the ileum becomes exposed to nutrients in the luminal content which would otherwise have been largely absorbed more proximally, including vitamin A. A series of studies utilizing an ileostomy model of SGS in mice and zebrafish have shown the critical effects of mechanoluminal flow on the structural adaptation process after SBR, including loss of structural adaptation and cellular proliferation in the distal bowel when isolated from the flow of luminal contents.43, 44 The importance of luminal nutrition or enteral feeding has also been highlighted by other groups.45, 46 As such, exposure to increased levels of dietary vitamin A is a likely mechanism driving structural adaptation, consistent with prior reports.21, 22, 32, 33 In the current study, we have identified novel regulatory networks and transcriptional changes downstream of RA signaling, which likely drive regional reprogramming as well.

In addition to vitamin A, the ileum also becomes exposed to higher levels of dietary fatty acids after SBR, and this likely constitutes an additional driving force for adaptation. In line with this, our lab has previously demonstrated that a high fat diet enhances villus growth following SBR, but not sham surgery, though enhanced structural adaptation did not correlate with enhanced functional adaptation (weight gain), despite increases in fatty acid transporters such as CD36. At the same time, this study did not examine the transcripts found to be most significant in the current dataset (Table 1), and it may be that ligands designed to induce those specific transcripts could have a different effect. Regardless, these findings again demonstrate the potential for incongruence between structural and functional adaptation, indicating the adaptation process is inherently multifaceted, and different aspects likely rely on different stimuli. For example, 2 of the upregulated transcripts in SBR—Fabp1 and Fabp2 (Table 1)—are regulated independently of one another, with Fabp1 activated by peroxisome proliferator-activated receptor alpha (or PPARα), in response to dietary fatty acids, and Fabp2 hormonally by PYY.

The PPAR signaling pathway was implicated by our interactome analysis. We investigated this further since PPARs form heterodimers with RXRs to activate PPAR response elements (PPREs) in the induction regions of many genes involved in lipid metabolism,49, 50 including Creb3l3. Interestingly, it is thought the ratio of RA binding protein to fatty acid gene expression determines the functional outcome of RA signaling, and further, treatment of mice with a PPARα agonist induced villus growth by facilitating cell differentiation, similar to the adaptation phenotype observed after SBR in which there are elongated villi with a preponderance of mature enterocytes. Ultimately, we observed neither changes in Pparα expression (1.5% decrease 3 days after SBR via qPCR, P = .96, not shown), nor significant changes in Pparδ expression (55.2% decrease 3 days after SBR, P = .16, not shown). Less is known about PPARδ, though it is thought to interact with corepressors and function as an inhibitor of PPARα. Given these findings, we suspect that either (1) PPAR signaling is important to SBR adaptation, but was not captured in our analysis temporally, or (2) minimal to no transcriptional change in these specific genes is needed to drive a significant biological effect.

Finally, another key finding from this study was the abrogation of villus zonation patterns during adaptive challenge, with the expansion of villus tip transcript Creb3l3 and ADA along the lengths of villi. The upper villus is typically responsible for fatty acid absorption, which is less metabolically demanding, and so preferential fat absorption/chylomicron secretion at the upper villus parallels the decreased bioavailability of oxygen in this area. Indeed, we have previously demonstrated that adaptation to SBR is associated with relative hypoxia, and so the metabolic incentive to prioritize fatty acid absorption is at least 2-fold after SBR: starvation and relative hypoxia. Mechanisms underlying redistribution of villus zonation during physiologic challenge warrant further investigation.

In summary, our analysis has revealed a significant shift in metabolic machinery and regional identity at the enterocyte level following SBR. Moving forward, additional studies are warranted to better delineate causal factors driving changes between sham and SBR enterocytes, in conjunction with or independent of RA signaling. This is especially true considering prior studies that demonstrated proliferative and morphometric effects of circulating factors on structural “jejunalization” of ileum following SBR,56, 57 which implicates nonluminal stimuli in driving structural adaptation. Similar studies exploring molecular changes in response to circulating factors would provide further insight. Discerning the stimuli for functional proximalization of ileum following SBR will prove critical, as it will provide insight toward targeted therapeutic approaches, via the enteral or parenteral route, for patients suffering from SGS. Targeted therapeutic approaches combining aspects of adaptation from a structural and molecular avenue could induce heightened adaptive responses, yielding better patient outcomes.

Conclusions

Here, we have characterized the transcriptome of adapted intestinal epithelium at the single-cell level following massive SBR, a laboratory model for SGS, using scRNA-seq. Our analysis revealed the emergence of unique enterocyte gene expression patterns between sham and SBR mice, which distinguished themselves on the basis of proximal vs distal SI patterning and cell identity, including critical absorptive and anti-inflammatory/immunological features. Pathways driving these changes, such as RA signaling, deserve further investigation, as they underlie the functional aspects of adaptation to SGS, which allow progressive tolerance of enteral feeding and weaning from PN.

Materials and Methods

Mice

A total of 50% proximal SBR was performed on male C57/B6 mice at 8–12 weeks of age, according to standard protocol (Table 4). Briefly, the SI was extruded via a midline laparotomy and the ileocecal valve identified. The SI was transected 12 cm proximal to the ileocecal valve, and ∼2 cm distal to the ligament of Treitz. The intervening SI was removed, the mesentery ligated with 3-0 silk suture, and the proximal and distal ends approximated with interrupted 9-0 nylon stitches. Sham surgery, consisting of distal transection and anastomosis only, was performed as control. Peritoneum and skin were approximated in separate layers, animals were resuscitated with a subcutaneous bolus of normal saline (repeated on postoperative day 1) and co-housed in a 33°C incubator until the end of the 7 day experiment. For longer studies (70 days), mice were moved to room temperature at day 7. Liquid diet (PMI Micro-Stabilized Rodent Liquid Diet LD 101; TestDiet, St. Louis MO) was initiated 24 hours before surgery, withheld the morning of surgery, and subsequently provided on postoperative day 1 until the end of the experiment. Food and water were available ad libitum, and animals were housed under 12-hour light/dark cycles with corn cob bedding and nestlet enrichment. All surgical and animal care procedures were approved by the Washington University Institutional Animal Care and Use Committee, and meet Animal Research: Reporting of In Vivo Experiments standards.

Table 4

Materials

Antibody list
Antibody	Species	Company	Catalog No	Lot No	Dilution
Primary
Selenoprotein P (SEPP1)	rabbit	Abcam (Cambridge, MA)	ab193193	GR192155-11	1:200 IHC 1:100 IF
Fatty acid binding protein 6 (FABP6)	rabbit	Abcam	ab91184	GR3210318	1:50 IHC 10 μg/mL IF
Fatty acid binding protein 1 (FABP1)	rabbit	Abcam	ab222517	GR3184395-2	1:3000 IHC 1:50 IF
Adenosine deaminase (ADA)	rabbit	Abcam	ab175310	GR3259503-1	1:100 IHC
Apolipoprotein C3 (APOC3)	rabbit	Abcam	ab55984	GR137329-12	1:1000 WB
Glyceraldehyde 3-phosphate dehydrogenase (GAPDH)	rabbit	Cell Signaling Technologies (Danvers, MA)	5174S	6	1:10000 WB
Mucin 2 (MUC2)	rabbit	Santa Cruz Biotechnology (Dallas, TX)	sc-15334	D0815	1:200 IHC
Secondary
Alexa Fluor anti-rabbit 555	goat	Invitrogen	A21428	1858260	1:500 IF
Horseradish peroxidase (HRP)	goat	Cell Signaling Technologies	7074S	27	1:2500 WB
Biotin-SP-conjugated AffiniPure Anti-Rabbit IgG	goat	Jackson ImmunoResearch Laboratories, Inc (West Grove, PA)	111-065-045	128361	1:400 IHC
Streptavidin-HRP		Jackson ImmunoResearch Laboratories, Inc	016-030-084	129316	1:400 IHC

IF, immunofluorescence; IHC, immunohistochemistry; WB, Western blot

Tissue Isolation and Processing

At day 7 after surgery, epithelium was isolated from a 1-cm segment of SI 3 cm distal to the anastomosis, similar to previously published protocols.58, 59 Briefly, the SI was flushed with ice cold sterile saline, filleted lengthwise, and placed in a conical tube containing ice cold 30 mM EDTA in phosphate-buffered saline (PBS). After 15 minutes on ice without agitation, the SI segment was transferred to a fresh conical tube containing 30-mM EDTA in PBS, briefly shaken, and placed in a 37°C water bath for 15 minutes. Subsequently, the tube was shaken aggressively by hand for 2 minutes. Subepithelial tissue floated to the top and was removed, and epithelium was pelleted by centrifugation. Epithelium was then re-suspended in a 0.3-U/mL dispase solution (07923; Stem Cell Technologies; Cambridge, MA) and incubated at 37°C for 15 minutes, shaking every 2 minutes. The solution was then quenched with media containing fetal bovine serum to a final concentration of 5%, pipetted several times, and sequentially passed through 100-, 70-, and 40-μm filters. Single-cell suspensions were confirmed by microscopy, pelleted by centrifugation, and resuspended in 200-μL ice cold PBS. Then, 800-μL ice cold 100% methanol was added, dropwise, with gentle mixing between drops. Samples were immediately stored in 80% methanol in PBS at –80°C, according to Alles et al, for later processing. For Western blotting analysis, whole epithelium was isolated using 30-mM EDTA, pelleted by centrifugation, and lysed in sodium dodecyl sulfate sample buffer (50-mmol/L Tris-HCL, pH 6.8, 2% sodium dodecyl sulfate, 10% glycerol, and 5% 2-mercaptoethanol). Lysate was heated to 100°C and stored at –20°C prior to processing. Protein concentration was measured using the RC DC (reducing agent and detergent compatible) Protein Assay Kit II (5000122; Bio-Rad; Hercules, CA). For qPCR experiments, RNA was isolated from homogenized whole SI using the standard Trizol method. RNA concentration was measured using a NanoDrop Spectrophotometer (ND-1000; NanoDrop Technologies, Wilmington, DE). A total of 1-μg RNA was converted to complementary DNA (cDNA) using qScript cDNA Synthesis Kit (95047, Quanta Bio, Beverly, MA), according to the manufacturer’s instructions, and stored at –20°C until use.

scRNA-seq Library Preparation

For single-cell library preparation on the 10x Genomics Chromium platform, we used the Chromium Single 3′ Library & Gel Bead Kit v2 (PN-120237), Chromium Single Cell 3′ Chip kit v2 (PN-120236), and Chromium i7 Multiplex Kit (PN-120262), according to the manufacturer’s instructions in the Chromium Single Cell 3′ Reagents Kits V2 User Guide. Methanol-fixed cells from sham (n = 3) and SBR (n = 3) animals were pooled for the first batch. Methanol-fixed cells from sham (n = 2) vs SBR (n = 1) were processed individually in a separate experiment for a final sample size of sham (n = 5) and SBR (n = 4). Just before cell capture, methanol-fixed cells were placed on ice, then spun at 3000 rpm for 5 minutes at 4°C, followed by resuspension and rehydration in PBS, as previously described. Resulting cDNA libraries were quantified on an Agilent Tapestation and sequenced on an Illumina HiSeq 2500.

scRNA-seq Analysis

The Cell Ranger v2.1.0 pipeline was used to align reads to the mm10 genome build, and generate a digital gene expression (DGE) matrix: (https://.support.10xgenomics.com/single-cell-gene-expression/software/downloads/latest). For initial filtering of these DGE matrices, we first excluded cells with a low number (<200) of unique detected genes. We then excluded cells for which the total number of unique molecules (UMIs) (after log10 transformation) was not within 3 standard deviations of the mean. This was followed by the exclusion of outlying cells with an unusually high or low number of UMIs/genes given their number of reads by fitting a Loess curve (span = 0.5, degree = 2) to the number of UMIs/genes with number of reads as predictor (after log10 transformation), removing cells with a residual more than 3 standard deviations the mean. Finally, we excluded cells in which the proportion of the UMI count attributable to mitochondrial genes was >25%. Raw and processed data files are available via GEO: accession number GSE130113. After filtering and normalization of the DGE, the R package Seurat (Version 3) was used to cluster and analyze the single-cell transcriptomes. Independent biological replicates from the sham and SBR surgeries were integrated by Canonical Correlation Analysis, identifying common sources of variation to align the datasets, reducing batch effects. Highly variable genes were identified and used as input for dimensionality reduction via canonical correlation analysis. The resulting Canonical Correlation Vectors and the correlated genes were examined to determine the number of components to include in downstream analysis, followed by clustering and visualization via UMAP.

Quadratic Programming Analysis to Assess Cell Identity and State

QP, previously described in Treutlein et al and successfully modified and used by our group in Biddy et al, was used to score cell identity. Here, we created a reference of SI epithelial cell types, collected previously. The R Package QuadProg was used for QP to generate cell identity scores, modifying our earlier approach by modeling cell type classification as a multivariate linear regression, solved for fractional cell types or identities. This approach enables cell identities to be determined. In addition, it enables more subtle changes in cell identity to be quantified.

Immunohistochemistry

Ileal tissue adjacent to the region collected for single-cell preparation was fixed in 10% neutral buffered formalin, paraffin embedded, and sectioned at a thickness of 5 μm. Deparaffinization and immunolocalization were performed as previously described. Briefly, slides were deparaffinized in xylene, rehydrated in sequential ethanol baths, and prepared in 3% hydrogen peroxide in methanol. Antigen retrieval was performed using 1× Diva Decloaking Solution (DV2004; Biocare Medical, Pacheco, CA), and blocking was performed using the Avidin-Biotin kit (AB972L; Biocare Medical). Primary antibodies were diluted in Da Vinci green (PD900L; Biocare Medical) and incubated overnight at 4°C. Slides were rinsed in phosphate-buffered saline + TWEEN 20 (PBST), incubated in biotin-labeled secondary IgG diluted in PBST, rinsed in PBST, incubated in streptavidin-horseradish peroxidase diluted in PBST, developed in DAB (D9015; Sigma, St. Louis, MO), counterstained with hematoxylin and bluing agent, run in successive dilutions of ethanol, and xylene, and cover-slipped using MM 24 mounting medium (100109; Surgipath, Richmond, IL). Of note, samples used for confirmatory staining were from a different litter of mice than those used for scRNA-seq analysis. This was done to validate consistency of results across cage and littermates, which has been previously reported as a confounding variable in murine gastrointestinal research. At least 3 sham and 3 SBR samples were analyzed after surgery; 20× images representative of the sample were obtained by a blinded investigator using a Nikon Eclipse 80i with Ds-Ri2 camera and NIS Elements V4.3 software (Nikon Instruments, Inc, Melville, NY).

Western Blotting

A total of 20 μg of each sample and 15-μL Novex Sharp Pre-stained Protein Standard (57318; Invitrogen, Carlsbad, CA) was loaded onto a 18% polyacrylamide gel. Western blotting was performed on a nitrocellulose membrane (IB301001; Invitrogen) after a dry transfer using the iBLOT Gel Transfer Device (IB1001; Invitrogen). The membrane was blocked with 5% bovine serum albumin in PBST for 1 hour at room temperature, followed by overnight incubation in primary antibody at 4°C. The membrane was washed in PBST 3 times in 10-minute increments, incubated for 1 hour in secondary antibody at room temperature, washed in PBST, and developed using GE Healthcare Amersham ECL Western Blotting Detection Reagents (6883S; GE Healthcare, Chicago, IL). Image Lab Software (Bio-Rad) was used to detect and quantify proteins.

Quantitative PCR

cDNA was amplified using TaqMan Gene Expression Master Mix (4369016; Applied Biosystems, Foster City, CA) and the specified primer probe on the Applied Biosystems 7500 Fast Real-Time PCR system. Primer probes were Actb (Mm02619580_g1, endogenous control), Creb3l3 (Mm00520279_m1), Ephx2 (Mm01313813_m1), and Rbp2 (Mm00436300_m1), all from Applied Biosystems.

Immunofluorescence

Ileal tissue adjacent to the region collected for single cell-preparation was fixed overnight in 4% paraformaldehyde and then overnight in 30% sucrose before being embedded in O.C.T. Compound (23-730-571; Fischer Healthcare, Houston, TX), sectioned at 5 μm, and stored at –80°C until use. Slides were washed in PBS and blocked in 5% goat serum with 0.3% Triton X-100 (T8787; Sigma) prior to incubation with primary antibody in antibody staining solution (1% goat serum and 0.3% Triton X-100) at 4°C overnight. Slides were rinsed in PBS for 5 minutes 3 times, and secondary antibody was applied in antibody staining solution for 1 hour at room temperature. A total of 300-nM DAPI (D1306; Invitrogen, Lot 1802085) was applied for 1 minute, and slides were rinsed in PBS for 5 minutes 3 times before a cover slip was applied using ProLong Gold Antifade Mountant (P10144; Invitrogen). The 40× images representative of the sample were obtained by a blinded investigator using a Nikon Eclipse 80i with Ds-Ri2 camera and NIS Elements V4.3 software (Nikon Instruments, Inc). Images were subsequently analyzed using the FIJI Distribution of ImageJ for FABP6 and SEPP1, or unbiased computational quantification (FABP1; see “Quantitative analysis of RNA-FISH images”). Representative images were chosen.

ImageJ Quantitative Analysis

Briefly, for SEPP1, regions of interest were drawn tracing the edges of all villi in the image. Each image was thresholded to the same intensity and the percent area within regions of interest expressing fluorescence was collected for each technical replicate. For FABP6, the full area of each villus was chosen as a region of interest. The mean intensity of that region was calculated to represent the fluorescence of the diffuse FABP6 staining. For percent villus occupancy, ImageJ was used to measure the areas of Creb3l3 expression or ADA expression, and this was divided by the total villus length. For Creb3l3, there were 12 sham villi (3 biological replicates) and 11 SBR villi (3 biological replicates). For ADA, there were 12 sham villi (3 biological replicates) and 16 SBR villi (4 biological replicates).

RNA Fluorescent In Situ Hybridization

Tissues were fixed and sectioned as previously described in methods for immunohistochemistry and immunofluorescence. RNA-FISH was performed using the RNAscope Multiplex Fluorescent v2 kit (323100; Advanced Cell Diagnostics, Newark, CA), following the protocol for Fixed Frozen Tissue. Briefly, tissue was pretreated with RNAscope Hydrogen Peroxide (322335; Advanced Cell Diagnostics), target retrieval was performed for 5 minutes, and tissue was treated with RNAscope Protease III (322337; Advanced Cell Diagnostics). Then, specified probes were hybridized using the RNAscope HybEZ II Oven (321710/321720; Advanced Cell Diagnostics). Probes were then amplified and the HRP signal was developed using TSA Fluorescein Plus Evaluation Kit (NEL741E001KT; Perkin Elmer, Waltham, MA). Finally, slides were counterstained with DAPI (323108; Advanced Cell Diagnostics) and mounted with ProLong Gold antifade reagent (P10144; Life Technologies Corporation, Eugene, OR). Imaging was performed on an Olympus FV1200 Confocal Microscope, in which multiple 60× images were taken per sample, of which 2–3 samples were used for each timepoint (day 7 vs day 70) in each treatment (sham vs SBR). Images were subsequently analyzed using unbiased computational quantification and representative images were chosen.

Quantitative Analysis of RNA-FISH Images

RNA-FISH images were processed with a custom python script to quantify gene expression level at single-cell resolution. Individual cell segmentation was achieved based on the nuclear signal. First, DAPI images were transformed into binary images by thresholding the fluorescent signal. The threshold values were determined by the Otsu method. Binarized nuclei images were processed by the watershed segmentation method to completely separate individual objects. The images were subjected to 2-step quality check: filtering of objects and filtering of images. First, inappropriately sized objects were filtered to remove noise and cell multiplets. Then, images with a large number of inappropriate objects were removed. The intensity of the fluorescent signal per individual cell area was then quantified. Fluorescent signals per image were averaged to obtain mean signals per sample and treatment.

Gene Coexpression and Interactome Analysis

To further reveal relationships among differentially expressed genes following SBR, we constructed gene coexpression networks using weighted gene correlation network analysis, adapted for single-cell analysis as in https://hms-dbmi.github.io/scw/WGCNA.html. The analysis was performed using the R package WGCNA. In brief, differentially expressed genes between SBR and sham samples were identified via Seurat analysis, as previous (using the function FindMarkers). Genes that were upregulated in SBR treatment were used for network construction. Correlations of each gene pair among all significantly differentially expressed genes were then calculated. A gene-gene correlation matrix was used to construct an adjacency matrix by raising the correlations to a soft-threshold power, from which Topological Overlap Matrix (TOM) was further computed to remove spurious correlations. With the TOM matrix, the algorithm identifies modules/clusters of genes via clustering using Ward’s method and Dynamic Branch Cut methods. Here, to identify the most significant connections, we select the top 5% of the most differentially expressed genes for visualization. The network was visualized using CytoScape based on the TOM matrix.

Statistical Analysis

Statistical comparison of the QP-generated identity scores between the 2 groups, sham and SBR, was performed using an unpaired Student’s t test. For Western blot, qPCR, and percent villus occupancy studies, statistical analyses were performed using Prism 7.0 (GraphPad Software, San Diego, CA). Differences in protein expression and mRNA expression between sham and SBR were assessed using unpaired Student’s t tests. Graphs with error bars represent mean ± SD. P < .05 was considered significant. For immunofluorescence quantitative analysis using FIJI, SEPP1 and FABP6 replicate values were averaged between sham and SBR treatment identities. A Wilcoxon rank sum test was performed to determine statistical significance. For quantitative analysis of RNA-FISH images, and FABP1 immunofluorescence, statistics were calculated using the Wilcoxon rank sum test. These processes were run with Python 3.6.7 and its libraries: scikit-image 0.13.1, numpy 1.14.3, pandas 0.24.1, oiffile 2019.1.1, matplotlib 3.0.3, seaborn 0.8.1, and jupyter 1.0.0.