JAK/STAT signaling promotes the emergence of unique cell states in ulcerative colitis

Abstract

The intestinal epithelium ensures uptake of vital nutrients and acts as a barrier between luminal contents and the underlying immune system. In inflammatory bowel diseases, such as ulcerative colitis (UC), this barrier is compromised, and patients experience debilitating symptoms. Here, we perform single-cell RNA profiling of epithelial cells and outline patterns of cell fate decisions in healthy individuals and UC patients. We demonstrate that patterns of hierarchical behavior are altered in UC patients and identify unique cellular states associated with Janus kinase/signal transducer and activator of transcription (JAK/STAT) activation in ulcerated and non-ulcerated areas of the colonic epithelium. These transcriptional changes could be recapitulated in human colonic organoids, wherein cytokine-mediated activation of JAK/STAT led to the emergence of cell populations with augmented regenerative properties. Altogether, our findings indicate that intricate relationships between epithelial and cytokine signaling regulate cell fate during epithelial tissue regeneration in humans and have important implications for the understanding of UC biology.

Figure 1

Characterization of healthy human colonic epithelium samples (A) Representative images of live colonoscopy (top) and histological analysis (bottom) confirming ulceration in the epithelium of UC patients (bottom). Histological samples were stained for the gastrointestinal epithelium marker KRT20 and counterstained with hematoxylin. Scale bar: 100 μm. (B) UMAP plot of 11,035 colonic epithelial cells profiled from healthy samples following annotation based on expression of cell lineage markers. (C) UMAP plots highlighting expression of known lineage markers in the cell type clusters identified in (B). (D) Cell type trajectory graph based on over-clustering. Nodes represent clusters colored by cell type of origin, while connecting edges indicate the similarity degrees. (E) RNA velocity measurements projected onto the UMAP embedding shown in (B). (F) Immunostaining for TFF3 (blue), PCNA (red), and β-Catenin (green) in the healthy colonic epithelium. Arrow indicates a TFF3/PCNA double-positive cell. Scale bar: 100 μm. (G–I) Spatial distribution of RNA molecules for MUC2 (red), HMGB2 (blue), KRT20 (orange), and CEACAM1 (green). DNA was counterstained with DAPI (gray). (G) and (H) show the entire crypts length. Scale bar: 50 μm. (I) shows a zoomed-in region of interest of crypt-top area. CT colonocyte (arrow) and goblet cell (arrowhead) are indicated. Scale bar: 10 μm (H).

Figure 2

Altered cell fate patterns in the ulcerated epithelium samples (A) Immunostaining for CEACAM1 (green), PCNA (red), and β-catenin (yellow) human tissues. Scale bar: 100 μm. (B) UMAP depicting 4,351 colonic epithelial cells from patients with UC following cell type annotation as in (1B). (C and D) Trajectory graph based on over-clustering and RNA velocity analysis of epithelial cells from ulcerated samples. Arrows (D) point to velocity vectors going from differentiated cell clusters toward progenitor clusters. (E) (Left) representative image of immunohistochemistry staining for KRT20 and KI67 in healthy and ulcerated human colonic samples. Arrows indicate KRT20/KI67 double-positive cells. (Right) bar plots showing percentage of KI67 single-positive and KRT20/KI67 double-positive cells in the indicate areas. Each dot represents an individual donor (average of at least 5 crypts each). H, healthy controls; UC, ulcerated samples. ^∗ for p < 0.05 by Student’s t test. Scale bar: 100 μm. (F) Density plots showing of KRT20 expression (y axis) and proliferation score (x axis) in the healthy (left) and ulcerated (right) datasets. (G) Representative flow cytometry analysis of KRT20-tdTomato human organoid reporter line. A wild-type line was used as a control for tdTomato levels. Schematic of the reporter gene is shown above. (H) Quantitative reverse-transcription PCR (RT-qPCR) analysis of the indicated genes in tdTomato⁻ and tdTomato⁺ sorted from KRT20-tdTomato reporter line. Expression levels of genes of interest were normalized to GAPDH. Bars show fold change relative to KRT20-tdTomato^- cells ±SEM; each dot represents an independent experiment. ^∗∗∗ for p < 0.005 by ordinary 2-way AVOVA with Sidak’s correction for multiple comparison. (I) (Left) representative phase-contrast images of human colon organoids formed from sorted tdTomato⁻ and tdTomato⁺ (bottom) cells 10 days after seeding. Scale bars: 275 μm. (Right) quantification of organoid-forming efficiency and size of organoids formed from purified KRT20-tdTomato^- and tdTomato⁺ cells. For efficiency, plot shows mean ± SEM; each dot represents an independent experiment (average of four replica per experiment). For size, plot was generated using the Tukey method. The number of organoids analyzed within one representative experiment is shown below the bars.

Figure 3

Emergence of new cell states in the ulcerated epithelium (A) UMAP of epithelial cells from the combined datasets from healthy individuals, and patient-derived healthy margin and ulcerated samples. Colors indicate different cell cluster types. Arrows indicate the inflammation-associated cell clusters. (B) UMAP of the combined datasets from healthy, healthy margin, and ulcerated epithelial cells colored by cell type and dataset of origin. (C) Bar plots showing cell fractions within the indicated lineages in the healthy, healthy margin, and ulcerated datasets, including the IA cell states. Colonocytes, BEST4+, and goblet cell lineages include TA and CT populations. Significance tested with two-sided Mann-Whitney U test. Each dot represents an individual donor. Bars represent the mean ± SEM. UC-HM, UC patient-derived healthy margin; UC-UL, patient-derived ulcerated samples. (D) (Left) Venn diagrams indicating the number of shared upregulated and downregulated genes between the stem cell clusters from the ulcerated and healthy margin samples (when compared to the healthy state). (Right) bar plots showing top 10 terms from GO term overrepresentation analyses of the shared differentially expressed genes. Numbers in each bar indicate the fraction of genes within that GO term found in the shared differentially expressed genes (DEGs). (E) Analysis of CD74 expression in the indicated samples. (Top) UMAP plots with heatmaps of normalized expression of the CD74 gene in the integrated dataset split by sample type of origin. (Bottom) Representative immunohistochemistry images for CD74 in tissue biopsies. Scale bars: 100 μm. (F) Similarity correlation matrix between healthy, ulcerated, and healthy margin samples. The last was split into CD74-positive and negative fractions (cutoff value based on the third quartile of CD74 expression in the healthy dataset). Similarity among samples calculated based on first 50 principal components using Pearson correlation. Dendrogram lines on the left shows hierarchical clustering of the sample types. (G) Trajectory graph of cell type clusters in the integrated datasets. Pie charts indicate proportions of cells in each cell-cycle phase. Outer ring colors match cell type clusters in (A). (H and I) Violin plots showing normalized enrichment score of the indicated signatures and expression levels of the indicated lineage-specific genes (CA1 and MUC2) in the stem cell, IA colonocyte, and IA goblet clusters from the ulcerated samples. Significance was tested with two-sided Mann-Whitney U test.

Figure 4

Inflammation-associated cell states are associated with JAK/STAT signaling (A) Heatmap of scaled transcription factor activity from the annotated DoRothEA regulons. Top three differentially active regulons for each of the cell type clusters are depicted. (B) Heatmap of scaled pathway activity for ulcerated cell populations. Activity scores have been calculated using genes differentially expressed in the ulcerated dataset, as compared to the healthy dataset. The DEGs have been found between pseudo-bulk profiles of cell types from the two datasets when at least two biological replicates per condition were available. In the ulcerated dataset, IA Colonocytes and IA goblet cells have been combined with TA colonocytes and TA goblet cells, respectively. (C) UMAP plot showing enrichment for STAT3 transcriptional regulon in the integrated dataset split by sample type. Arrows indicate location of IA clusters in the ulcerated samples. (D) Venn diagram of the intersection between shared upregulated genes in IA colonocytes and IA goblet cells (compared to their respective lineages) in the ulcerated dataset. (E and F) UMAPs with a heatmap of the normalized thirteen IA signature genes (E) and REG1A expression in as in (C). (G) Representative immunofluorescent images for REG1A (green), TFF3 (red), and b-catenin or E-Cadherin (white) in healthy and ulcerated tissue samples. Arrowheads indicate single REG1A-positive cells. Arrows indicate TFF3/REG1A double-positive cells. Scale bars: 100 μm.

Figure 5

IL-22 activates JAK/STAT in human colonic epithelium organoids (A) Bright-field images of control and IL-22-treated human colonic organoids for 3 days. Scale bars: 200 μm. (B) Bar plot showing reduction of tdTomato+ cells (KRT20-tdTomato reporter line) following treatment with IL-22. Each dot represents an independent experiment (average of 3–6 replica wells per experiment). ^∗p < 0.05 by paired Student’s t test. (C) RT-qPCR analysis showing relative expression levels of indicated genes following IL-22 treatment. Plot shows fold change relative to untreated organoids ± SEM. Each dot represents an independent experiment. ^∗p < 0.5 by two-way ANOVA. (D) Representative images of human colonic organoids with IL-22 and/or JAK inhibitor tofacitinib for 3 days. Scale bars, 200 μm. (E) RT-qPCR analysis of showing relative expression of the IA-associated genes following treatment with IL-22 and/or tofacitinib as in (D). Plot depicts fold change relative to untreated control organoids ±SEM; each dot represents an independent experiment. ^∗∗∗∗ for p < 0.001 by paired two-way ANOVA followed by Dunnett’s multiple comparison test. (F) Schematics of organoid formation assay setup. Bottom part shows representative images of organoids at 168 h formed following treatment with with IL-22 at the indicated concentrations during the initial 72 h. (G) Time course RT-qPCR analysis confirming upregulation of REG1A upon IL-22 treatment. Pink shade indicates exposure window to IL-22. Each dot represents technical replicas of a single representative experiment. (H) Quantification of organoid forming efficiency and organoid size following IL-22 treatment as indicated in (F). For efficiency, plot shows mean ± SD. Each dot represents a replica well of three independent experiments. For size, plot was generated using the Tukey method. The number of organoids analyzed within one representative experiment is shown below the bars. a.u, arbitrary unit. ^∗∗, ^∗∗∗, and ^∗∗∗ for p < 0.01, 0.005, and 0.001, respectively, by two-way ANOVA.

Figure 6

IL-22 treatment recapitulates the emergence of IA cellular states in human colonic organoids (A) UMAP plots of control and IL-22-treated organoid datasets colored by condition (left) and cell type annotation (right). (B and C) UMAP plots showing enrichment for STAT3 transcriptional regulon (C) and the shared IA signature (C) in control and IL-22-treated organoids. (D) UMAP plots with a heatmap of normalized expression REG1A in control and IL-22-treated organoids. (E and F) Violin plots showing normalized enrichment score of the indicated signatures and expression levels of the indicated lineage-specific genes (CA2 and ZG16) in the cell type clusters identified in organoids treated with IL-22. (G) Representative flow cytometry analysis of REG1A-mNeon human organoid reporter line (clone#3) following IL-22 treatment for 3 days. Simplified schematic of the reporter is shown above. (H) Quantification of REG1A-expressing cells in human colonic organoids treated with IL-22. Bar plots showing percentage of mNeon-positive cells (mean ± SD) in the REG1A-mNeon reporter line (clone#3; left) and percentage of high REG1A cells (normalized expression in the third quantile) in the scRNA-seq dataset (right). Each dot represents an independent experiment. ^∗∗ by paired Student’s t test. (I) Representative images of REG1A-mNeon reporter line (clone#3) treated with IL-22 and/or tofacitinib for 3 days. Images show merged phase contrast and GFP channels. The secreted REG1A-mNeon fusion protein accumulates in the lumen of the organoids. Scale bars: 200 μm. (J) RT-qPCR analysis showing relative expression levels of the indicated genes in the fluorescence-activated cell sorting (FACS)-purified cell populations from the REG1A-mNeon reporter line (clone#3). Fold change relative to mNeon− cells ±SEM. Each dot represents an independent experiment. ^∗p < 0.05, ^∗∗p < 0.01, ^∗∗∗p < 0.001, ^∗∗∗∗p < 0.0001; two-way ANOVA with Tukey’s multiple comparison test. (K) Organoid-forming assay of purified mNeon-negative and mNeon-positive cells from IL-22-treated organoids (clone#3). mNeon-negative cells from untreated organoids were used as control. Left shows representative images of organoids formed from each purified cell fraction. Scale bar: 300 μm. Right shows quantification of formation efficiency and organoid size, as described in (2I). Efficiency plot shows mean ± SEM. Each dot represents an independent experiment. Size is shown as Tukey plot in log2 scale. The number of organoids analyzed within one representative experiment is shown below the bars. ^∗∗ for p < 0.01 by one-way ANOVA.

Figure 7

Effect of tofacitinib (JAK inhibitor) in mouse DSS-induced colitis (A) Time course analysis of body weight changes of mice in all experimental groups. Gray shade (from day 0 to day 5) indicates presence of DSS in drinking water. Arrow heads indicate doses of vehicle or tofacitinib. Individual dots represent the relative weigh of each animal compared to day 0. ^∗∗ for p < 0.01 at day 10 by two-way ANOVA followed by Šídák’s multiple comparisons test. (B) Swiss roll sections of colon from control and DSS groups as indicated in (A). Sections were stained with hematoxylin & eosin (H&E). Scale bar: 600 μm. (C) Assessment of tissue damage in the histological preparations. Left plot shows length of epithelial erosions from H&E staining. Right plot shows total length of the colon for comparison. Individual dots represent the values for individual animals at day 10.