Construction of competitive endogenous RNA network reveals regulatory role of long non‐coding RNAs in type 2 diabetes mellitus

Abstract Increasing epidemic of type 2 diabetes mellitus (T2DM) and its comorbidities makes it urgent to understand the pathogenesis and regulatory mechanism. However, little is known about the regulatory role of lncRNAs in diabetes. Here, we constructed a T2DM‐related competitive endogenous RNA (ceRNA) network (DMCN) to explore biological function of lncRNAs during the development of diabetes mellitus. This network contained 351 nodes including 98 mRNAs, 86 microRNAs and 167 lncRNAs. Functional analysis showed that the mRNAs in DMCN were annotated into some diabetes‐related pathways. Furthermore, mTOR‐centred subnetwork was extracted and ncRNA‐involved mTOR pathway was established. Finally, we validated that NEAT1 was potentially communicated with mTOR signalling target protein mLST8 via the association with miR‐181b. These findings provide significant insight into lncRNA regulatory network in T2DM.


Introduction
The increasing epidemic of T2DM and its comorbidities urge our understanding of the pathogenesis and regulatory mechanisms of T2DM [1]. The recent discovery of functional non-coding RNAs (ncRNAs) as specific gene expression regulators suggests that manipulating ncRNAs could be a novel therapeutic approach for combating metabolic disorders such as diabetes mellitus [2]. Numerous studies demonstrated that microRNAs have been implicated in diabetes mellitus [3]. However, less is known about the regulatory role of lncRNAs in diabetes. Recent studies provide evidence that endogenous RNAs influence each other's levels by competing for a limited pool of microRNAs [4]. ceRNA theory proposed that all types of RNA transcripts communicate through a new 'language' mediated by microRNA binding sites or microRNA response elements (MREs) [4].
Recently, various ceRNA networks have been constructed to gain a global view of regulatory interaction networks of ceRNAs in certain diseases. Yang et al. developed Starbase v 2.0 to systemically identify the RNA-RNA and protein-RNA interaction networks from 108 CLIP-Seq data sets generated by tumour samples [5]. The ceRNA theory has been applied to elucidate the pathogenesis of numerous diseases such as cancer, muscular dystrophy and neurodegenerative disease [6][7][8]. However, little is known about potential role of ceRNAs in metabolic disease, especially diabetes mellitus.
In this study, we constructed a T2DM-related ceRNA network (DMCN) to explore the biological functions of lncRNAs in the development of diabetes mellitus. After functional analysis, we found that the mRNAs in DMCN were annotated into some diabetes-related pathways. Then, mTOR-centred subnetwork was extracted and ncRNAinvolved mTOR pathway was established. Finally, we validated that NEAT1 was potentially communicated with mTOR signalling target protein mLST8 via the association with miR-181b.

Materials and methods
Collection of human T2DM-related genes (T-HOD) [9] and Online Mendelian Inheritance in Man (OMIM)'s Morbid Map [10]. A total of 632 human T2DM-related genes were obtained.

Identification of microRNA and T2DM-related gene pairs
The DIANA-TarBase (V6.0) [11] and miRTarBase (V4.5) [12] databases collected the microRNA-mRNA interaction pairs, which were manually curated and experimentally verified. From these two databases, we downloaded human microRNA-mRNA pairs that were verified by qPCR and Western blot analysis. After integrating process, we obtained 43,497 non-redundant microRNA-mRNA pairs. By integrating the collected human T2DM-related genes, we then screened T2DM-related microRNA-mRNA pairs from these non-redundant microRNA-mRNA pairs.

Identification of microRNA and T2DM-related lncRNA pairs
The human microRNA-lncRNA interaction pairs were identified from starBase (V2.0) [5], which were supported by ClipSeq experiments. A total of 10,213 experimentally verified microRNA-lncRNA pairs were identified. By combining microRNAs with T2DM-related target genes, we identified microRNA-lncRNA interaction pairs associated with T2DM.

Construction of human T2DM-related ceRNA network (DMCN)
The microRNA-lncRNA-gene interactions were integrated if microRNA-mRNA and microRNA-lncRNA pairs shared common microRNAs. By combining the collected T2DM-related microRNA-mRNA pairs with T2DM-related microRNA-lncRNA pairs, we identified 543 pairs of T2DM-related microRNA-lncRNA-gene interactions and then constructed the human T2DM-related ceRNA network.

Topological features of T2DM-related ceRNA network
The degree is the number of edges that connect to a certain node. K is the number of edges that connect to node i. It can be defined as follows: Deg (i) = K For a graph G: = (V,E) with n nodes, the betweenness B i stands for a node i as where S st is the number of shortest paths from s to t, and S st (i) represents the number of shortest paths from s to t that pass through a node i. The measure is normalized by the number of pairs of nodes except i, which is (n-1)(n-2).

Function annotation of ncRNAs to mTOR signalling pathway
To investigate biological functions of T2DM-related genes in the network, we used DAVID 6.7 for gene ontology (GO) functional annotation of T2DM genes and KEGG pathway for pathway enrichment [13]. After GO enrichment by DAVID tool, the GO BP terms with FDR value <0.1 were chosen as statistically significant terms. KEGG (Kyoto Encyclopedia of Genes and Genomes, http://www.genome.jp/kegg/) is a collection of online databases. T2DM genes were subjected to KEGG pathway analysis, and FDR values <0.1 were considered as statistically significant pathways.
We manually draw figure of mTOR signalling according to KEGG pathway database and annotated mirRNAs and lncRNAs to this pathway according to their interactions with genes involved in mTOR signalling pathway. The mirRNAs and their target genes were downloaded from miRTarBase (http://mirtarbase.mbc.nctu.edu.tw), which is an experimentally validated microRNA-target interactions database [14].

Cell culture
Human liver cell line HL7702 was purchased from the Type Culture Collection of the Chinese Academy of Sciences (Shanghai, China). Cells were grown in DMEM or RPMI-1640 supplemented with 10% foetal bovine serum, 100 U/ml penicillin and 50 mg/ml streptomycin at 37°C with 5% CO 2 .

Quantitative RT-PCR
Total RNA was isolated, and complementary DNA (cDNA) was generated with the use of a High Capacity cDNA Reverse Transcription Kit (ABI, Applied Biosystems, Carlsbad, CA, USA). cDNA from miRNAs was then generated with the use of TaqMan miRNA Reverse Transcription Kit (GenePharma, Shanghai, China). Quantitative reverse transcriptase PCR (RT-PCR) was conducted using TaqMan gene expression assays for NEAT1 MALAT1 and GAPDH or TaqMan Universal Master Mix II with TaqMan microRNA assay for miR-181b, miR-144-3p and U6.

Dual-luciferase reporter assay
Two Luc-NEAT1 constructs were generated including potential binding site (NEAT1-WT) or mutant binding site (NEAT1-Mut). HEK293 cells were seeded in 12-well plates, cultured for 24 hrs and then cotransfected with NEAT1-WT or NEAT1-Mut plasmid and miR-181b mimics or negative control sequences using Lipofectamine 2000 (Invitrogen, Carlsbad, CA, USA). The luciferase activity was measured 48 hrs after transfection using Dual-Luciferase reporter assay system (Promega, Madison, WI, USA).

3205
(Abcam, Cambridge, UK) overnight at 4°C. The membranes were then incubated with infrared (IR) fluorescent dye-labelled secondary antibody (Alexa Fluor; Molecular Probes, Eugene, OR, USA) for 1 hr. The bands were analysed using an IR Imaging System (LI-COR Biosciences, Lincoln, NE, USA), and band density was quantified using Odyssey 3.0 software (LI-COR Biotechnology, Lincoln, NE, USA) and normalized to b-actin.

RNA-binding protein immunoprecipitation assay
RNA-binding protein immunoprecipitation (RIP) assay was performed using EZ-Magna RIP Kit (Millipore, Billerica, MA, USA) according to the manufacturer's instruction. The enrichment of miRNAs binding AGO2, a key component of the microRNA-containing RISC complex, was used as the positive control, and U6 as the negative control. Briefly, cells were lysed with RIP lysis buffer, followed by incubated with RIP buffer containing magnetic beads conjugated with human anti-Ago2 antibody (Millipore). Proteinase K was used to digest the protein, and then immunoprecipitated RNA was isolated. Purified RNA was subjected to qRT-PCR analysis.

Statistical analysis
All experiments were repeated three times. Data were shown as mean AE S.D. and analysed with a two-tailed Student's t-test. P ≤ 0.05 was considered statistically significant.

T2DM-related ceRNA network (DMCN)
We constructed DMCN that included 98 genes, 86 microRNAs and 167 lncRNAs (Fig. S1). All 351 nodes formed a giant component (Fig. 1). A node's degree is defined as the total number of edges connecting it in the network (Fig. S2). Among the top 10 microRNAs with the highest degree, hsa-mir-21 was at top of the list (degree = 35), indicating 35 edges (35 targets) which directly connect to miR-21. As for gene VEGFA (degree = 8) and lncRNA MIR22HG (degree = 18), eight edges and 18 edges are directly connected to VEGFA and MIR22HG, respectively. Node with the higher degree means that this node is highly connected and likely to be 'hub' node. Yellow, red and blue nodes represented type 2 diabetes-related genes, microRNAs and lncRNAs, respectively (Fig. 1). The edges between the nodes represented the interactions between these RNAs. The size of each node represented the degree. The top 10 nodes of genes, microRNAs and lncRNAs with highest degree were shown in Table 1. The highest degree nodes were VEGFA/ESR1, hsa-mir-21 and lncRNA MIR22HG, respectively (Table 1 and Fig. 2A-C).
We also calculated a network feature-node betweenness, defined as the fraction of shortest paths between node pairs that pass through a given node. The betweenness of a microRNA represents its role in connecting mRNAs and lncRNAs to each other (Fig. S2). In our DMCN, the higher betweenness value of a microRNA means that it regulates more potential ceRNA pairs. Our results showed that miR-21, miR-155, miR-30a had not only higher degrees, but also higher betweenness values (Table 1), indicating their role not only in connecting other nodes, but also in the interactions with mRNAs and lncRNAs. However, miR-632 and miR-99a had the highest betweenness value (betweenness value = 1) but relative low degree (Table 1 and Fig. 2D), suggesting their role in connecting genes and lncRNAs rather than being a 'hub' node in a local network.

Gene ontology and KEGG pathway analysis of T2DM genes
To explore biological functions of lncRNAs during the development of diabetes mellitus, we annotated DMCN into GO and KEGG pathways. The GO BP (biological process) analysis showed that the most significant enriched terms were bidirectional regulation of apoptosis and cell death, and others were related to response to external and internal stimulus and cellular metabolism (Fig. 3A). The top 20 most significant enriched pathways included pathways in cancer, insulin signalling pathway, type 2 diabetes and mTOR signalling (Fig. 3B). Among these 20 pathways, half were related to cancer, and others were mainly involved in immune response, energy and substance metabolism, neuronal function.

The components of mTOR subnetwork
In global diabetes-related ce-network, we focused on mTOR subnetwork due to the central role of mTOR in energy metabolism. The components of mTOR subnetwork were shown in Figure 4A. miR-144 and miR-181b and their predicted or validated binding sites with mTOR were shown in Figure 4B and C [14,15]. We detected the interactions between miR-181b-MALAT1 and miR-144-MALAT1 ( Fig. 4D and E). The two interactions sites between miR-181b and NEAT1 were shown in Figure 4F and G.

The ncRNAs-involved mTOR signalling pathway
The mTOR subnetwork extracted from DMCN consisted of mTOR gene, three microRNAs and 15 lncRNAs (Fig. 4A). We focused on miR-181b and miR-144 and lncRNAs MALAT1 and NEAT1 due to their potential relationship with diabetes. To annotate these ncRNAs to mTOR signalling pathway, we mapped MALAT1 and NEAT1 to mTOR signalling pathway according to a study using proteomics with CHART-enriched material (CHARTMS) to identify proteins associated with NEAT1 and MALAT1 in vivo [16]. The results showed that SREBP1, ATG13, LAMTOR3, EIF4G were associated with both MALAT1 and NEAT1. LAMTOR3 is a  member of trimeric protein complex called Ragulator, which is a Rag-interacting complex essential for amino acid signalling to mTORC1 [17] (Fig. 5). In addition, Ragulator-Rag complex serves as a docking site for mTORC1 on lysosomal membrane, and lysosomes play roles in autophagy, energy metabolism and intracellular signalling [18]. Downstream of mTORC1, SREBPs are activated through S6K1 to induce the expression of genes involved in fatty acids and sterols metabolism [19]. Another mTOR signalling downstream effector, ATG13, functions as a positive regulator of ULK1 (an autophagy initiating kinase) to suppress initiation of autophagy during nutrient abundance [18]. EIF4G acts as a scaffolding protein for translation initiation complex, and insulin signalling regulates translation initiation by regulating eIF4G1 [20].

NEAT1 and miR-181b levels in HL7720 cells
After transfection of miR181b mimics into HL7720 cells, we found that miR-181b level increased (Fig. 6A left) while NEAT1 level decreased (Fig. 6B). On the contrary, after transfection of NEAT-siRNA into HL7720 cells, NEAT1 level decreased (Fig. 6A right) while the level of miR-181b was up-regulated compared to cells transfected with NC siRNA (Fig. 6C). RIP assay showed that the relative enrichment of AGO2-immunoprecipitated NEAT1 and miR-181b decreased (Fig. 6D). These results suggest that NAET1 may be one target of miR-181b.

Discussion
Accumulating evidence has demonstrated regulatory role of ncRNAs in metabolic disease, and their potential as novel targets for drug development [21,22]. In this study, we constructed a ceRNA network to identify regulatory role of ncRNAs in diabetes.
DMCN contained 98 genes, 167 lncRNAs and 86 microRNAs. As shown in Table 1, the most connected genes were VEGFA and ESR1, consistent with that allelic variants of VEGFA gene are associated with microvascular complications of diabetes and atherosclerosis [23] and that ESR1 gene was one of the most highly interconnected T2DM network 'hub' genes [24]. LncRNA MIR22HG was at top of the degree list. A study reported that MIR22HG was significantly up-regulated in the endothelial cell in response to hypoxia both in vitro and in vivo [25]. Another study reported that MIR22HG responded to chemical stresses in induced pluripotent stem cells, indicating its role in response to external stimulus [26]. OIP5-AS1 is highly expressed in the nervous system and is important for neurogenesis during development [27]. However, the relationship between other lncRNAs and diabetes mellitus has not been reported.
Degree and betweenness are two major parameters to characterize DMCN. Among the top 10 microRNAs with the highest degree, hsa-mir-21 was at top of the list (degree = 35). miR-21 was reported to be positively correlated with BMI [28]. Moreover, miR-21 has been shown to contribute to cytokine-mediated b-cell dysfunction in type 1 diabetes [29]. miR-155 (degree = 23) played a role in linking adipose tissue dysfunction and the development of obesity-associated disorders including type 2 diabetes [30]. We also calculated a network feature-node betweenness, defined as the fraction of shortest paths between node pairs that pass through a given node [31]. The betweenness of a microRNA represents its role in other nodes that connect to each other. Our results showed that miR-21, miR-155, miR-30a had not only higher degrees, but also higher betweenness values (Table 1), confirming the strong correlation between degree and betweenness. However, miR-632 and miR-99a had the highest betweenness value (betweenness value = 1) but relative low degree, suggesting their role in connecting genes and lncRNAs rather than being a 'hub' node in a local network.
To investigate the functions of T2DM-related genes, we annotated these genes into GO and KEGG analysis. We found that the most significant terms were involved in the regulation of cell death. Apoptosis of b-cells has been demonstrated to be involved in autoimmune T1DM and T2DM [32]. Moreover, apoptosis can lead to pancreatic bcell loss and endothelial cell injury [33], suggesting its significant role in the pathogenesis in DM and its complications.
Due to its central role in energy metabolism, we extracted mTOR subnetwork from the global DMCN and aimed to investigate the regulatory role of ncRNAs in mTOR pathway. As shown in Figure 5A, 15 lncRNAs including MALAT1, NEAT1 and three microRNAs hsa-mir-144-3p, hsa-mir-181b-5p and hsa-mir-16-5p were grouped in this subnetwork. We only focused on lncRNA MALAT1 and NEAT1 considering their potential role in diabetes. mTOR is an atypical serine/threonine protein kinase that integrates various extracellular signals into a  variety of anabolic processes, cell growth and autophagy [34]. Dysregulation in mTOR signalling is implicated in various diseases such as obesity, T2DM, cancer and ageing [35].
Recently, various studies showed that microRNAs affect mTOR pathway by targeting upstream regulators including IGF-R, PI3K and Akt [36]. In our study, we found that miR-144 and miR-181b were involved in the regulation of mTOR signalling by targeting both upstream and downstream effectors. In silico analysis has revealed two miR-144 binding sites in the mTOR 3 0 untranslated region, and in vitro assays showed that mTOR is a direct target of mir-144 [37]. Increased circulating level of miR-144 has been found to inhibit insulin receptor substrate 1(IRS-1) to impair insulin signalling in T2DM [38]. Insulin growth factors (IGF-1 and IGF-2) and their receptors play a significant role in insulin resistance and T2DM, and miR-181b overexpression inhibited PI3K/AKT signalling through targeting IGF-1R [21]. lncRNA metastasis-associated lung adenocarcinoma transcript 1 (MALAT1) was associated with diabetes-induced microvascular dysfunction, as up-regulation of MALAT1 could activate p38/MAPK signalling and regulate retinal endothelial cell function under diabetic condition [39]. A recent study indicated that MALAT1 regulated hyperglycaemia induced inflammatory process in endothelial cells [40]. NEAT1 (nuclear enriched abundant transcript 1) is a nuclear lncRNA that was shown to participate in miR-140 induced adipogenesis [41].
Then, we annotated these ncRNAs to mTOR signalling pathway. We further validated two ceRNA pairs (MALAT1-miR-144-mTOR and NEAT1-miR-181b-mTOR) involved in mTOR pathway and explored their significance in diabetes mellitus. NEAT1 is positively co-expressed with mLST8, and mLST8 expression was decreased by mir-181b mimics, indicating potential communication between NEAT1 and mTOR signaling pathway medicated by microRNA-target sites. These findings establish a novel connection between lncRNA and microRNA in diabetes mellitus.
In summary, we constructed a T2DM-related ceRNA network (DMCN) to explore the biological functions of lncRNAs during the development of diabetes mellitus. A mTOR-centred ceRNA subnetwork was extracted, and we experimentally validated that NEAT1 was potentially communicated with mTOR signalling pathway target protein mLST8 via the association with miR-181b. These findings provide significant insight into lncRNA regulatory network in T2DM.

Conflict of interest
None.

Supporting information
Additional Supporting Information may be found online in the supporting information tab for this article: