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

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 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.

Materials and methods

Collection of human T2DM‐related genes

We downloaded human T2DM‐related genes from The Text‐mined Hypertension, Obesity and Diabetes candidate gene database (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 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% CO2.

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).

Western blot analysis

Proteins from cells were separated by sodium dodecyl sulphate‐polyacrylamide gel electrophoresis (10–12%). After blocking, the membranes were probed with antibody for mLST8 (ab25974) and β‐actin (ab8227) (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 β‐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 ± S.D. and analysed with a two‐tailed Student's t‐test. P ≤ 0.05 was considered statistically significant.

Results

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).

Figure 1

ceRNA network related to type 2 diabetes mellitus. Nodes in yellow, red and blue represented genes, microRNAs and lncRNAs, respectively.

Table 1

The top 10 degree list of miRNAs, genes and lncRNAs

miRNA	degree	Gene	degree	LncRNA	degree
hsa‐miR‐21	35	VEGFA	8	MIR22HG	18
hsa‐miR‐155	23	ESR1	8	HCG18	16
hsa‐miR‐30a	23	TP53	7	OIP5‐AS1	14
hsa‐miR‐31	19	EP300	6	TUG1	14
hsa‐miR‐183	17	FOXO1	5	TRAF3IP2‐AS1	13
hsa‐miR‐26b	17	KRAS	4	TP73‐AS1	12
hsa‐miR‐15a	16	HMGA2	3	ZNF518A	10
hsa‐miR‐182	16	HIF1A	3	HOXD‐AS1	10
hsa‐miR‐29a	16	PPARG	3	TTC28‐AS1	9
hsa‐miR‐15b	15	EDNRA	3	RP11‐65F13.2	8

Figure 2

(A, B and C) showed the degree distribution of genes, lncRNAs and microRNAs, respectively. Horizontal axis represented the degree values, and vertical axis represented the number of nodes with a certain degree value. (D). Horizontal axis represented the betweenness value of microRNAs.

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.

Figure 3

(A) Gene ontology (GO) analysis of T2DM‐related genes. (B) KEGG pathway analysis of T2DM‐related genes. Only top 20 terms or pathways were presented.

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.

Figure 4

(A) mTOR subnetwork. (B) Predicted interaction of mir‐144 and mTOR (by miRanda). (C) Predicted interaction of mir‐181b‐5p and mTOR (by miRanda). (D) Interactions between mir‐181b and MALAT1. (E) Interactions between mir‐144 and MALAT1. (F) The first interaction site between mir‐181b and NEAT1. (G) The second interaction site between mir‐181b and NEAT1. All data were downloaded from Starbase (v2.0) and miRTarBase (v4.5).

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.

Figure 5

The ncRNAs involved in mTOR signalling pathway. MicroRNAs were labelled in red, and lncRNAs were labelled in blue. MALAT1 and NEAT1 were mapped to mTOR signalling pathway according to a study using proteomics with CHART‐enriched material (CHARTMS) to identify proteins associated with NEAT1 and MALAT1 16. SREBP1, ATG13, LAMTOR3, EIF4G were associated with both MALAT1 and NEAT1.

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.

Figure 6

NEAT1 and miR‐181b levels in HL7702 cells. (A) Relative miR‐181b and NEAT1 levels in HL7702 cells transfected with negative control (NC) mimics versus miR‐181b mimics, or negative control (NC) siRNA versus NEAT1 siRNA. (B) The expression of NEAT1 in HL7702 cells transfected with miR‐181b mimics. (C) The expression of miR‐181b in HL7702 cells transfected with NEAT1‐siRNA. For comparison, expression level in cells transfected with NC mimics or NC siRNA was set to 1. All experiments were repeated three times. Each bar represents the mean ± S.D. of three independent experiments. *P < 0.05. (D) RIP assay showing the relative enrichment of AGO2‐immunoprecipitated NEAT1 and mir‐181b. anti‐Ago2: samples incubated with Ago antibody; Input: cell lysate.

The interaction between NEAT1 and miR‐181b by mTOR signalling pathway

To investigate whether mTOR signalling pathway plays an important role in ceRNA of NEAT1 and miR‐181b, we searched the TargetScan database 7.1 and found that mLST8, a subunit of both mTOR complex 1 (mTORC1) and complex 2 (mTORC2), was a potential target of miR‐181b. As expected, miR‐181b level increased after transfection with miR‐181b mimics (Fig. 7A left), while NEAT1 level decreased after transfection with NEAT‐siRNA (Fig. 7A right). Western blot analysis showed that the expression of mLST8 decreased after transfection with miR‐181b mimics (Fig. 7B). Luciferase assay showed that miR‐181b mimics reduced mLST8‐WT luciferase activity (Fig. 7C). Furthermore, NEAT1‐siRNA decreased protein level of mLST8 (Fig. 7D). Similarly, NEAT1‐siRNA reduced mLST8‐WT luciferase activity (Fig. 7E). These results demonstrate that NEAT1‐miR‐181b‐mLST8 ceRNA pair regulates mTOR signalling pathway.

Figure 7

The interaction between NEAT1 and miR‐181b. (A) Relative miR‐181b and NEAT1 levels in cells transfected with negative control (NC) mimics versus miR‐181b mimics, or negative control (NC) siRNA versus NEAT1 siRNA. The expression of mLST8 protein in cells transfected with NC mimics versus miR‐181b mimics (B) or NC siRNA versus NEAT1 siRNA (D). The luciferase activity of mLST8 promoter in cells transfected with NC mimics versus miR‐181b mimics (C) or NC siRNA versus NEAT1‐siRNA (E). All experiments were repeated three times. Each bar represents the mean ± S.D. of three to six independent experiments. *P < 0.05.

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 β‐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 β‐cells has been demonstrated to be involved in autoimmune T1DM and T2DM 32. Moreover, apoptosis can lead to pancreatic β‐cell 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′ 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.

Figure S1 Schematic figure of construction of T2DM related ceRNA network (DMCN).

Figure S2 Example network with five nodes and five edges.

Click here for additional data file.