Transcriptome and methylome profiling in a rat model of irritable bowel syndrome induced by stress.

Irritable bowel syndrome (IBS) is a common gastrointestinal disorder that is associated with psychological stress. However, the full landscape of IBS‑related epigenetic factors remains unveiled and needs to be elucidated. The water‑avoidance stress (WAS) method was used to induce a rat IBS model. Each rat was monitored, and its defecation and behavior were recorded. Total colon RNA was isolated and subjected to Affymetrix GeneChip analysis. Reduced Representation Bisulfate Sequencing (RRBS) was applied to determine the genome‑wide methylation pattern in both IBS and control rats. Rats with IBS egested a significantly increased amount of dry and loose stools compared with the control animals, without significant changes in body weight. Compared with the control group, 309 genes were upregulated and 224 genes were downregulated in the colon of the IBS rats. Notch signaling and focal adhesion were increased in the differentially expressed genes (DEGs). A total of 541 genes had significant lower methylation level and 626 genes had significantly higher methylation level in their promoter regions. Adherens junction and leukocyte transendothelial migration were enriched in the differentially methylated genes (DMGs). Few genes were identified in common in both DEGs and DMGs, suggesting that gene expression was not altered by promoter methylation. Reverse transcription‑quantitative polymerase chain reaction validation revealed that the mRNA levels of SSX2IP, PARD3 and VCL were significantly downregulated in the IBS group, in accordance with hypermethylation of their promoters. In summary, the present study used a WAS‑induced IBS rat model to provide transcriptome and methylome profiling. Most DEGs were associated with Notch signaling and focal adhesion, and only a few were altered by promoter methylation. The present results demonstrated that psychological stress could influence the integrity of the intestinal mucosa barrier and regulate inflammatory response.

Introduction

Irritable bowel syndrome (IBS) is one of the most commonly diagnosed gastrointestinal (GI) disorders, with a 10-15% global incidence (1,2). IBS patients usually have chronic abdominal pain/discomfort and bowel habit change, such as diarrhea and/or constipation (3,4). According to different symptoms, IBS is divided into the following three subtypes: IBS-D (frequent occurrence of diarrhea), IBS-C (frequent occurrence of constipation), and IBS-M (switch between diarrhea and constipation) (3). IBS is frequently associated with psychiatric disorders, which can significantly affect patients' quality of life (2,5,6). It is estimated that the direct cost of IBS management is >$1,000 per patient per year in the USA, without taking into account the substantial economic burden and the loss of productivity (7).

The etiology and pathophysiological mechanisms of IBS are largely unknown (8-10). Increasing evidence indicates that inflammation and immune function are dysregulated in IBS (4,11-13). Many psychiatric disorders, such as depression and dysphoria, are closely associated with IBS, and stress is also noticeably associated with IBS (14-16). Emotional and physical stress has been linked to the first onset, the exacerbation, and the treatment outcome of IBS patients (17,18). It has been reported that stressful events at young ages (including childhood trauma, maternal neglect, sexual or emotional abuse) can be significant risk factors that contribute to IBS in adulthood (19).

Acute stress events in adulthood can induce an immediate IBS-like phenotype, and this has been widely used to generate animal models of IBS (20). One of the most widely used paradigms in IBS research is the water-avoidance stress (WAS) method (21). It is commonly accepted that the WAS method can induce immediate gut motility change and exaggerate visceral pain (22-25). A WAS stimulation of 1 h is normally sufficient to reduce colorectal distension (CRD) in a rat model by inducing delayed visceral hypersensitivity (26).

Even though a number of IBS-related factors have been identified in animal and patients studies (1,12,13,24,25), a full landscape of IBS-associated genes at the molecular level is needed to get a comprehensive understanding of IBS pathogenesis. A previous small bowel mucosal gene expression study has reported that both transcriptome expression and immune-related mechanisms are altered in IBS-D patients (27). Another DNA microarray analysis has demonstrated that IBS is associated with abnormal serotonin metabolism and vitamin D deficiency (28). Currently, most omics studies only report the observational results of IBS-related differentially expressed genes (DEGs). The investigation of epigenetic factors, such as DNA methylation, remains absent.

In order to better understand the pathogenesis of IBS and the symptom etiology, a stress-induced IBS animal model was used to study the transcriptome and methylome profiles, of IBS.

Materials and methods

Ethical approval

The protocols complied with regulations regarding the care and use of laboratory animals published by National Institutes of Health. The present study was approved by the Ethics Committee of Beijing Friendship Hospital (Beijing, China). A total of 20 specific-pathogen-free (SPF) male Wistar-Kyoto rats weight 180-200 g and aged 8 weeks were purchased from Vital River Laboratory Animal Technology Co., Ltd. (Beijing, China). All the animals were maintained at 22±2°C with a 12-h light/dark cycle and free access to standard laboratory food and water, while environmental noise was kept to a minimum range.

Animal experiments

A total of 10 rats were randomized into the control group and the IBS group (n=5 per group). WAS paradigm was adopted to induce the IBS model (29). Every morning from 8:00-9:00 am, each rat in the IBS group was subjected to 1 h stress test as following: each animal was placed on a small platform (10×10 cm) that was raised 5 mm above the water level, and this test was repeated once a day for 10 consecutive days. The control animals were treated the same way except that each control rat was kept in a single cage alone without stress. All behavior experiments were performed in triplicate.

Gastrointestinal syndrome and behavior test

Each rat was weighed daily at 8:00 am. From 8:00-9:00 am, the defecation number of each rat in both the IBS and control groups was recorded from day 1 to day 10. Feces with regular ellipsoidal shape were defined as dry stools, whereas the wet amorphous feces without solid pieces were defined as loose stools. All rats were sacrificed on day 10, and the colon tissues were collected and divided into two parts, one for mRNA microarray analysis, and the other for Reduced Representation Bisulfate Sequencing (RRBS) methylation analysis.

RNA extraction and reverse transcription-quantitative polymerase chain reaction (RT-qPCR)

Total RNA extraction

Colon tissues were thoroughly ground with liquid nitrogen, 100 mg tissue sample added to 1 ml fresh TRIzol (Thermo Fisher Scientific, Inc., Waltham, MA, USA). A total of 200 µl chloroform was added and agitated for 15 sec. This was left at room temperature for 3 min to allow for phase separation. The sample was centrifuged at 12,800 × g in a microcentrifuge at 4°C for 15 min. The upper 3/4 of the aqueous phase was collected. The aqueous phase was combined with an equal volume of isopropanol. The sample was left for 10-20 min at -20°C for precipitation. Then centrifuged at 12,800 x g in a microcentrifuge at 4°C for 15 min. The pellet was washed with 1 ml 75% pre-cold ethanol. The sample was then centrifuged at 12,800 × g in a microcentrifuge at 4°C for 10 min. The pellet was resuspended with 20-50 µl DEPC-H2O (65°C for 15 min). The RNA concentration and A260/A280, A260/A230 value was measured in a NanoDrop 1000.

RT-PCR

A total of 2 µg RNA sample was added to 10 µl DEPC water for reverse transcription. The RT reaction mixture was prepared using the High Capacity cDNA Reverse Transcription kits (Thermo Fisher Scientific, Inc.; cat. no. 4368813) as follows (total per sample 10 µl): 10X RT buffer 2.0 µl, 25X dNTP mix (100 mM) 0.8 µl, 10X RT Random primers 2.0 µl, Multiscribe reverse transcriptase 1.0 µl, RNase inhibitor 1.0 µl and Nuclease-free H2O 3.2 µl. A total of 10 µl master mix into the tube containing 10 µl RNA sample. The sample was briefly centrifuged at 300 × g, at 4°C for 2 min and then placed into a BioRad C1000 Thermal Cycler and run program to synthesize cDNA as follows: 25°C for 10 min, 37°C for 120 min, 85°C for 5 min, 4°C hold. Total cDNA was stored at -20°C until ready for qPCR analysis.

qPCR

Primer sequences used for amplification of the SSX2IP promoter were: 5′-ATGGGAGATTGGATGACTGT-3' (forward), and 5′-CCCGACTGGCTATGATATTTTG-3′ (reverse); for the PARD3 promoter were: 5′-GGAGATGGCCGCATGAAAGTT-3′ (forward), and 5′-CTCCAAGCGATGTACCTGTAT-3′ (reverse); and for the VCL promoter were: 5′-CGCTGGCCTCCATAGACTC-3′ (forward), and 5′-CCTAGCATCTTGCAGGTTCCTA-3′ (reverse). The reaction mix was prepared as follows: Prepare sense/anti-sense primer solution in final concentration 1 µM. Mix qPCR reaction system (total per sample 12 µl): Fast SYBR™-Green Master Mix (Thermo Fisher Scientific, Inc.) 6 µl, Primer mix (F+R, 4 µM) 2.5 µl x2 and cDNA 1 µl. The sample was centrifuged at 300 x g, 4°C, 2 min. qPCR was performed using the following program: Stage 1 pre-denaturation 1 cycle at 95°C for 2 min; stage 2 thermal cycle 40 cycles at 94°C for 15 sec, 56°C 20 sec, 72°C 30 sec; Stage 3 melting stage (Dissociation stage). The RNA copy numbers were analysed by the 2-ΔΔCq method (30) via QuantStudio 7 Flex (Thermo Fisher Scientific, Inc.) and compare it with housekeeping gene.

mRNA microarray analysis

Total colon tissue RNA was isolated with the TRIzol reagent (Thermo Fisher Scientific, Inc.), according to the manufacturer's protocol. Affymetrix GeneChip Rat Genome 230 Array 2.0 (Affymetrix; Thermo Fisher Scientific, Inc.) was applied in the present study and all associated procedures were conducted by GeneChem Co., Ltd. (Shanghai, China), according to the standard operating procedure. All CEL data were corrected by RMA method. Then, log2 transformation along with quantile normalization were applied subsequently. Expressional data were analyzed in R software using the Limma package to identify differentially expressed genes.

RRBS methylation analysis

Colon tissue DNA was extracted with the TIANamp Genomic DNA kit (Tiangen, Beijing, China) according to the manufacturer's protocol. A total of 500 µg DNA from each sample was used for library preparation. Total DNA was digested with MspI, end-repaired and A-tail was added with Klenow fragment and ligated with 5mC-methylated paired end sequencing adapters (Illumina, Inc., San Diego, CA, USA). All enzymes were purchased from New England Biolabs, Inc. (Ipswich, MA, USA), unless otherwise mentioned. Fragments were purified by electrophoresis, and converted with the EZ DNA Methylation kit (Zymo Research, Irvine, CA, USA). On an Illumina HiSeq2000 platform, each library was subjected to sequencing by Annoroad Co. (Shanghai, China), using standard operating procedures. Rat genomic sequences (Rnor6.0) were downloaded from the Ensembl database for RRBS reads mapping. Adapter sequences were mapped to Rnor6.0 using Bismark, and the remained sequences were removed. Using the methyl kit package, methylation data were analyzed in R software to identify the differentially methylated region.

Statistical analysis

Statistical tests and data visualization were performed using R 3.3.1 (www.r-project.org). All statistical tests were two-tailed and False Discovery Rate (FDR) was adjusted. P<0.05 was considered to indicate a significant difference. Packages of dplyr, plyr, and reshape2 were used for data sorting and restructuring. GEOquery, annotate, DOSE and ReactomePA were used for annotating and enriching genes. VennDiagram, pheatmap, and ggplot2 were used for visualization of results. All R packages were obtained from Bioconductor (https://www.bioconductor.org) and CRAN webset (https://cran.r-project.org/web/packages).

Results

WAS paradigm results in IBS-like phenotype in rats

The present results demonstrated that the average number of bowel movement in the IBS group was significantly increased compared with the control group from day 1 to day 10 (Fig. 1A). On day 10, compared with the control group, the number of defecation in the IBS group was increased by nearly 10 times (10.67±2.12 vs. 1.33±0.21; P=0.001; Fig. 1B). When stratified by the status of stools, both numbers of dry and loose stools in the IBS group were larger compared with the control group (Fig. 1C). However, the total body weight in both groups was not statistically different from day 1 to day 10 (Fig. 1D and E).

Figure 1

Effect of WAS stimulation on animal bowel movements. (A) Stool counts from day 1 to day 10 in the animals of the control and IBS groups. (B) Defecation counts of IBS and control animals on day 10. (C) Total counts of dry and loose stools from day 1 to day 10. (D) Animal total body weight measurements from day 1 to day 10. (E) Animal total body weights on day 10. Data were presented as mean ± standard deviation (n=5 per group). All behavior experiments were performed in triplicate. **P<0.01, with comparisons indicated by lines. WAS, water-avoidance stress; IBS, irritable bowel syndrome; NS, not significant.

Gene expression changes in IBS group rats

Results from the colon tissue mRNA microarray demonstrated that 309 genes were upregulated and 224 genes were downregulated in the IBS group compared with the control group. The mRNA levels in the IBS group were changed by 30% compared with the control group (Fig. 2A). Gene Ontology (GO) enrichment analysis revealed that 30 GO biological process terms were enriched, and some of these enriched terms had important roles in stimulus, cell communication, phosphorus metabolic process (Fig. 2B). Kyoto Encyclopedia of Genes and Genomes (KEGG) and BIOCARTA enrichment analysis revealed that 22 terms were enriched, including notch signaling, neuroactive ligand-receptor interaction and focal adhesion (Fig. 2C).

Figure 2

Gene expression analysis of IBS and control rats. (A) Heatmap of representative differentially expressed genes between IBS and control groups (>30% fold-change; FDR adjusted P<0.05). (B) Bubble plot of representative GO enrichment results of differentially expressed genes. (C) Bubble plot of representative KEGG and BIOCARTA enrichment results of differentially expressed genes. The bubble size reflects the number of genes enriched in a certain term, and the bubble color reflects the P-value. IBS, irritable bowel syndrome; FDR, false discovery rate; GO, gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Gene methylation changes in IBS group rats

RRBS was applied to determine the genome-wide methylation patterns in both IBS and control rats. The distribution of CpGs in different functional components of the genome is illustrated in Fig. 3A; nearly 97% of GpGs were located in non-coding areas, while 9% of GpGs were in promoter regions. Compared with the control group, the IBS rats had significantly lower CpG methylation levels, especially in the 5'UTR and promoter regions (Fig. 3B). CpG was detected in the promoter regions of all genes. A total of 541 genes had significantly lower methylation levels in their promoter regions in the IBS group compared with the control group (fold-change >2; P<0.05). A total of 626 genes had significantly higher methylation in their promoter region in the IBS group compared with the control group (fold-change >2; P<0.05). In addition, 17 terms had KEGG enrichment, and some of these enriched terms were involved in bacterial invasion of epithelial cells, adherens junction and leukocyte transendothelial migration (Fig. 3C). Schematics showing the changed promoter methylation levels in genes involved in adherens junctions (hsa04530) and leukocyte transendothelial migration (hsa04670) are presented in Fig. 3D and E.

Figure 3

Gene methylation analysis of IBS and control rats. (A) Distribution of CpGs in different functional components of the genome. (B) CpG change (%) in the IBS group compared with the control group of different genome functional components. (C) Bubble plot of representative KEGG enrichment results of genes with differentially methylated promoters. (D) Promoter methylation levels in the adherens junction pathway. (E) Promoter methylation levels in the leukocyte transendothelial migration pathway. Up-methylated genes were marked as red; down-methylated genes were marked as in blue; other genes were marked as green. IBS, irritable bowel syndrome; KEGG, Kyoto Encyclopedia of Genes and Genomes; UTR, untranslated region; CGI, CpG islands.

Gene expression and methylation co-analysis in the IBS rat model

The expression levels and methylation levels of all genes were aligned according to their locations on the chromosomes (Fig. 4A). The density plots of gene expression and methylation are illustrated in Fig. 4B. The peak locations of the curves indicated that the overall gene expression was promoted and the overall gene methylation was inhibited in the IBS group compared with the control group. Gene set enrichment analysis (GSEA) revealed that genes that located in the Chr3P26 region were significantly hypermethylated (Fig. 4C). In addition, the gene expression levels were negatively correlated with their promoter methylation levels (R=-0.250; P=0.008; Fig. 4D), which partially validated the biological effect of the hypermethylation in this region.

Figure 4

Gene expression and methylation landscapes in the IBS rat model. (A) Ratios of mRNA changes (upper panel), and log2 transformed ratios of promoter methylation (lower panel) between IBS and control groups. Each gene was aligned with respect to their own chromosome position. Each dot represents abundance changes for one gene, the corresponding mRNA and the corresponding promoter methylation level, respectively. (B) Upper panel: mRNA expression fold-change density histogram. The vertical dashed line represents equal level of expression between two groups (fold-change=1). Lower panel: Promoter methylation fold-change (log2 transformed) density histogram. The vertical dashed line represents equal level of methylation between two groups (log2 fold-change=0). (C) GSEA analysis indicated that genes located in the Chr3P26 region were largely differentially expressed between IBS and control groups (NES=1.363; P=0.009). (D) Gene expression levels were negatively correlated with their promoter methylation levels in the Chr3P26 region (R=-0.250; P=0.008). IBS, irritable bowel syndrome; GSEA, gene set enrichment analysis.

A Venn diagram was constructed to illustrate the overlap of differentially expressed or methylated genes (Fig. 5A). Only 3 genes were identified as hypomethylated and upregulated, and 5 genes were hypermethylated and downregulated in the IBS group, suggesting that the majority of differentially expressed genes were not changed by promoter methylation.

Figure 5

Gene expression and methylation co-analysis and RT-qPCR validation. (A) Venn diagram illustrating the overlap of differentially expressed genes and differentially methylated genes. (B) Heatmap of mRNA expression levels of representative differentially methylated genes (from the KEGG hsa04530 pathway) between IBS and control groups. (C) Promoter methylation change and expression fold-change of representative genes (from the KEGG hsa04530 pathway). (D) RT-qPCR validation of gene expression in KEGG hsa04530. Data were presented as mean ± standard deviation (n=5 per group). *P<0.05. RT-qPCR, reverse transcription-quantitative polymerase chain reaction; KEGG, Kyoto Encyclopedia of Genes and Genomes; IBS, irritable bowel syndrome; CREBBP, CREB binding protein; SSX2IP, SSX family member 2 interacting protein; PARD3, par-3 family cell polarity regulator; VCL, vinculin; TCF7, transcription factor 7; WASF1, WAS protein family member 1; NS, not significant.

In order to evaluate the biological effect of methylation change, the differentially methylated genes [CREB binding protein (CREBBP), SSX family member 2 interacting protein (SSX2IP), par-3 family cell polarity regulator (PARD3), vinculin (VCL), transcription factor 7 (TCF7), and WAS protein family member 1 (WASF1)] in the KEGG pathway 'adherens junction' (hsa04530) were extracted. The expression levels of those genes in each sample are presented in Fig. 5B. The relationship between the promoter methylation and the gene expression levels is illustrated in Fig. 5C. SSX2IP, PARD3 and VCL were all hypermethylated and downregulated, even though the expression levels did not reach a significant level (FDR adjusted P>0.05). However, in the RT-qPCR analysis, both SSX2IP and PARD3 were significantly downregulated in the IBS group compared with the control group (Fig. 5D). These data indicated that DNA hypermethylation only mildly affected gene expression in this IBS model, and this mild effect could be easily overlooked by the strict multiple testing correction procedure.

Discussion

The etiology and pathophysiological mechanisms of IBS remain unclear. The brain-gut axis, immune regulation defection and risk gene mutations are all suggested to be associated with IBS. Further studies of IBS pathogenesis could help us to identify high risk populations and develop new clinical therapies to treat IBS patients. A full landscape of IBS-related genes at the molecular level is needed to get a comprehensive understanding of IBS pathogenesis. Stress-altered signaling pathways are involved in the pathogenesis of IBS (31); therefore, stress-induced IBS animal models are well accepted and widely used (1,10,18,24). The present study established a classical WAS-induced IBS rat model, which exhibited diarrhea-like bowel habits (21,24), symptoms that are often seen in IBS-D patients. Then, expression and methylome profiling was performed in the WAS-induced IBS rat model. As a fundamental mechanism research, the contribution of the present study may have limited clinical application; however, it proved that psychological stress resulted in the onset of IBS and suggested that antidepressants could be used in IBS therapy.

Transcriptome profiling revealed that 309 genes were upregulated and 224 genes were downregulated in the colon of IBS rats. Compared with the control group, IBS could induce significant changes in multiple genes, and most of the changed genes were involved in several pathophysiological functions, including cell communication, phosphorus metabolic process, Notch signaling, neuroactive ligand-receptor interaction and focal adhesion. Previous studies have also demonstrated that genes involved in these pathways were altered in the colon tissue of IBS rats compared with control rats (4,13,24,26,29). Attention should be paid to focal adhesion, which was strongly associated with gut mucosal barrier function and could be a possible inducer of IBS (32,33).

Few previous studies focused on the role of DNA methylation in IBS. Mahurkar et al (34) identified several DNA methylation markers from the peripheral blood mononuclear cells (PBMCs) in IBS patients; however, DNA methylation in IBS colon tissue has not been reported. The present methylome profiling revealed that in the colon mucosa of IBS rats, 626 genes were up-methylated and 541 genes were down-methylated in their promoter regions. Most of these IBS-related genes were associated with bacterial invasion of epithelial cells, adherens junction and leukocyte transendothelial migration. The genes involved in the epithelial cell junctions and gut mucosal barrier were significantly affected (32,33). However, the DMGs identified in the present study did not have a significant overlap with those identified previously (34). It is possible that WAS stimulation had different effects on PBMCs and intestinal epithelial cells.

The present transcriptome-methylome co-analysis revealed that the gene methylation levels were downregulated and the overall gene expression levels were upregulated in IBS. In addition, GSEA analysis demonstrated that genes located in the Chr3P26 region were significantly hypermethylated, and the expression levels of these genes were negatively correlated with their promoter methylation levels, which partially validated the biological effect of the hypermethylation in this region. By contrast, only 8 genes had synchronized significant changes in methylation and expression (hypomethylated and overexpressed, or hypermethylated and underexpressed), suggesting that most of the IBS-related genes were not changed by promoter methylation.

The biological effect of methylation in IBS was also investigated by RT-qPCR in the present study. Three genes that were enriched in the KEGG pathway 'adherens junction' were selected for validation, because they were synchronized but had no significant change in expression with their methylation level. Both SSX2IP and PARD3 were significantly down-regulated in the IBS group compared with the control group. VCL was downregulated but did not reach a significant level. It appeared that in this IBS rat model, the promoter methylation had a complicated role in gene expression. Other factors, such as transcription factors and histone modifications, could also have affected gene expression (35,36) by interfering with methylation-related changes.

In conclusion, the transcriptome and methylome profiling in the WAS-induced IBS model employed in the present study provided useful information on IBS-related genes at both the expression level and methylation level. The results revealed that psychological stress influenced the integrity of intestinal mucosa barrier and regulated the inflammatory response at the transcriptome and methylome level. The present molecular profiling of an IBS rat model could be beneficial to better understand the pathogenesis of IBS, and to potentially develop novel therapeutic targets to treat IBS.