Dynamic co-expression network analysis of lncRNAs and mRNAs associated with venous congestion

Introduction

Despite advances in pharmacological and mechanical therapies, acute decompensated heart failure (ADHF) remains one of the most frequently encountered cardiovascular diseases and its post-discharge mortality rates can reach 8.1 and 20.5% within 3 months and 1 year, respectively (1–3). Over the past decades, it has become evident that acute heart failure (HF) and renal failure (RF) are associated with systemic inflammatory activation (4). It has become increasingly recognized that venous congestion and volume overload are involved in the cardiorenal syndrome (4,5). Certain animal and human studies indicate that congestion may lead to venous endothelial activation with peripheral synthesis and release of proinflammatory mediators, including cyclooxygenase-2, inducible nitric oxide synthase, interleukin-6 (IL-6) and tumor necrosis factor-α (6–8). These suggest that the venous endothelium is in a proinflammatory state in times of venous congestion (4). To investigate global mRNA changes during peripheral venous congestion, Colombo et al (5) simulated peripheral venous congestion and analyzed mRNA by using the Affymetrix HG-U133 Plus 2.0 microarray.

Long noncoding RNAs (lncRNAs) have become a research hotspot for numerous diseases. A number of lncRNAs have been demonstrated to have important and diverse functions (9,10). LncRNA-associated dysfunction has been demonstrated to be important in cancer (11), cardiovascular diseases (12), and neurodegeneration diseases (11). Particularly, it is becoming evident that lncRNA may be involved in cardiovascular diseases. For example, the myocardial infarction-associated transcript lncRNA is associated with myocardial infarction (13). Another study identified 15 lncRNAs modulated in the heart of mice subjected to aortic constriction (14). However, global analysis of lncRNA associated with peripheral venous congestion is required and the potential underlying regulatory mechanisms remain unclear, due to the limited RNA sequencing (RNA-Seq) data. Thus, the present study re-annotated an Affymetrix microarray associated with peripheral venous congestion, then constructed a dynamic lncRNA-mRNA co-expression network (5). Following functional analysis of this network, it was demonstrated, although the genes (mRNAs) and lncRNAs were different, certain venous congestion-associated GO terms, including ion channel activity, were identified. The current study also identified certain lncRNA dynamically regulated pathways, including dilated cardiomyopathy and the calcium signaling pathway, in which the involvement of lncRNAs persistently occurred from normal and peripheral venous congestion conditions. To the best of our knowledge, the present study was the first to analyze the dynamic lncRNA-associated mechanism of peripheral venous congestion and provide insights into the understanding of the functional mechanism of peripheral venous congestion and lncRNAs.

Materials and methods

Microarray data

The microarray data set GSE38783 was accessed using the Gene Expression Omnibus database (ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE38783). This data set was obtained from 12 healthy subjects using the Affymetrix HG-U133 Plus 2.0 microarray (Affymetrix, Inc., Santa Clara, CA, USA) (5). Venous arm pressure was increased to 30 mmHg above the baseline level by inflating a tourniquet cuff around the dominant arm (test arm). Then endothelial cells were obtained from blood samples from the test and control arm (lacking an inflated cuff) before and after 75 min of venous congestion (5).

Functional re-annotation of lncRNAs

To re-annotate micro-array data, a non-coding RNA function annotation server (ncFANs) was used to re-annotate the probes of the HG-U133 Plus 2.0 array as described on the ncFANs website (15). Then each probe was converted into gene Ensembl Gene IDs. There were 3495 lncRNAs re-annotated. When a gene matched more than one probe, the expression value of this mRNA or lncRNA was computed by taking the average expression value of all the corresponding probes.

Construction of dynamic lncRNA-mRNA co-expression network

Pearson's correlation coefficient (PCC) was calculated between expressed values of each lncRNA-mRNA pair in normal samples and venous congestion samples. The co-expressed lncRNA-mRNA pairs with PCC>0.99 or <−0.99 and P<0.01 were selected. To construct the dynamic lncRNA-mRNA co-expression network, two lncRNA-mRNA co-expression networks were initially constructed based on lncRNA-mRNA co-expression associations in normal samples and venous congestion samples. Then, the different sets of the above two networks were calculated and the following two novel networks identified: i) 'Lost' network, in which lncRNA-mRNA co-expression pairs only appeared in normal samples and not in the venous congestion samples; and ii) 'obtained' network in which co-expression pairs only appeared in venous congestion samples and not in normal samples. Finally, the 'lost' and 'obtained' networks were combined together to obtain the final dynamic lncRNA-mRNA co-expression network (the edges of 'lost' and 'obtained' network were added). The procedure is illustrated in Fig. 1.

Figure 1 Flow chart of construction of the dynamic lncRNA-mRNA co-expression network. The 'normal' and 'venous congestion' samples represented the arms before and after the induced hypertension treatment, respectively. PCC, Pearson's correlation coefficient; lncRNA, long non-coding RNA.

Functional enrichment

To determine potential functions of lncRNAs, co-expressed genes for the lncRNAs were obtained then GO enrichment was implemented to identify the affected GO terms of lncRNA co-expressed genes using the Database for Annotation, Visualization and Integrated Discovery 6.7 (16). P<0.01 was considered to indicate a statistically significant difference. Pathway annotation used the Kyoto Encyclopedia of Genes and Genomes (KEGG) search and color tools.

Results

Co-expression networks based on normal and venous congestion samples

To construct the dynamic lncRNA-mRNA co-expression network, two lncRNA-mRNA co-expression networks based on lncRNA-mRNA co-expression associations in normal samples and venous congestion samples were initially constructed. In the 'venous congestion' co-expression network, there were 56 lncRNA nodes, 61 mRNA nodes and 64 edges, whereas the 'normal' network contained 364 lncRNA-gene co-expressing associations composed of 175 lncRNAs and 275 mRNAs (Fig. 2A and B).

Figure 2 lncRNA-mRNA co-expression network in (A) venous congestion samples and (B) normal samples. lncRNA, long noncoding RNA.

'Lost' and 'obtained' network

The different sets of the above two networks were then calculated and the following two novel networks were created: i) The 'lost' network contained 350 edges composed by 172 lncRNA and 266 mRNA nodes and ii) the 'obtained' network contained 56 edges composed of 52 lncRNA and 54 mRNA nodes (Fig. 3A and B). To investigate the 'lost' and 'obtained' functions, GO enrichment was performed on the two networks (Tables I and II). Certain terms (including ion channel activity) were significantly enriched (P<0.01), although the mRNAs in the two networks differed.

Figure 3 (A) 'Obtained' lncRNA-mRNA co-expression network and (B) the 'lost' lncRNA-mRNA co-expression network. lncRNA, long noncoding RNA.

Table I GO terms enriched by 'obtained' genes.

Table I

GO terms enriched by 'obtained' genes.

GO term	P-value
GO:0005216 ion channel activity	0.002033
GO:0022838 substrate specific channel activity	0.002323
GO:0015267 channel activity	0.002700
GO:0022803 passive transmembrane transporter activity	0.002728
GO:0005261 cation channel activity	0.003950
GO:0022836 gated channel activity	0.006033
GO:0046873 metal ion transmembrane transporter activity	0.007346

[i] GO, gene ontology.

Table II GO terms enriched by 'lost' genes.

Table II

GO terms enriched by 'lost' genes.

GO terms	P-value
GO:0015267 channel activity	0.001265
GO:0022803 passive transmembrane transporter activity	0.001296
GO:0043565 sequence-specific DNA binding	0.001758
GO:0003700 transcription factor activity	0.002589
GO:0005261 cation channel activity	0.002731
GO:0046873 metal ion transmembrane transporter activity	0.003682
GO:0005216 ion channel activity	0.005109
GO:0043176 amine binding	0.005961
GO:0022838 substrate specific channel activity	0.006561
GO:0000287 magnesium ion binding	0.007641

[i] GO, gene ontology.

Dynamic lncRNA-mRNA co-expression network

The 'lost' and 'obtained' networks were combined to obtain the final dynamic lncRNA-mRNA co-expression network (Fig. 4). This network contained 524 nodes (211 lncRNA nodes and 313 mRNA nodes) and 406 edges. The giant component was formed of 172 nodes. The top 10 highest degree nodes are demonstrated in Tables III and IV. The current study demonstrated that certain lncRNAs exhibited co-expressed mRNAs in the normal and venous congestion conditions, indicating that these lncRNAs may dynamically regulate certain biological process during the development of venous congestion. Certain mRNAs also had co-expressed lncRNA partners in both conditions, indicating that they were dynamically regulated by different lncRNAs during the development of venous congestion.

Figure 4 Dynamic lncRNA-mRNA co-expression network. The triangle nodes denote lncRNA and rectangle nodes denote mRNAs. The red and green edge denotes co-expression associations between a lncRNA and an mRNA were lost in normal conditions or 'newly' obtained in venous congestion conditions, respectively. lncRNA, long noncoding RNA.

Table III Top 10 long non-coding RNA nodes with the highest degree.

Table III

Top 10 long non-coding RNA nodes with the highest degree.

Ensembl ID	Name	Degree
ENSG00000196273	LINC00523	14
ENSG00000225946	RP11-395B7.2	11
ENSG00000267649	CTD-2587H24.10	11
ENSG00000233101	HOXB-AS3	8
ENSG00000260394	LA16c-313D11.9	8
ENSG00000224184	AC096559.1	8
ENSG00000270168	LA16c-380H5.4	7
ENSG00000272933	RP11-47A8.5	7
ENSG00000228008	CTD-2330K9.3	7
ENSG00000259905	PWRN1	7

[i] LINC00523, long intergenic non-protein coding RNA 523; HOB-AS3, HOXB cluster antisense RNA 3; PWRN1, Prader-Willi region non-protein coding RNA 1.

Table IV Top 10 mRNA nodes with the highest degree.

Table IV

Top 10 mRNA nodes with the highest degree.

Ensembl ID	Name	Degree
ENSG00000176884	GRIN1	7
ENSG00000183250	C21orf67	5
ENSG00000167011	NAT16	5
ENSG00000108342	CSF3	4
ENSG00000086205	FOLH1	4
ENSG00000186868	MAPT	4
ENSG00000123989	CHPF	4
ENSG00000178498	DTX3	3
ENSG00000141150	RASL10B	3
ENSG00000154252	GAL3ST2	3

[i] GRIN1, glutamate ionotropic receptor N-methyl D-aspartate 1; C21orf67, long intergenic non-protein coding RNA 1547; NAT16, N-acetyltransferase 16; CSF3, colony stimulating factor 3; FOLH1, folate hydrolase (prostate-specific membrane antigen) 1; MAPT, microtubule associated protein tau; CHPF, chondroitin polymerizing factor; DTX3, deltex 3, E3 ubiquitin ligase; RASL10B, RAS like family 10 member B; GAL3ST2, galactose-3-O-sulfotransferase 2.

Detection of dynamically alternative pathways potentially regulated by lncRNAs

To investigate the dynamic alternative biological functions of lncRNAs during the development of venous congestion, lncRNA co-expressing mRNAs in the 'lost' and 'obtained' networks were mapped into KEGG pathways. A dynamic alternative pathway network was constructed (Fig. 5A). The present study demonstrated that i) There were certain pathways directly associated with venous congestion (hsa05414, dilated cardiomyopathy; these pathways are illustrated with red labels in Fig. 5A); and ii) certain pathways (hsa05414, dilated cardiomyopathy; hsa04020, calcium signaling pathway), were annotated by the 'lost' and 'obtained' mRNAs. For detailed investigation of the dynamic changes, the 'lost' and 'obtained' mRNAs were mapped into proteins in the calcium signaling pathway (Fig. 5B; pink rectangle and blue pentagram represent proteins which were potentially regulated by lncRNAs in normal and venous congestion samples).

Figure 5 (A) Dynamically alternative pathway network. (B) Calcium signaling pathway; the pink rectangle and blue pentagram denoted proteins which were annotated by 'lost' and 'obtained' mRNAs, meaning that they were dynamically regulated lncRNAs in normal and venous congestion condition, respectively. lncRNA, long non-coding RNA; NCX: Na+/Ca2+ exchanger-3; PMCA, plasma membrane calcium ATPase; GPCR, G protein-coupled receptor; ORAI, calcium release-activated calcium modulator; CAV, certain voltage-activated Ca2+ channels; ROC, nicotinic acetylcholine receptor α-7; Gq, guanine nucleotide-binding protein G(q) subunit α; PTK, receptor tyrosine-protein kinase erbB-2; TCR, T cell receptor α-chain V region; BCR, immunoglobulin heavy chain; CD38, ADP-ribosyl cyclase 1; Gs, guanine nucleotide-binding protein G(s) subunit α; ADCY, adenylate cyclase 1; PLC, phosphatidylinositol phospholipase C; SPHK, sphingosine kinase; PKA, protein kinase A; STIM, stromal interaction molecule; PLN, phospholamban; SERCA, Ca2+ transporting ATPase, sarcoplasmic/endoplasmic reticulum; RYR, ryanodine receptor; IP3R, inositol 1,4,5-triphosphate receptor type; CALM, calmodulin; VDAC, voltage-dependent anion channel protein; TnC, troponin C, slow skeletal and cardiac muscles; MLCK, myosin-light-chain kinase; PHK, phosphorylase kinase; CaN, serine/threonine-protein phosphatase 2B; CAMK, calcium/calmodulin-dependent protein kinase; NOS, nitric-oxide synthase; PDE1, phosphodiesterase 1; FAK2, focal adhesion kinase 2; IP33K, 1D-myo-inositol-triphosphate 3-kinase; PKC, protein kinase C.

Discussion

In recent years, the significant functional molecular mechanisms of lncRNAs have been recognized, particularly in cardiovascular diseases. RNA-Seq is currently a prevalent technique to profile lncRNAs. However, due to relatively high cost of RNA-Seq techniques, publicly available RNA-Seq data are limited compared with microarray expression profile data. Furthermore, RNA-Seq data lack high sample numbers, in contrast with microarray expression profile data that often contains dozens to hundreds of pair-matched samples. Thus, the current study adopted a re-annotation method to identify lncRNAs associated with venous congestion (15).

Furthermore, increasing evidence has demonstrated that lncRNAs may be important for regulating gene expression. The functions of lncRNAs are predominantly performed by their secondary structures, however, this is difficult to decipher (17). Due to considerable challenges in investigating the lncRNA functions, a co-expression-based method was developed, in which lncRNA functions were predicted based on the functions of their co-expressed protein-coding genes (17), as genes with similar expression patterns under multiple conditions tend to be involved in the same pathways (18). The co-expressed protein-coding genes are potentially regulated by the corresponding lncRNAs.

Thus, in the current study, an Affymetrix microarray associated with peripheral venous congestion was re-annotated and normal and venous congestion co-expression networks were constructed (Fig. 2A and B). Based on these networks, another two networks were created; the 'lost' and 'obtained' networks. The 'lost' network represented the lncRNA-mRNA co-expression associations that occurred in the normal condition but disappeared in the venous congestion condition. By contrast, the 'obtained' network represented the novel lncRNA-mRNA co-expression associations that occurred in the venous congestion conditions, but were not present in the normal conditions. Then, GO enrichment was performed on the protein-coding mRNAs of the two networks. It was assumed that the GO terms of 'obtained' network represented 'novel' biological processes that were potentially influenced by lncRNAs, whereas the GO terms of the 'lost' network represented dysregulated biological processes that maintained the normal functions of the human body. Notably, numerous terms closely associated with venous congestion, including ion channel activity (19), were significantly enriched, suggesting that lncRNAs may be involved in similar biological processes by regulating the same mRNAs or different mRNAs in the same functional term (Tables I and II).

To decipher the dynamically regulated process during the development of venous congestion, the 'lost' and 'obtained' networks were combined to form the dynamic lncRNA-mRNA co-expression network (Fig. 4). The nodes with the highest degree are demonstrated in Tables III and IV. The most-connected mRNAs nodes were glutamate ionotropic receptor N-methyl D-aspartate 1 (GRIN1). It has been previously reported that the release of proinflammatory cytokine IL-1β in venous congestion condition enhanced GRIN1 responses through IL-1 receptor activation (20).

The current study also demonstrated that certain lncRNAs were 'dynamically' co-expressed with different mRNAs, and certain mRNAs were 'dynamically' co-expressed with different lncRNAs. To further investigate the dynamic biological pathways, a dynamic alternative pathway network was constructed (Fig. 5A). This demonstrated that there were certain pathways (for example, the pathway of 'dilated cardiomyopathy') directly associated with venous congestion or cardiovascular diseases (21). Furthermore, certain pathways, for example, hsa04020: calcium signaling pathway, were annotated by the 'lost' and 'obtained' mRNAs, suggesting that lncRNA dynamically regulate this pathway. The 'lost' and 'obtained' mRNAs and their corresponding lncRNAs were annotated to the calcium signaling pathway (Fig. 5B). G protein-coupled receptors (GPCRs) were potentially regulated in the normal and venous congestion conditions, but by different lncRNAs (LINC01210 in venous congestion and RP11-435O5.5 in normal samples; Fig. 5B). A previous novel study demonstrated the importance of GPCR signaling partners in endothelial cell nitric-oxide synthase (eNOS) function and the fundamental implications for vascular disorders involving dysregulated eNOS (22). In the normal samples, lncRNA-LINC00523 were co-expressed with PMCA, indicating the potential regulated association of LINC00523 with PMCA. However, in the venous congestion, this association was changed. LINC00523 were co-expressed with GPCRs. PMCA is a calcium pump, with a pivotal function in maintaining intracellular calcium homeostasis (23,24). Therefore, this observation indicated that, although further experimental validation is required, the dynamic regulated switch of LINC00523 from PMCA to GPCR may be involved in the damage to steady state intracellular calcium levels.

In the current study, the underlying role of lncRNAs in the development of venous congestion was investigated by the construction of a dynamic lncRNA-mRNA co-expression network. Functional analysis identified several dynamic alternative pathways, including the calcium signaling pathway. Notably, it was suggested that the dynamically regulated switch of LINC00523 from PMCA to GPCR may be involved in damage to the intracellular calcium homeostasis, leading to the development of venous congestion. In brief, the current study demonstrated a potential novel mechanism of lncRNA function during venous congestion.

Acknowledgments

This work was supported in part by the Health and Family Planning Commission of Heilongjiang Province [Grant Nos. 2012-357] and Health and Family Planning Commission of Heilongjiang Province [Grant Nos. 2013-110].

References