lncRNA expression profiles and associated ceRNA network analyses in epicardial adipose tissue of patients with coronary artery disease

Epicardial adipose tissue (EAT) contributes to the pathophysiological process of coronary artery disease (CAD). The expression profiles of long non-coding RNAs (lncRNA) in EAT of patients with CAD have not been well characterized. We conducted high-throughput RNA sequencing to analyze the expression profiles of lncRNA in EAT of patients with CAD compared to patients without CAD. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were executed to investigate the principal functions of the significantly dysregulated mRNAs. We confirmed a dysregulated intergenic lncRNA (lincRNA) (LINC00968) by real-time quantitative PCR (RT-qPCR). Subsequently, we constructed a ceRNA network associated with LINC00968, which included 49 mRNAs. Compared with the control group, lncRNAs and genes of EAT in CAD were characterized as metabolic active and pro-inflammatory profiles. The sequencing analysis detected 2539 known and 1719 novel lncRNAs. Then, we depicted both lncRNA and gene signatures of EAT in CAD, featuring dysregulation of genes involved in metabolism, nuclear receptor transcriptional activity, antigen presentation, chemokine signaling, and inflammation. Finally, we identified a ceRNA network as candidate modulator in EAT and its potential role in CAD. We showed the expression profiles of specific EAT lncRNA and mRNA in CAD, and a selected non-coding associated ceRNA regulatory network, which taken together, may contribute to a better understanding of CAD mechanism and provide potential therapeutic targets. Trial registration Chinese Clinical Trial Registry, No. ChiCTR1900024782.


Epicardial adipose tissue (EAT) contributes to the pathophysiological process of coronary artery disease (CAD). The expression profiles of long non-coding RNAs (lncRNA) in EAT of patients with CAD have not been well characterized. We conducted high-throughput RNA sequencing to analyze the expression profiles of lncRNA in EAT of patients with CAD compared to patients without CAD. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses
were executed to investigate the principal functions of the significantly dysregulated mRNAs. We confirmed a dysregulated intergenic lncRNA (lincRNA) (LINC00968) by real-time quantitative PCR (RT-qPCR). Subsequently, we constructed a ceRNA network associated with LINC00968, which included 49 mRNAs. Compared with the control group, lncRNAs and genes of EAT in CAD were characterized as metabolic active and pro-inflammatory profiles. The sequencing analysis detected 2539 known and 1719 novel lncRNAs. Then, we depicted both lncRNA and gene signatures of EAT in CAD, featuring dysregulation of genes involved in metabolism, nuclear receptor transcriptional activity, antigen presentation, chemokine signaling, and inflammation. Finally, we identified a ceRNA network as candidate modulator in EAT and its potential role in CAD. We showed the expression profiles of specific EAT lncRNA and mRNA in CAD, and a selected non-coding associated ceRNA regulatory network, which taken together, may contribute to a better understanding of CAD mechanism and provide potential therapeutic targets.

Trial registration Chinese Clinical Trial Registry, No. ChiCTR1900024782.
Cardiovascular disease is the leading cause of death in the world. A large number of studies have confirmed that obesity is closely related to CAD. Recent studies have found in addition to the accumulation of peripheral fat, cardiac fat distribution is also a significant risk factor for cardiovascular disease 1,2 . EAT is visceral adipose tissue located between the heart and the visceral pericardium 3 . EAT is anatomically adjacent to the coronary arteries and myocardium, and directly contact with the coronary arteries without fascia separation. EAT can secrete a large amount of adipokines, affecting the physiological functions and pathophysiological processes of the myocardium and coronary arteries 4,5 . Various clinical studies have confirmed that EAT is closely related to the occurrence and development of coronary heart disease, but little research has been conducted on how EAT mediates the occurrence and development of CAD 6 .
LncRNA is a type of RNA with more than 200 nucleotides in length and lacks a complete open reading frame with little or no protein-coding ability 7 . lncRNA can regulate gene expression through different mechanisms of action, such as RNA-RNA base pairing, RNA-protein interaction, RNA-DNA interaction, etc. 8,9 . lncRNA has been reported to play an essential role in the regulation of cardiovascular physiological and pathophysiological processes such as heart development, heart failure, and CAD 10 . CeRNAs were described as groups of non-coding

Results
Study population. We enrolled 19 CAD patients who underwent coronary artery bypass graft surgery (CAD group), and 24 non-CAD patients who underwent cardiac valve replacement with negative clinical and evidence of CAD (CTRL group), of which 5 pairs of subjects' EAT samples were used for high-throughput RNA Sequencing, the other samples were used for verifying the sequencing results. Compared with CTRL, CAD patients has higher proportion in hypertension and hyperlipidemia, while age, male sex, body mass index (BMI), abdominal circumference (AC), blood lipids, and left heart function showed no difference between two groups. Characteristics of the study population were provided in Table 1. We selected 10 subjects from the two groups matched in age, gender, BMI, abdominal circumference (AC), and RNA-seq were performed. The main clinical characteristics of those 10 subjects were summarized in Table 2. The other samples were used as an expansion for verification. Patients in both groups were routinely treated with heparin anticoagulation.
Identification of differentially expressed lncRNA and mRNA. To identify lncRNAs and mRNAs expression patterns that characterize CAD, we used the RNA-seq method. We found that 844 EAT-specific  Table 3. To show the distribution of all differentially expressed RNA in the dimensions of -log (FDR) and logFC, we used the volcano plot ( Fig. 1).   www.nature.com/scientificreports/

GO and pathway analysis of differentially expressed genes (DEGs).
To provide a framework for the interpretation of our results, we then functionally clustered significant biological pathways by using GO and KEGG pathway analysis 12 . GO analysis showed changes in the biological processes (BP) of DEGs were enriched in the extracellular structure organization, collagen metabolic process, etc., as was shown in Fig Figure 2C showed that the changes in molecular functions (MF) of DEGs were manifested in the functions of actin binding, platelet-derived growth factor binding, etc. Figure 2D showed the results of KEGG pathway analysis.
Construction and analysis of the LINC00968-related ceRNA network. LncRNAs can act as ceR-NAs to interfere with the function of specific miRNAs, thereby releasing the mRNAs affected by the miRNAs which can be predicted by softwares such as miRanda etc. We built the ceRNA network on the basis of the lncRNA and mRNA expression profiles in patients with CAD. Among those DElncRNAs, 28 lncRNAs were identified as lincRNAs and we found that LINC00968 level was more than 3-fold greater in CAD patients confirmed by RT-qPCR (Fig. 3A). Through logistic analyses, we found that LINC00968 level was an independent predictor of CAD when adjusted for possible confounding factors including sex, systolic blood pressure, hypertension and hyperlipidemia (Table 4). Furthermore, the ROC curve was built to assess the sensitivity and specificity of these lincRNAs as predictors for CAD, and the areas under the curve (AUC) were 0.8158 (p < 0.05) (Fig. 3B). ROC curve of LINC00968 relative expression with a cut off value of 2.27 achieved a sensitivity of 57.1%, and specificity of 100.0%. Then, we selected LINC00968 to build the network. The LINC00968-related ceRNA network was predicted based on the potential number of MREs shared by mRNAs and LINC00968, which includes 49 mRNAs. The network was shown in Fig. 4.

PPI network construction and key mRNAs verification.
To further identify key mRNAs in the LINC00968-related ceRNA network, mRNAs of the ceRNA network were submitted to the STRING database. PPI network was constructed with a comprehensive score greater than 0.4 (Fig. 5A), and the most important module was obtained using MCODE in Cytoscape. The key mRNAs were selected from the PPI network using

Discussion
EAT is a type of visceral adipose tissue located between the myocardium and the pericardium. The lineage tracing study confirms that EAT is differentiated from epicardial cells 13 . In addition to acting as a buffer zone around the heart, EAT also acts as an energy reservoir for myocardial tissue to provide free fatty acids for cardiomyocyte metabolism 14 . Under pathological condition, EAT exhibits as metabolically active endocrine and paracrine organ that secretes a large amount of anti-inflammatory and pro-inflammatory adipokines, affecting the dynamic balance of glycolipid metabolism, cardiac structure and function [15][16][17] . EAT is directly nourished by coronary arteries, and there is no fascia isolation between EAT and myocardial tissue, which provides an important anatomical basis for the interaction among them 18 . Several evidence-based medical studies have confirmed that EAT is closely related to the occurrence and development of CAD 19,20 , and may be associated with the increase in the incidence of long-term adverse events. However, thus far, most of these studies are observational. Mechanism research and prospective studies are sparse. It is, therefore, necessary to understand whether changes in EAT physiology are the cause of atherosclerosis or just an epiphenomenon of metabolic syndrome. Here, we provide novel information that describes the molecular identity (mRNAs and lncRNAs) of EAT in CAD and non-CAD by performing a comprehensive analysis of the expression profiles. Through GO and KEGG analyses, we further analyzed the pathways and functions in which the DEGs are involved. We present evidence that, compared with CTRL, EAT of CAD shows enhanced inflammation and promotion of lipid metabolic disorder.
Intergenic lncRNAs (lincRNAs) are defined as non-coding RNAs longer than 200 nucleotides without overlapping protein-coding genes, which are distinguished from the general lncRNA transcripts, as a lot of lncRNAs share sequence with coding loci 21 . LincRNAs are important regulators of gene expression and often exhibit remarkable tissue-specificity 22 . One of the potential regulative mechanisms of lincRNAs is to act as ceRNAs. miRNAs can bind to mRNAs through MREs, causing mRNAs degradation or inhibiting its translation. When lincRNAs and mRNAs share the same MREs, they form a relationship of competition for the same type of miR-NAs, thus preventing mRNAs degradation, and in such manner lincRNA indirectly regulates the expression level of mRNA 23,24 . The ceRNA mechanism shows a novel way of studying lncRNA-mRNA-miRNA cross talk and has received wide attention. Some studies have showed and highlighted the core function of ceRNA network across CVD. lncRNA HIF1A-AS1 acts as a ceRNA of miR-204 and elevates SOCS2 expression, which contributes to ventricular remodeling after myocardial ischemia/reperfusion injury 25 . Another research by Liang et al. 26 identified the lncRNA 2810403D21Rik/Mirf can act as a ceRNA of miRNA26a in ischemic myocardial injury. The lncRNA TNK-AS1 regulates smooth muscle cell proliferation and migration via acting as a ceRNA of miR-150-5p 27 . The uncovering of the ceRNA mechanism may provide a new layer for exploring the physiopathology process and shed new light on cardiovascular research. www.nature.com/scientificreports/ In this study, we observed an adipose tissue-specific lincRNA of high expression, namely LINC00968, which appeared to be significantly dysregulated in CAD. The ROC analysis suggests that LINC00968 may be used as a potential target for further study. Although LINC00968 was previously reported to be connected to the function of endothelial cells 28 , the mechanism of its increase in an adipose tissue-specific manner is still unclear. Here, we applied the bioinformatics method to construct the LINC00968-related ceRNA regulatory network, which provides new ideas for exploring the biological functions of LINC00968. From the constructed LIN00968 related ceRNA network, we found several mRNAs which have been previously reported to play a role in atherosclerosis, such as CXCL12 29 , TAZ 30 , PDLIM5 31 , GNB3 32 , LMNA 33 , TRPM6 34 . Similarly, we used bioinformatics to predict the top five key mRNAs in the LINC00968-related ceRNA network, including GNG13, GNB3, CXCL12, ADCY3 and LMNA, among which CXCL12 and LMNA were significantly dysregulated in CAD according to verification results. LINC00968 may interact with specific miRNAs and cause the dysregulation of those genes. Interestingly, in the network, we found that lncRNA H19 was significantly down regulated, while LncRNA H19 is a significant CAD-associated lncRNA 35 . These results might help us to study EAT's role in CAD from a new perspective. In the future, we will conduct an in-depth study of the regulatory mechanisms underlying the LINC00968 ceRNA network.  www.nature.com/scientificreports/  www.nature.com/scientificreports/

Limitations
The main limitation of the study could be the limited amount of total adipose tissue obtained from human due to ethical restrictions. In addition, the present study stopped short of validating the functions of LINC00968 and further exploring the related mechanisms because of the limitations in methodology of culturing human epicardial adipose tissue cells. Further studies are required to gain insight into the pathogenesis.

Conclusion
In conclusion, our study has identified expression profiles of lncRNAs in EAT that potentially predict the progression of CAD by using RNA-seq. Among them, LINC00968 was identified and its related ceRNA network may play a pivotal role in the process of CAD pathophysiology.

Materials and methods
Patients, samples, and ethics. 43 non-diabetic caucasian patients were recruited at Xiangya Hospital of Central South University; 24 underwent cardiac valve surgery (CTRL group; no evidence of CAD), and 19 underwent coronary artery bypass graft surgery (CAD group). Inclusion criteria of CAD group: CAD patients undergoing coronary artery bypass surgery (confirmed by coronary angiography); inclusion criteria of CTRL group: patients undergoing valve replacement for valvular diseases or surgeries for congenital heart diseases, with age, gender, and BMI matched with CAD group and with no evidence of coronary heart disease or negative findings in coronary angiography. Exclusion criteria include history of PCI, history of thoracotomy, diabetes mellitus, tumor, severe infection, renal insufficiency (glomerular filtration rate < 60 mL/min), active period of autoimmune disease, and severe obesity BMI > 35 kg/m 2 . Specimens were collected in the way as previously described with informed consent signed 36 . EAT samples were collected from the atrioventricular groove near the right coronary artery. EAT samples were immediately snap-frozen and stored using liquid nitrogen for RNAseq and qRT-PCR analysis. The research design was approved by the Medical Ethics Committee of the Xiangya Hospital of Central South University. Hypertension was defined as: systolic pressure ≥ 140 mmHg or diastolic pressure ≥ 90 mmHg without the use of antihypertensive drugs, or blood pressure within the normal range with self-reported use of antihypertensive drugs. Diabetes mellitus was identified if the random plasma glucose concentration was ≥ 11.1 mmol/L, or the fasting plasma glucose concentration was ≥ 7.0 mmol/L, or the postprandial plasma glucose concentration (2 h after meals) was ≥ 11.1 mmol/L in oral glucose tolerance test. Hyperlipidemia was diagnosed if the total cholesterol was ≥ 6.22 mmol/L (240 mg/dL), and/or triglyceride was ≥ 2.26 mmol/L (200 mg/dL), and/or low-density lipoprotein cholesterol was ≥ 4.14 mmol/L (160 mmol/L).
Total RNA extraction and lncRNA library construction. Total RNA was extracted from EAT samples using TRIzol (Invitrogen, Carlsbad, CA, USA) according to the manufacturer's instructions. The RNA quality and quantity were evaluated using NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, MA, USA) and Agilent 2100 bioanalyzer (Thermo Fisher Scientific, MA, USA).
The Ribo-Zero Magnetic Kit (Epicentre) was used to process approximately 1 μg of total RNA in each sample to deplete rRNA. The First Strand Master Mix (Invitrogen) was added to fragment the retrieved RNA. The synthesized cDNA was end-repaired and then 3′adenylated. The ends of these 3′adenylated cDNA fragments were connected using adaptors. PCR Primer Cocktail and PCR Master Mix were used for PCR amplification to enrich cDNA fragments. The PCR product was then purified using Ampure XP beads. Two methods were used for the qualification and quantification of the final library: using Agilent 2100 bioanalyzer to check the size distribution of fragments, and using real-time quantitative PCR (RT-qPCR) (TaqMan Probe) for the library quantification. Qualified libraries were sequenced on the Hiseq 4000 or Hiseq X-ten platform (BGI, Shenzhen, China).

Prediction of lncRNA target genes.
To identify lncRNA cis-regulated target genes, we screened for protein-coding genes located within 10 kb upstream or 20 kb downstream of a differentially expressed lncRNA. The Pearson and Spearman correlation coefficient was calculated. To identify an mRNA as a lncRNA trans-regulated target gene, the correlation coefficient between a lncRNA and mRNA pair r ≥ 0.6.
Gene ontology (GO) and pathway analysis. Quantitative signal data for each gene from CAD and control EAT samples were compared to obtain absolute fold change (FC) values, which were then log2 FC transformed. Differential expression of mRNAs and lncRNAs was defined by a |log2 FC| ≥ 1 and Q value < 0.05. To predict the potential biological roles of the identified lncRNAs, we performed enrichment analyses on differentially expressed mRNAs and lncRNAs targeted genes using GOseq and KOBAS (2.0), respectively. Those enriched GO terms or pathways with a Q value < 0.05 were considered significant.