Identification of ceRNA (lncRNA-miRNA-mRNA) Regulatory Network in Myocardial Fibrosis After Acute Myocardial Infarction

Introduction

Acute myocardial infarction (AMI) is a life-threatening cardiovascular disease with high mortality and morbidity globally, featured by myocardial necrosis caused by acute or persistent ischemia and hypoxia of cardiac muscle.¹ Myocardial fibrosis (MF) is a main pathological changes appear in the heart after AMI, characterized by increasing cardiac fibroblasts activity, increasing the risk of heart failure.² Thus, preventing and reversing cardiac fibrosis may have potential therapeutic effects on improving heart function and even reducing mortality risk.

With the rapid development of RNA-sequencing technology, the roles of non-coding RNAs (ncRNAs) in cardiovascular disease has been extensively studied. Chistiakov et al have highlighted the important role of miRNA in cardiac function under physiological conditions and cardiac dysfunction under pathological conditions.³ Wang et al concluded the regulatory functions or diagnostic potential of miRNAs, siRNAs, lncRNAs and circRNAs in cardiovascular diseases.⁴ Similarly, a review indicated that ncRNAs, including lncRNAs, miRNAs and circRNAs, are emerging key regulators of gene expression of several physiological and pathological processes in AMI, as well as lncRNA/circRNA-miRNA-mediated interaction (ie, competitive endogenous RNA (ceRNA) mechanism).⁵ Thum reviewed the role of ncRNAs in MF and suggested that ncRNAs can serve as powerful therapeutic targets to treat MF.⁶ PLAUR was found to play a vital role in transducing inter-cellular signals from endothelial cells and fibroblast cells to intra-cellular pathways of myocardial cells, leading to gene expression involved in cellular response to hypoxia.⁷ Zhuo et al indicated that lncRNA SNHG8 may serve as diagnostic or prognostic biomarkers of AMI.⁸ It was reported that miR-223 enhanced cell proliferation, migration, and differentiation in cardiac fibrosis, thus mediated cardiac fibrosis after MI partially via the involvement of RASA1.⁹ Based on the above-mentioned findings, therefore, we performed this current work.

In current study, differentially expressed mRNAs (DEmRNAs), miRNAs (DEmiRNAs) and lncRNAs (DElncRNAs) in AMI and MF after AMI were acquired, respectively. Then, shared DEmRNAs/DEmiRNAs/DElncRNAs in AMI and MF after AMI were further determined. Finally, a ceRNA (DElncRNA-DEmiRNA-DEmRNA) regulatory network was constructed to explore the potential molecular mechanisms of MF after AMI.

Materials and Methods

Subjects and Samples

The cohort subjected to RNA-Seq comprised 3 patients with AMI, 3 patients with MF after AMI and 3 healthy individuals. The inclusion criteria of AMI patients were as follows: (1) time of chest pain or distress >30min within 24h, and the level of the myocardial enzyme of creatine kinase (CK)-MB and cardiac troponin T (cTnT) was higher than the normal range; (2) patients had their first episode; (3) patients received no medical or surgical treatment prior to admission; (4) patients had blood samples before admission, at discharge, and 6 months after MI; and (5) patients had complete clinical data. The exclusion criteria of AMI patients were as follows: (1) patients with myocarditis and other diseases caused by chest pain or distress; (2) patients with a history of renal failure, advanced liver disease, malignant tumors, and other inflammatory diseases such as psoriasis and rheumatoid arthritis; (3) recurrent patients; (4) patients with incomplete clinical data; and (5) patients with missing blood samples before admission, at discharge, and 6 months after MI. Those patients diagnosed with MF (1 year after MI) were enrolled. All samples were collected after obtaining written informed consent from every participant. This study was approved by the ethics committee of Shijiazhuang People’s Hospital and performed in accordance with the Declaration of Helsinki. Peripheral whole blood (2.5 mL) drawn from each subject was collected in PAXgene^® RNA blood tubes. Specifically, the blood sample from patients with AMI was collected within 12 h after onset of AMI.

RNA Isolation and Sequencing

According to the manufacturer’s protocol, RNA was extracted with PAXgene blood RNA kit. With Agilent 2100 Bioanalyzer (Agilent RNA 6000 Nano Kit), the concentration, integrity and RNA integrity number (RIN) values of RNA were assessed. After total RNA DNase I treatment, removal of rRNA, RNA interruption, synthesis of reverse transcriptional one and two strands, end repair, cDNA with an “A” at the 3ʹ end, connection of cDNA 5ʹ adapter, digestion of two strands of cDNA, PCR reaction, and product recovery and library quality inspection, DNBSEQ platform (PE100 strategy) was used to perform RNA sequencing for mRNA and lncRNA. The Fastx-Toolkit was used to trim 5ʹ and 3ʹ segments of reads to remove bases with mass <20 and delete reads with N >10%. HISAT2 was applied to align the clean reads with the human reference genome (GRCh38). Expression of mRNAs and lncRNAs was normalized and outputted with StringTie. The Rfam was used for annotation analysis on measured small RNA. The mature miRNA and miRNA precursor sequences were downloaded from miRBase. The expression of miRNA was quantified with miRDeep2. Finally, a total of 19,951 mRNAs, 1313 miRNAs and 16,054 lncRNAs were obtained.

Identification of DEmRNAs, DEmiRNAs and DElncRNAs in AMI and MF After AMI

DEGseq2 (http://bioconductor.org/packages/DESeq2/) was applied to identify DEmRNAs, DEmiRNAs and DElncRNAs in AMI and MF after AMI with |log₂FC| ≥ 1 and p-value < 0.05. Hierarchical clustering analysis of DEmRNAs, DEmiRNAs and DElncRNAs was performed with R (https://www.r-project.org/) package “pheatmap”. Then, shared DEmRNAs in AMI and MF after AMI were further identified with Venny 2.1.0 (https://bioinfogp.cnb.csic.es/tools/venny/). CPDB (http://cpdb.molgen.mpg.de/CPDB) was used to perform GO and KEGG enrichment analysis for up-regulated and down-regulated shared DEmRNAs in AMI and MF after AMI. A value of p < 0.05 was considered to be represented statistically significant.

CeRNA (DElncRNA-DEmiRNA-DEmRNA) Regulatory Network

Six bioinformatic algorithms (PITA, RNA22, miRmap, microT, miRanda and PicTar) were used to predict the targeted DEmRNAs of DEmiRNAs. Then, the DEmiRNA-DEmRNA pairs recorded by ≥4 algorithms were remained for network construction. DIANA-LncBase v2.0 was applied to predict DElncRNA-DEmiRNA interaction pairs. Based on ceRNA theory, the ceRNA (DElncRNA-DEmiRNA-DEmRNA) regulatory network was constructed with Cytoscape, by using a combination of lncRNA-miRNA pairs and miRNA-mRNA pairs.

Quantitative Real-Time Polymerase Chain Reaction (qRT-PCR) Validation

Thirteen blood samples were obtained from 5 patients with AMI (3 for RNA sequence), 3 patients with MF after AMI (2 for RNA sequence) and 5 healthy controls (3 for RNA sequence). Total RNA was isolated with the Trizol reagent (Invitrogen, USA) following manufacturer’s protocol. The qRT-PCR reactions were performed in ABI 7300 Real-time PCR Detection System with SuperReal PreMix Plus (SYBR GREEN) (Invitrogen, USA). Relative gene expression was analyzed by 2^−ΔΔCT method. The human GAPDH and ACTB were used as endogenous controls for mRNA and lncRNA expression in analysis. The human U6 was used as endogenous controls for miRNA in analysis. The PCR primers were displayed in Table S1.

Validation in the GEO Database

To further validate our findings, GSE66360, including miRNA and mRNA expression profiles and consisting of 49 patients experiencing acute myocardial infarction and 50 healthy cohorts, was downloaded from GEO database. Then, the expression levels of selected mRNAs and lncRNAs were validated with GSE66360.

Results

Clinical Characteristics of Enrolled Individuals for RNA Sequencing and qRT-PCR Validation

A total of 14 individuals (5 patients with AMI, 4 patients with MF after AMI, and 5 healthy controls) were enrolled in RNA sequencing and qRT-PCR validation. The clinical information of these individuals was presented in Table S2.

Identification of DEmRNAs, DEmiRNAs and DElncRNAs in AMI and MF After AMI

Compared to normal controls, 1990 DEmRNAs (1047 up- and 943 down-regulated) were identified in AMI, 1477 DEmRNAs (290 up- and 1187 down-regulated) were identified in MF after AMI, respectively (Figure 1A and B). Of these, PFKFB2 and SLC4A10, TVP23C-CDRT4 and UBAP1L was the most up- and down-regulated DEmRNA in AMI and MF after AMI, respectively (Table 1). A total of 509 shared DEmRNAs (84 up- and 425 down-regulated) in AMI and MF after AMI were acquired (Figure S1A).

Table 1 Top 10 Up- and Down-Regulated DEmRNAs

Figure 1 The heatmap of the DEmRNAs in AMI (A) and MF after AMI (B), DEmiRNAs in AMI (C) and MF after AMI (D), and DElncRNAs in AMI (E) and MF after AMI (F). Row and column represented DEmRNAs/DEmiRNAs/DElncRNAs and samples, respectively. The color scale represented the expression levels.

Figure 2 Continue.

Figure 2 Significantly enriched GO terms and KEGG pathways of up-regulated (A–D) and down-regulated (E–H) shared DEmRNAs in AMI and MF after AMI. (A and E) BP, biological process; (B and F) CC, cellular component; (C and G) MF, molecular function; (D and H) KEGG pathways. The x-axis shows counts of DEmRNAs enriched in GO terms or KEGG pathways and the y-axis shows GO terms or KEGG pathways. The color scale represented -lg p-value.

Compared to normal controls, 88 DEmiRNAs (22 up- and 66 down-regulated) were identified in AMI, 36 DEmiRNAs (16 up- and 20 down-regulated) were identified in MF after AMI, respectively (Figure 1C and D). Of these, hsa-miR-223-3p and hsa-miR-6802-3p, hsa-miR-1299 and hsa-miR-197-3p was the most up- and down-regulated DEmiRNA in AMI and MF after AMI, respectively (Table 2). A total of 22 shared DEmiRNAs (6 up- and 16 down-regulated) in AMI and MF after AMI were acquired (Figure S1B).

Table 2 Top 10 Up- and Down-Regulated DEmiRNAs

Compared to normal controls, 473 DElncRNAs (315 up- and 158 down-regulated) were identified in AMI, 502 DElncRNAs (212 up- and 290 down-regulated) were identified in MF after AMI, respectively (Figure 1E and F). Of these, DLEU1 and FP236383.3, SBF2-AS1 and CTBP1-DT was the most up- and down-regulated DElncRNA in AMI and MF after AMI, respectively (Table 3). A total of 90 shared DElncRNAs (39 up- and 51 down-regulated) in AMI and MF after AMI were acquired (Figure S1C).

Table 3 Top 10 Up- and Down-Regulated DElncRNAs

Functional Annotation of Shared DEmRNAs in AMI and MF After AMI

For up-regulated shared DEmRNAs group, the result indicated that these genes were enriched in response to stress (p = 9.94E-06), cell activation (p = 1.25E-04), enzyme regulator activity (p = 2.19E-04), Pantothenate and CoA biosynthesis (p = 3.53E-03) and Leishmaniasis (p = 4.09E-03) (Figure 2A–D). These, down-regulated genes were enriched in lymphocyte costimulation (p = 2.45E-05), leukocyte activation (p = 1.28E-04), protein binding (p = 3.69E-05), Th17 cell differentiation (p = 1.85E-03), NF-kappa B signaling pathway (p = 5.95E-03) (Figure 2E–H).

CeRNA (DElncRNA-DEmiRNA-DEmRNA) Regulatory Network

A total of 509 shared DEmRNAs, 22 shared DEmiRNAs and 90 shared DElncRNAs in AMI and MF after AMI were used to construct the ceRNA regulatory network. By overlapping DElncRNA-DEmiRNA interaction pairs and DEmiRNA-DEmRNA interaction pairs, a ceRNA network was acquired. The ceRNA network contained 9 lncRNA-miRNA pairs and 9 miRNA-mRNA pairs, including 9 lncRNAs, 3 miRNAs, and 8 mRNAs (Figure 3).

Figure 3 CeRNA (DElncRNA-DEmiRNA-DEmRNA) regulatory network. The inverted triangles, rectangles and ellipses indicate DElncRNAs, DEmiRNAs and DEmRNAs, respectively. Red and green color represents up-regulation and down-regulation, respectively.

QRT-PCR Validation and Validation in GEO Database

Three DEmRNAs (BCL9L, CENPB and JAK2), three DEmiRNAs (hsa-miR-135a-5p, hsa-miR-144-3p and hsa-miR-197-3p) and two DElncRNAs (IL21R-AS1 and LEF1-AS1) were selected randomly to use for qRT-PCR validation. Based on our RNA-seq results, JAK2 and hsa-miR-135a-5p were up-regulated while BCL9L, CENPB, hsa-miR-144-3p, hsa-miR-197-3p, IL21R-AS1 and LEF1-AS1 were down-regulated in AMI and MF after AMI. In AMI, all these candidate markers generally exhibited the same pattern as that in our RNA-seq results; while in MF after AMI, except for CENPB, JAK2 and hsa-miR-197-3p, the expression of the others in the qRT-PCR results exhibited the same pattern as that in our RNA-seq results (Figure 4). As shown in Figure 5, the expression of BCL9L, CENPB, LEF1-AS1, LINC00664 and SNHG22 were found to exhibit the same pattern as that in our RNA-seq results.

Figure 4 The qRT-PCR results of the DEmRNAs, DEmiRNAs and DElncRNAs in AMI (A) and MF after AMI (B). The x-axis represents the DEmRNAs/DEmiRNAs/DElncRNAs and the y-axis represents log2 (fold change). *Indicates p < 0.05, **Indicates that p < 0.01. Abbreviation: FC, fold change.

Figure 5 Validation of mRNAs and lncRNAs in GEO. The x-axes represent normal control and AMI groups. The y-axes represent gene expression levels. (A) BCL9L, (B) CENPB, (C) LEF1-AS1, (D) LINC00664, (E) SNHG22.

Discussion

MF is a common hallmark in various heart diseases and predominant cause of morbidity and mortality in the world. The potential biomarkers for cardiovascular disease have been explored widely. It has been reported that ADMA is an independent predictor of cardiovascular morbidity and mortality.^10,11 Zinellu et al suggested that plasma ADMA was involved in carotid narrowing after carotid endarterectomy intervention.¹² Bivona et al indicated that serum Gal 3 levels rise immediately after AMI and decrease significantly within 5 days after the acute event.¹³ In addition, it has been concluded that H-FABP is not a reliable marker for AMI diagnosis.¹⁴ Despite available drug treatment strategy for MF, there is still an urgent need for development of new biomarkers to prevent and reverse MF after AMI.

Many lncRNAs have been reported to be involved in the development of various diseases. In this study, abundant lncRNAs were identified to be associated with AMI and MF after AMI, such as IL21R-AS1 and LEF1-AS1. Cai et al indicated that IL21R-AS1 was a candidate biomarker for coronary artery disease.¹⁵ Overexpressed LEF1-AS1 was detected in coronary artery atherosclerosis, which regulated vascular smooth muscle cell proliferation and migration.¹⁶ In agreement with previous studies, we observed IL21R-AS1 and LEF1-AS1 were significantly dyregulated and validated in the qRT-PCR as well in this analysis, which reminds us to focus on the roles of IL21R-AS1 and LEF1-AS1 in AMI and MF after AMI.

High LINC00664 expression was associated with shorter survival in endometrial and clear cell carcinomas patients, whereas, it is converse for high-grade serous carcinomas patients.¹⁷ Zhang et al determined an eight-lncRNA signature with high prognostic value for colorectal cancer, including LINC00664.¹⁸ Chen et al indicated LINC00664 was a prognostic marker for ovarian carcinoma.¹⁹ Chouvarine et al suggested that hypoxia alone could drive cardiac inflammation and identified miR-197-3p as a hypoxia-regulated miRNA.²⁰ McManus et al examined associations of miRNA levels with 6 cardiometabolic (CM) traits and determined that miR-197-3p was significantly associated with body mass index, triglycerides, and systolic and diastolic blood pressures.²¹ The role of JAK/STAT signaling pathway in the heart has been extensively studied. Negoro et al demonstrated that JAK/STAT signaling pathway was activated in the rat model of AMI.²² Mascareno et al suggested that activation of the JAK/STAT signaling pathway contributed to the pathogenesis of myocardial ischemia.²³ Chen et al indicated that activation of JAK2/STAT3 pathway can prevent apoptosis after AMI.²⁴ Zhang et al revealed inhibition of JAK2/STAT3 pathway can attenuate myocardial hypertrophy induced by transverse aortic constriction.²⁵ Chen et al implied that JAK/STAT pathway made contribution to left-atrial-fibrogenesis.²⁶ Eid et al reported that activation of JAK2/STAT3 signaling contributed to reduce MI-induced left ventricle injury.²⁷ In this study, JAK2 was a target of hsa-miR-197-3p, and LINC00664 may function as a ceRNA to capture hsa-miR-197-3p.

Altered expression of GAS6-AS1 has been correlated with poor outcome in various tumors. Han et al indicated that decreased GAS6-AS1 was related to poor prognosis in non-small cell lung cancer patients.²⁸ Zhang et al demonstrated that high GAS6-AS1 expression was associated with advanced TNM stage in gastric cancer.²⁹ Ai et al suggested that up-regulated GAS6-AS1 was negatively correlated with the overall survival of hepatocellular carcinoma patients.³⁰ Similarly, SNHG22, located on chromosome 18q21.1, has been associated with multiple cancers as well. Zhang et al indicted that high expressed SNHG22 predicted poor prognosis in patients with epithelial ovarian carcinoma.³¹ Fang et al found that SNHG22 exerted carcinogenic effects in triple-negative breast cancer.³² Gao et al suggested that SNHG22 overexpression was correlated with worse overall survival in papillary thyroid cancer patients.³³ However, to our best knowledge, the role of GAS6-AS1 and SNHG22 in the heart has not been reported yet.

Recent studies indicated that miR-135a-5p may play important roles in a variety of tumors, including gliomas, lung cancer and head and neck squamous cell carcinoma.^34–36 In addition, Wu et al suggested that miR-135a inhibited cardiac fibrosis and it may be a therapeutic target for the amelioration of cardiac fibrosis.³⁷ Xu et al demonstrated that miR-135a-5p negatively regulated hypoxia-induced cardiomyocyte injury in H9c2 cells.³⁸ CENPB, a cell cycle-related gene, was reported to be differentially regulated in the heart tissue in serum-response factor null embryonic stem cells.³⁹ BCL9L is a nuclear co-activator in the Wnt/β-catenin pathway, which was shown to be associated with various cancers. BCL9L knockdown was reported to inhibit pancreatic cancer progression and liver metastasis.⁴⁰ Huge et al indicated that high expression of BCL9L was associated with poor overall survival in hepatocellular carcinoma patients.⁴¹ In addition, Cantù et al suggested a causative link between mutations in BCL9L and human congenital heart disease.⁴² Both CENPB and BCL9L were targeted by hsa-miR-135a-5p in this analysis. We speculated that GAS6-AS1 and SNHG22 may participate in development of MF after AMI via acting as ceRNAs of CENPB and BCL9L by sponging hsa-miR-135a-5p.

Conclusions

In conclusion, we speculated that LINC00664/hsa-miR-197-3p/JAK2 and GAS6-AS1/SNHG22/hsa-miR-135a-5p/CENPB/BCL9L interaction pairs may serve as potential biomarkers in MF after AMI. Due to small sample size included in this study, further research with a larger cohort is essential to be performed to confirm the conclusion. In addition, the specific mechanism of the lncRNA-miRNA-mRNA axes in MF after AMI need to be validated in vivo and in vitro.