Pediatric Crohn disease patients exhibit specific ileal transcriptome and microbiome signature

Abstract

Interactions between the host and gut microbial community likely contribute to Crohn disease (CD) pathogenesis; however, direct evidence for these interactions at the onset of disease is lacking. Here, we characterized the global pattern of ileal gene expression and the ileal microbial community in 359 treatment-naive pediatric patients with CD, patients with ulcerative colitis (UC), and control individuals. We identified core gene expression profiles and microbial communities in the affected CD ilea that are preserved in the unaffected ilea of patients with colon-only CD but not present in those with UC or control individuals; therefore, this signature is specific to CD and independent of clinical inflammation. An abnormal increase of antimicrobial dual oxidase (DUOX2) expression was detected in association with an expansion of Proteobacteria in both UC and CD, while expression of lipoprotein APOA1 gene was downregulated and associated with CD-specific alterations in Firmicutes. The increased DUOX2 and decreased APOA1 gene expression signature favored oxidative stress and Th1 polarization and was maximally altered in patients with more severe mucosal injury. A regression model that included APOA1 gene expression and microbial abundance more accurately predicted month 6 steroid-free remission than a model using clinical factors alone. These CD-specific host and microbe profiles identify the ileum as the primary inductive site for all forms of CD and may direct prognostic and therapeutic approaches.

Figure 7. Host gene expression and microbial shifts across the spectrum of ileal IBD.

Progressive induction of an ileal DUOX2 host gene coexpression signature in association with expansion of Proteobacteria taxa was observed in UC, cCD, and iCD, with the greatest change relative to Ctls detected in iCD-DU (far right). By comparison, alteration of an ileal APOA1 host gene coexpression signature in association with reduction in Firmicutes taxa was specific to all forms of CD, whereby the majority of the molecular signature of iCD was present in cCD and hence was largely independent of the degree of local inflammation. This emphasizes the central role of the ileum in the pathogenesis of all forms of CD. Maximal overall alteration of these microbial shifts and host responses favored oxidative stress and Th1 polarization and was associated with the most severe tissue injury manifested as ileal deep ulcers (far right).

Figure 6. Covariation of the ileal microbial community structure with ileal gene expression and clinical subgroup and severity.

(A) The biplot depicts covariation of the ileal microbial community structure with clinical group, clinical disease activity (PCDAI), and ileal gene expression using NMDS rotated by CD and based on significant associations obtained from MaAsLin. Patients with IBD and healthy Ctls are plotted as orange circles (Ctl), green circles (UC), red circles (cCD), purple circles (iCD-noDU), or purple triangles (iCD-DU). Arrows indicate the direction of covariation for APOA1 or DUOXA2 gene expression and clinical severity (PCDAI). Microbial taxa are labeled with capitalized lettering. The stress measurement tests for the goodness of fit of the biplot 2-dimensional depiction of the multidimensional data, with stress < 0.2 regarded as a good fit. (B) The bar plot graphs show effect sizes (r coefficients) on the x axes for significant associations between microbial taxa and the indicated ileal gene expression, clinical disease activity (PCDAI), and clinical groups obtained from MaAsLin, while controlling for age, gender, BMI, antibiotic exposure, and NOD2, FUT2, and ATG16L1 risk allele carriage (P < 0.05, q < 0.25 were considered significant). Representative genes are from the APOA1 (APOA1 and CXCL9) and DUOX2 (DUOXA2, MUC4, and LCT) gene coexpression signatures. The comprehensive presentation of these data with P values, q values, and taxon abundance can be found in Supplemental Excel file 16.

Figure 5. The ileal microbial community in patients with IBD and Ctls.

(A) Univariate analysis (LDA Effect Size) was performed to test for association between ileal microbial composition and disease state in 240 CD (48 cCD, 178 iCD, and 14 CD with no indicated ileal involvement), 163 non-IBD Ctl, and 56 UC subjects from the RISK cohort that had not received antibiotics prior to endoscopy. The cladogram illustrates the output of this univariate analysis by demonstrating the relationship between the different bacterial taxa and disease state. Colored nodes from the center to the periphery represent marked phylum (p), class (c), order (o), family (f), genus (g), and species (s) differences detected between groups for Ctls (green), CD (red), and UC (blue). Only phylum, class, order, and family levels are indicated on the right side of the cladogram. (B) Fold change for each taxon was calculated by dividing the mean abundance in the cases (cCD [48 patients] or iCD [178 patients]) by that of the Ctls (154 patients) and is shown for microbiota with differential abundance in CD compared with Ctl by univariate analysis. The cCD group was further subdivided to cCD with abnormal histological features (25 patients) and cCD with normal histology (18 patients).

Figure 4. Gene expression signature and biologic pathways associated with mucosal ulceration.

(A) Image of ileal deep ulcers and Venn diagram showing 345 genes from the core iCD gene signature, which overlap with differentially expressed genes (Audic Claverie method with Benjamini-Hochberg FDR correction [0.05], fold change ≥ 1.5), between iCD-DU and iCD-noDU. A heat map of averaged gene expression and that after hierarchical clustering for the iCD-DU gene list for the indicated clinical subgroups. Specific upregulated (red) and downregulated (blue) genes are listed. (B) IPA functional annotation enrichment analyses to detect upstream regulators and biologic functions associated with the mucosal ulceration (iCD-DU) gene expression signature are shown. The activation z score for biologic function enrichment (P value range: 6E-8 to 3E-26) and upstream regulator enrichment (P value range: 3E-10 to 1E-56) depicts the degree of activation or suppression of a given regulator or function. (C) Immune cell–type enrichment of upregulated genes for iCD-DU (238 of 345 genes) and upregulated genes for iCD-noDU (524 of 936 genes) was determined using the Immunological Genome Project data series as a reference through ToppGene (20). Ileal enrichment for a given immune cell class (e.g., granulocytes [GN]) is illustrated by colored bars on the x axis, with the significance for each individual cell subtype within the class shown as the –log₁₀ (P value) on the y axis. pDC, plasmacytoid DC; MF, macrophages; MO, monocytes.

Figure 3. The ileal APOA1 gene coexpression signature is specific to all forms of CD.

(A) The Venn diagram shows the overlap of 179 genes differentially expressed 2-fold between patients with cCD and UC and genes within the DUOX2 and APOA1 coexpression modules. (B) A heat map of 93 of the above-mentioned 179 genes after averaging and hierarchical clustering of gene expression, with a fold change of 2.5 between the cCD and UC training cohorts for the indicated clinical subgroups. The training cohort included 17 cCD with abnormal histological features, 17 cCD with normal histology, and 45 UC. The independent validation cohort included 14 cCD with abnormal histological features, 11 cCD with normal histology, and 28 UC. The blue arrow indicates APOA1, and the red arrow indicates DUOX2. (C) The results of functional annotation enrichment analyses (IPA, Ingenuity Systems) using the DUOX2 and APOA1 gene coexpression signatures. Upregulated transcription factors and target genes are shown in orange and red, respectively, while downregulated transcription factors and target genes are shown in blue and green, respectively. Biologic functions associated with these groups of genes are also shown.

Figure 2. Functional annotation enrichment analysis of the APOA1 and DUOX2 gene coexpression signatures.

(A) Representative ileal DUOX2, APOA1, and 4HNE immunohistochemistry (original magnification, ×40) for Ctl (n = 3) and iCD (n = 7). Arrows indicate apical DUOX2 and intracellular APOA1 enterocyte staining and lamina propria 4HNE staining for lipid peroxidation. (B) Functional annotation enrichment analyses of the DUOX2 and APOA1 gene coexpression signatures using ToppCluster (21) and Cytoscape (27). Upregulated transcription factors and biologic functions are shown as red circles and orange boxes, respectively, while downregulated transcription factors and biologic functions are shown as blue circles and boxes, respectively.

Figure 1. APOA1 and DUOX2 gene coexpression signatures define pathogenic ileal modules.

(A) Study groups included Ctl, UC, cCD, and iCD. Ileal 16S DNA sequencing and mRNA sequencing were used to define microbial profiles and gene expression signatures, respectively. (B) Venn diagram shows 1,281 genes (core iCD signature) that were differentially expressed (fold change ≥ 1.5) between 2 independent iCD and Ctl groups. (C) Heat maps of average gene expression for the clinical subgroups for DUOX2 or APOA1 gene coexpression signatures (Pearson correlation 0.98<|r|<1). Venn diagram indicates DUOX2 and APOA1 gene coexpression signature overlap. Selected upregulated (red) and downregulated (blue) genes are listed, with the graphs showing the average DUOX2, LCT, APOA1, and IFNG gene expression across groups. Differences between patient subgroups were tested using Kruskal-Wallis with Dunn’s multiple comparison test. *P < 0.05, **P < 0.01, ***P < 0.001 vs. Ctl. cCD-mic, cCD with abnormal histological features; cCD N histo, cCD with normal histology.