Intestinal epithelial cells related lncRNA and mRNA expression profiles in dextran sulphate sodium‐induced colitis

Abstract Intestinal epithelial barrier damage caused by intestinal epithelial cells (IECs) dysfunction plays a crucial role in the pathogenesis and development of inflammatory bowel disease (IBD). Recently, some studies have suggested the emerging role of long non‐coding RNAs (lncRNAs) in IBD. The aim of this study was to reveal lncRNAs and mRNA expression profiles in IECs from a mouse model of colitis and to expand our understanding in the intestinal epithelial barrier regulation. IECs from the colons of wild‐type mice and dextran sulphate sodium (DSS)‐induced mice were isolated for high‐throughput RNA‐sequencing. A total of 254 up‐regulated and 1013 down‐regulated mRNAs and 542 up‐regulated and 766 down‐regulated lncRNAs were detected in the DSS group compared with the Control group. Four mRNAs and six lncRNAs were validated by real‐time quantitative PCR. Function analysis showed that dysregulated mRNAs participated in TLR7 signalling pathway, IL‐1 receptor activity, BMP receptor binding and IL‐17 signalling pathway. Furthermore, the possibility of indirect interactions between differentially expressed mRNAs and lncRNAs was illustrated by the competing endogenous RNA (ceRNA) network. LncRNA ENSMUST00000128026 was predicted to bind to mmu‐miR‐6899‐3p, regulating Dnmbp expression. LncRNA NONMMUT143162.1 was predicted to competitively bind to mmu‐miR‐6899‐3p, regulating Tnip3 expression. Finally, the protein‐protein interaction (PPI) network analysis was constructed with 311 nodes and 563 edges. And the highest connectivity degrees were Mmp9, Fpr2 and Ccl3. These results provide novel insights into the functions of lncRNAs and mRNAs involved in the regulation of the intestinal epithelial barrier.


| INTRODUC TI ON
Inflammatory bowel disease (IBD) is a group of chronic, non-specific inflammatory conditions of the gastrointestinal tract that mainly include Crohn's disease (CD) and ulcerative colitis (UC). 1,2 During the last decades, the incidence and prevalence of IBD are increasing in the world, due to many factors such as environmental exposures, better detection techniques (eg colonoscopy), advances in healthcare infrastructure and so on. [3][4][5] However, the precise pathogenesis of IBD is still unknown. Recent evidence implicates that the disruption of intestinal epithelial barrier contributes to many intestinal diseases, including IBD. [6][7][8] The intestinal epithelial barrier consists of different types of intestinal epithelial cells (IECs) with intact tight junction (TJ). In the intestinal lumen, this barrier acts as the first physical and immunological protective wall against the toxins and pathogenic organisms. 9 Excessive death of IEC, altered expression and distribution of TJ proteins can cause the intestinal epithelial damage, which is characteristic of IBD. [10][11][12][13] For example, increased apoptosis was reported in the intestinal epithelium of both patients with UC or CD. 10,[14][15][16][17][18] In addition, clinical studies showed that changes in the expression and distribution of TJ proteins led to an altered TJ structure and barrier dysfunction in the active Crohn's disease. 13 Long non-coding RNAs (lncRNAs) are RNA transcripts more than 200 bp in length, but lack of protein-coding capacity. 19 Accumulating evidence has shown that lncRNAs act as important regulators in a variety of physiological and pathological processes, such as chromatin modification, transcriptional regulation, post-transcriptional regulation and so on. 19,20 Abnormal expression of lncRNA is closely related to the development of various diseases.
Increasing evidence suggests that lncRNA plays an important role in the pathogenesis of IBD. 21 For example, lncRNA interferon-γ-antisense 1 (IFNG-AS1), which is up-regulated in the intestinal mucosa of patients with actively inflamed IBD, has been reported to be a mediator of an inflammatory signalling cascade in IBD pathophysiology. 22,23 Interestingly, several studies have shown that some lncRNAs are involved in the modulation of intestinal epithelial barrier function. For instance, the apoptosis of IECs was reported to be regulated by lncRNA BC012900. 24 Overexpression of lncRNA H19 decreased the expression of TJ protein, zonula occludin 1 (ZO-1). 25 However, studies investigating the role of lncRNAs in IBD are still limited.
To date, there have been several studies that used high-throughput RNA-seq analysis or microarray technology to determine the coding or non-coding gene differences both in mouse colitis models and in human IBD samples. [26][27][28][29][30][31] However, in these previous studies, the whole colon tissue was utilized for high-throughput RNA-seq or microarray analysis. Given the complexity of the cell types of the whole colon, such as IECs, intestinal intraepithelial lymphocytes (IELs) and lamina propria lymphocytes (LPLs), it is necessary to purify and collect different cell populations for the characterization of lncRNAs.
In this study, we focused on IECs and sought to determine the dysregulated lncRNAs and mRNAs that may be involved in the regulation of intestinal epithelial barrier in IBD. Therefore, we estab-

| Animal model
Mouse model of DSS-induced colitis was established according to the protocol in our previous study. [32][33][34] Male, 6-to 8-week-old C57BL/6 mice were obtained from the Experiment Animal Center (Army Medical University, Chongqing, China). A total of 10 mice per group were used to establish animal models. Control mice were provided with distilled water for 7 days. And mice in the DSS group were provided with drinking water containing 3% DSS for 7 days.
Then, mice were sacrificed for haematoxylin-eosin (H&E) staining and IEC isolation. All procedures were approved by the University Committee on the Use and Care of Animals of the Army Medical University.

| Isolation and purification of IECs
Intestinal epithelial cells isolation was performed as previously described. 33 Briefly, the colon was removed and placed in a tissue culture medium (RPMI 1640, with 10% foetal calf serum).
Then, the colon was cut into 5-mm pieces followed by extensively rinsed with ice-cold PBS containing 2% foetal calf serum.

| RNA isolation and high-throughput sequencing
Total RNA was isolated using RNeasy mini kit (Qiagen). Then, the RNA concentration and quality were determined by the Qubit ® 2.0

| Analysis of differentially expressed genes (DEGs)
The analysis for DE mRNA and DE lncRNA was performed using R package edgeR. Differentially expressed genes with |log2(FC)| value >1 and P < .05 were considered as significantly modulated and retained for further analysis. This choice is motivated by the decision to maximize the sensitivity of this analysis, in order to perform a massive screening and identify candidate genes to be validated with a wider sample population with real-time PCR analysis.

| GO and KEGG pathway analysis
The GO analysis (http://www.geneo ntolo gy.org) for biological processes, cellular components and molecular function and the KEGG pathway analysis (http://www.genome.ad.jp/kegg) were performed via enrich R package. P < .05 was considered to be a statistically significant enrichment.

| Quantitative real-time polymerase chain reaction (qRT-PCR) confirmation
Six lncRNAs and four mRNAs were selected for validation. The qRT-PCR was performed as previously described. 33 Total RNA was extracted from the IECs using RNAiso Plus (Takara), according to the manufacturer's instructions. The total RNA was reverse transcribed into complementary DNA (cDNA) using a SuperScript FirstStrand Synthesis System RT-PCR Kit (Invitrogen). The sequences of the primers used in the present study were listed in Table S1. The relative gene expression ratio was analysed using the 2 − ΔΔC t method. 35 Glyceraldehyde 3-phosphate dehydrogenase (GAPDH) was used as a reference for normalization. The qRT-PCR results were statistically tested with a Student's t test. A P < .05 was considered statistically significant.

| Profiles of DE mRNAs and lncRNAs
C57BL/6J wild-type (WT) mice were treated with 3% DSS for 7 days; changes in bodyweight, colon length, mortality and histology between the wild-type mice (Control group) and DSS-induced mice (DSS group) confirmed the successful induction of colitis. The DSS group exhibited significant weight loss after 3-4 days of DSS treatment ( Figure S1A). In addition, colon length was decreased in the DSS group ( Figure S1B). And there was a 23% rate of mortality in DSS group ( Figure S1C). As shown in Figure S1D, remarkable mucosal inflammatory infiltration also occurred in the colon mucosa from DSS group. Then, intestinal epithelial cells (IECs) were isolated from the colons of the two groups. The epithelial cell marker E-Cadherin was used to identify the purity of the IECs. As shown in Figure S1E, the purity of the IECs is about 96%. Then, a high-throughput sequencing F I G U R E 1 DE mRNAs in IECs from the colon of DSS-treated mice and Control mice. A, The number of overlapping mRNAs in Control group and DSS group was displayed by Venn diagram. B, Hierarchical clustering analysis of differentially expressed mRNAs between DSS group and Control group was showed by the heat map. Red colour represents up-regulated mRNAs, and the green colour represents downregulated mRNAs. Up-regulated and down-regulated mRNAs were showed by the scatter plot (C) and volcano plot (D). Red points indicate up-regulated mRNAs, and blue points indicate down-regulated mRNAs was used to detect the levels of mRNAs and lncRNAs in IECs from the Control group and DSS group.
A total of 28 673 mRNAs were detected in the Control group and DSS group. As shown in Figure 1A, 2129 were detected only in the Control group, 2464 were detected only in the DSS group, and 24 080 were detected in both groups. Then, using the following criteria: Pvalue <.05 and fold change >2, we identified 254 up-regulated mRNAs and 1013 down-regulated mRNAs in the DSS group compared with the Control group. A heat map was used to show the hierarchical clustering features of altered mRNAs ( Figure 1B). The variation of mRNA expression between the DSS group and the Control group was assessed by the scatter plot ( Figure 1C). The volcano plot was used to visualize the differential mRNA expression between two groups ( Figure 1D).  Figure 2B, the hierarchical clustering features of altered lncRNAs were presented in a heat map.
The variation of lncRNA expression between the DSS group and the Control group was assessed by the scatter plot ( Figure 2C). The volcano plot was used to visualize the differential lncRNA expression between two groups ( Figure 2D). All dysregulated lncRNAs between the DSS group and the Control group were displayed in Table S2.

| Features of DE lncRNAs
In order to detect the expression characteristics of DE lncR-NAs, the length distribution and chromosome distribution of the  Figure 3A). Based on the association with neighbouring protein-coding genes by genomic architecture, lncRNAs can be placed into one or more of five broad categories: sense, antisense, bidirectional, intronic or intergenic. 36,37 We then analysed the frequency distribution of DE lncRNAs in each category. As shown in Figure 3B, the fractions of exonic_sense, intronic_sense, exonic_antisense, in-tronic_antisense, bidirectional and intergenic lncRNAs were 29%, 12%, 12%, 1%, 6% and 42%, respectively. In addition, the chromosome distribution of DE lncRNAs was shown in Figure 3C

| GO annotation and KEGG pathway analysis of DE mRNAs
GO analysis and KEGG pathway analysis were performed on significantly dysregulated mRNAs in the DSS group and the Control group.
As shown in Figure

| Construction of the ceRNA network
With a threshold of P-value <.05 and |Log2 FC| > 8, 20 mRNAs that would be further studied were selected to construct the ceRNA network. The interactions of lncRNA-miRNA-mRNA were calculated by miRanda. Cytoscape software was used for network import and visualization. The selected mRNAs were listed in Table S3. For instance, Dnmbp was predicted to interact with 3 miR-NAs, including mmu-miR-7035-5p, mmu-miR-194-2-3p and mmu-miR-6899-3p. LncRNA ENSMUST00000128026 was predicted to act as a ceRNA and compete for binding to mmu-miR-6899-3p, thereby regulating the expression of Dnmbp. Tnip3 was predicted to interact with mmu-miR-351-5p and mmu-miR-6899-3p. LncRNA NONMMUT143162.1 was predicted to competitively bind to mmu-miR-6899-3p, regulating Tnip3 expression.

| The protein-protein interaction (PPI) network construction
Protein-protein interaction networks were constructed according to the information of the STRING database. The PPI networks analysis included 311 nodes and 563 edges. The highest connectivity degrees were Mmp9 (matrix metallopeptidase 9, degree = 27), Fpr2 (formyl peptide receptor 2, degree = 17) and Ccl3 (C-C motif chemokine ligand 3, degree = 17). The PPI networks of the top 20 nodes were showed in Figure S2. In addition, the top 20 core genes and their corresponding degree were shown in Table S4.
Moreover, MCODE was used to screen the modules of PPI networks ( Figure 7A-C). Then, GO analysis and KEGG pathway analysis were performed to analyse the functional annotation of the genes involved in the top three significant modules ( Figure 7D-G). As shown in Figure 7, the GO enrichment analysis and KEGG pathway analysis of the modules were involved in collagen fibril organization, extracellular structure organization, extracellular matrix organization, complex of collagen trimers, extracellular matrix component, fibrillar collagen trimer, extracellular matrix structural constituent, platelet-derived growth factor binding, growth factor binding, protein digestion and absorption, ECM-receptor interaction, proteoglycans in cancer and so on.

| D ISCUSS I ON
The intestinal epithelial barrier, which is primarily formed by a single-cell layer of IECs, is essential in the maintenance of intestinal homoeostasis. Integrated IECs and TJs between IECs are the main components of the barrier. 38 Both the increased level of IEC death and the disruption of TJ contributed to the pathogenesis and development of IBD. Recently, the importance of lncRNAs in IBD has become increasingly obvious. RNA-seq analysis or microarray technology was performed to identify the dysregulated mRNAs and lncRNAs using the whole colon RNA isolated from the mouse models of colitis or IBD patients. 28,30 However, due to the diversity of cell types in the whole colon, a gene may be significantly changed in one cell type, but in other cell types, the differential effects can be attenuated or eliminated completely. Thus, it is necessary to isolate and purify different cell types and conduct genetic analysis on them.
In the present study, focusing on the role of IECs in the pathogenesis and development of IBD, we investigated the mRNAs and lncRNAs expression profile for IECs which were isolated from the colons of wild-type mice and DSS-induced mice.
The RNA-seq analysis exhibited that 1267 mRNAs and 1308 ln-cRNAs were significantly differentially expressed in the DSS group compared with the Control group (Figures 1 and 2). We selected 6 lncRNAs and 4 mRNAs for validation purposes (Figure 4). The results of qRT-PCR analysis suggested that up-regulated and down-regulated lncRNAs and mRNAs were consistent with the data of RNA-seq analysis. A previous study, using whole colon sections for RNA-seq, reported 12 mRNAs (11 up-regulated and 1 down-regulated) that are differentially expressed in dextran sodium sulphate (DSS) and interleukin-10-deficient mice as well as in IBD patients. 30 Compared with these 12 mRNAs, 7 up-regulated mRNAs match with our data: Mmp3, Hcls1, Il1b, Lcn2, Plek, Steap4 and Ubd.
Subsequently, the GO analysis was performed to identify the biological function of DE mRNAs ( Figure 5A). In the present study, the enriched GO terms were significantly associated with TLR7 signalling pathway, IL-1 receptor activity, BMP receptor binding and so on. In addition to forming a physical barrier against foreign antigens, IECs also participate in the innate immune response since IECs express pattern recognition receptors, such as TLRs. TLR7, a member of TLRs family, is known as an important regulator of innate immunity.
A previous study has reported that TLR7 agonists Imiquimod could ameliorate DSS-induced acute colitis. 39 IL-1 cytokines are key mediators in immune regulation and inflammatory processes. For example, excess IL-1β expression was found in CD and other inflammatory conditions of the gut. 40,41 In addition, IL-1β has been shown to cause a functional opening of intestinal TJ barrier. 42 Al-Sadi et al 43 demonstrated that the IL-1β-induced increase in TJ permeability was mediated by an increase in MLCK expression and activity. It is known that not only the ratio between the proliferation and cell death, but also the balance between the proliferating progenitor cells and differentiating IECs, contributed to the maintenance of intestinal epithelial homoeostasis. 6 BMPs, which belong to the transforming growth factor-β superfamily, were known to regulate the epithelial proliferation and differentiation. 44 Conditional inactivation of the type I BMP receptor Bmpr1a in mice resulted in disturbed homoeostasis of intestinal epithelial regeneration. 45 The KEGG signalling pathway analysis showed that the most significant pathways were involved in IL-17 signalling pathway, inflammatory bowel disease (IBD), focal adhesion, cytokine-cytokine receptor interaction and so on ( Figure 5B). IL -17 cytokines were thought to induce the mucosal inflammation but also be involved in the restitution and repair of the intestinal mucosa after resolution of inflammation. 46 Increased expression of IL-17 family members was found in both human IBD and animal models of colitis. 47  Recent research has demonstrated that lncRNAs can act as miRNA sponges to protect mRNAs from miRNA inhibition. 54,55 Then, we constructed a ceRNA network to better understand the function of lncRNAs ( Figure 6). For example, Dnmbp, also known as Tuba, is a Cdc42-specific guanine-nucleotide-exchange factor. In IECs, Dnmbp can control the shaping of cell junctions through the local activation of Cdc42 and its effectors, suggesting that Dnmbp may contribute to the function of intestinal epithelial barrier. 56,57 Interestingly, our qRT-PCR analysis results showed the decrease of Dnmbp in the DSS group. According to the ceRNAs network, lncRNA ENSMUST00000128026 was predicted to bind to mmu-miR-6899-3p and further regulate Dnmbp expression. Therefore, lncRNA ENSMUST00000128026 may regulate the intestinal epithelial barrier function through affecting the Dnmbp expression.

F I G U R E 6
The ceRNA network construction. In the network, red represents up-regulation and green represents down-regulation. LncRNAs, miRNAs and mRNAs are presented as quadrilateral, arrowheads and hexagons, respectively. And grey edges indicate lncRNA-miRNA-mRNA interactions However, further research is needed to elucidate the function of ENSMUST00000128026 in IBD.
Tnip3 (TNFAIP3 interacting protein 3), also known as ABIN-3, can negatively regulate the nuclear factor-κB (NF-κB) activation in response to tumour necrosis factor (TNF) and lipopolysaccharide (LPS). 58 In addition to the pro-inflammatory and anti-inflammatory effects, NF-kB can also mediate wound healing in IECs during the inflammatory processes. As a negative regulator of NF-κB activation, Tnip3 may contribute to the control of wound healing in IECs. It has been reported that Tnip3 expression is significantly up-regulated in F I G U R E 7 Functional annotation of top three modules in the PPI network. A, Module 1, B, Module 2, C, Module 3. GO enrichment analysis (D-F) and KEGG pathway analysis (G) were performed for DE mRNAs in the modules human masticatory mucosa during wound healing. 59 To date, the role of Tnip3 in IBD is still unknown. Recently, a study showed that Tnip3 was one of the top 10 down-regulated genes when comparing between the long-duration and short-duration UC patients. 60 Interestingly, in the present study, we found that Tnip3 was predicted to be up-regulated when comparing between the DSS and control mice. These findings suggest that Tnip3 may play an important role in IBD. However, the association between Tnip3 and IBD remains to be elucidated. Furthermore, the ceRNAs network showed that lncRNA NONMMUT143162.1 was predicted to regulate Tnip3 expression by competitively binding to mmu-miR-6899-3p. Therefore, we will further illuminate the concrete mechanisms in future studies to understand the relationship between this novel RNA crosstalk and IBD.

CO N FLI C T O F I NTE R E S T
The authors declare that they have no conflict of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings in the current study are available from the corresponding author upon reasonable request.