Characteristics of circular RNA expression in lung tissues from mice with hypoxia‑induced pulmonary hypertension.

Pulmonary hypertension (PH) is a life‑threatening lung disease, characterized by an increase in pulmonary arterial pressure caused by vasoconstriction and vascular remodeling. The pathogenesis of PH is not fully understood, and there is a lack of potential biomarkers for the diagnosis and treatment of patients with PH. Non‑coding RNAs with a characteristic covalently closed loop structure, termed circular RNAs (circRNAs), are present in a number of pulmonary diseases. To the best of our knowledge, the present study is the first to use microarray analysis to determine the expression profile of circRNAs in lung tissues from mice with hypoxia‑induced PH. In total, 23 significantly upregulated and 41 significantly downregulated circRNAs were identified. Of these, 12 differentially expressed circRNAs were selected for further validation using reverse transcription‑quantitative polymerase chain reaction. Putative microRNAs (miRNAs) that bind to the dysregulated circRNAs were predicted. Subsequently, bioinformatics tools were used to construct circRNA‑miRNA‑mRNA networks for the two most promising circRNAs, namely mmu_circRNA_004592 and mmu_circRNA_018351. Gene Ontology and Kyoto Encyclopedia of Genes and Genomes analyses of target genes of the dysregulated circRNAs revealed that these dysregulated circRNAs may serve an important role in the pathogenesis of hypoxia‑induced PH. Therefore, these dysregulated circRNAs are candidate diagnostic biomarkers and potential therapeutic targets for PH.

Introduction

Pulmonary hypertension (PH) is a complex and life-threatening type of lung disease, which is defined by a mean pulmonary arterial pressure of ≥25 mmHg at rest as measured by right heart catheterization (1). Recently, PH was clinically classified into five groups: Pulmonary arterial hypertension, PH due to left heart disease, PH due to lung diseases and/or hypoxia, chronic thromboembolic PH, and PH with unclear and/or multifactorial mechanisms (2). Hypoxia-induced PH, which is the focus of the present study, results from hypoxic pulmonary vasoconstriction, polycythemia and pulmonary arterial remodeling (3). Pathological changes in patients with hypoxia-induced PH manifest as dysfunction of smooth muscle cells, fibroblasts and pulmonary artery endothelial cells, which are associated with progressive obliteration of the pulmonary arteries (4,5). Several pathways participate in the complex process of hypoxia-induced PH; however, the underlying mechanism is not well-defined. Furthermore, the current treatment strategies available are limited and the prognosis is poor for patients with PH (6). Thus, novel strategies for the diagnosis and treatment of patients with hypoxia-induced PH are urgently required.

Circular RNAs (circRNAs) are a novel class of the non-coding RNA family that feature a covalently closed loop structure without a 5′-terminal cap or 3′-terminal poly(A) tail (7). Initially, circRNAs were considered to be functionless byproducts of the normal splicing processing. With the advent of next-generation sequencing, it was revealed that circRNAs were widely distributed stable entities in various cell types and were expressed in species-, tissue-, disease- and developmental stage-specific manners (8,9). To date, the functions of circRNAs have not been fully elucidated. The most widely studied function of circRNAs is their ability to function as miRNA sponges. The representative circRNA, ciRS-7, has >70 conventional binding sites, which bind to miR-7 to decrease its expression levels and further increase the levels of miR-7-targeted transcripts (10). CircRNA located in the nucleus may participate in the regulation of transcription and alternative splicing (11). Furthermore, a number of circRNAs have been demonstrated to be translated into proteins. Legnini et al (12) firstly demonstrated that circ-ZNF609 could be translated into a small protein that specifically controlled the myoblast proliferation. These unique features of circRNAs have brought a new perspective to our understanding of the biological mechanisms underlying various diseases, including cardiovascular diseases, cancer, neurological disorders and diabetes, and the possibility of novel diagnostic biomarkers and innovative strategies for treatment (13–16). However, the role of circRNAs in hypoxia-induced PH remains unknown.

In the present study, a murine model of hypoxia-induced PH was established. To the best of our knowledge, circRNA microarray analysis was used for the first time to determine the differential expression profile of circRNAs in the lungs of mice with hypoxia-induced PH compared with the control mice. A total of 12 candidate circRNAs were selected for further validation by reverse transcription-quantitative polymerase chain reaction (RT-qPCR). Bioinformatics analyses identified the potential functions of the candidate circRNAs in the mechanism of hypoxia-induced PH. These candidate circRNAs have the potential to function as diagnostic biomarkers and targets for the treatment of patients with the disease.

Materials and methods

Animal models

A total of 20 male C57BL/6 mice (6–8 weeks old, weighing 17–20 g) were obtained from JSJ Laboratory Animal Co., Ltd. (Shanghai, China). The experimental protocols involving the animals were approved by the Institutional Animal Care and Use Committee of Huadong Hospital, Fudan University (Shanghai, China). All mice were housed in room temperature and relative humidity of ~60% environment with a 12-h light-dark cycle, and fed with normal chow and water. Following acclimatization for 1 week, mice were randomly assigned into two groups, including the hypoxia-PH group (n=10) and the control group (n=10). The hypoxia-PH mice were exposed to hypoxia (10% O2) in a normobaric chamber for 21 days (17). The control mice were housed in room air. The feeding conditions were the same for both groups.

PH measurements

The right ventricular systolic pressure (RVSP) and right ventricular hypertrophy (RVH) were considered to be indicators of PH. At the end of the treatment period, mice were anesthetized with 20% urethane (0.4 ml/100 g, intraperitoneal), and then a catheter was inserted into the right ventricle (RV) via the right jugular vein. RVSP was measured by a gauge pressure transducer (RM-6000; Nihon Kohden Corporation, Tokyo, Japan). Animals were then sacrificed and the hearts were collected. The RV was dissected from the left ventricle (LV) and septum, and the two ventricles were then weighed. The RVH was determined based on the ratio of the RV weight compared with the weight of the LV and septum [calculated as RV/(LV+S)]. Additionally, lung tissues were harvested and stored at −80°C for subsequent circRNA microarray analysis.

CircRNA microarray analysis

Lung samples from the hypoxia-PH group were randomly distributed into three subgroups for further analysis. The two subgroups had three lung samples from the hypoxia-PH group and the other subgroup has four lung samples from the hypoxia-PH group. The same grouping strategy was applied to the control group. Total RNA from each subgroup was isolated using TRIzol reagent (Invitrogen; Thermo Fisher Scientific, Inc., Waltham, MA, USA) according to the manufacturer’s protocol. The RNA concentration was quantified using a NanoDrop spectrophotometer (ND-1000; Thermo Fisher Scientific, Inc., Wilmington, DE, USA), while RNA integrity was determined using gel electrophoresis. Following RNA quality control, fluorescence labeling and microarray hybridization were performed according to the manufacturer’s protocol for the Arraystar Super RNA Labeling kit (Arraystar, Inc., Rockville, MD, USA). Briefly, linear RNAs were eliminated using RNase R (Epicentre; Illumina, Inc., San Diego, CA, USA), and the enriched circRNAs were subsequently amplified and transcribed into fluorescent complementary RNA (cRNA) by employing a random priming method (Arraystar Super RNA Labeling kit; Arraystar, Inc.). Labeled cRNA was then purified using the RNeasy Mini kit (Qiagen GmbH, Hilden, Germany). Labeled cRNA quality and quantity were determined using the NanoDrop ND-1000 spectrophotometer. Subsequently, the labeled cRNA was fragmented and hybridized onto an Arraystar mouse circRNA microarray (6x7K; Arraystar, Inc.). Subsequent to washing with gene expression wash buffer 1 (Agilent Technologies, Inc., Santa Clara, CA, USA) for 1 min at 37°C and then gene expression wash buffer 2 (Agilent Technologies, Inc.) for 1 min at 37°C, the hybridized arrays were scanned using an Axon GenePix 4000B microarray scanner (Molecular Devices, LLC, Sunnyvale, CA, USA).

CircRNA data analysis

GenePix Pro 6.0 software (Axon; Molecular Devices, LLC) was used to analyze the scanned images of circRNA microarray in order to extract the raw data. Normalization and data processing were accomplished using the R software package (version 3.4.1; R Foundation Inc., Vienna, Austria). The statistical significance of the expression of circRNAs between the two groups was estimated by fold-change filtering and an unpaired t-test. To increase the number of circRNAs detected, the screening criterion for differentially expressed circRNAs was a fold-change of ≥1.5 and a P-value of <0.05. Hierarchical cluster analyses, box plot, volcano plot and scatter plot were used to display the differential expression of circRNAs between the two groups.

RT-qPCR validation

The top 12 candidate dysregulated circRNAs (including six upregulated circRNAs and six downregulated circRNAs) were selected for validation by RT-qPCR. Briefly, total RNA was extracted from the two groups using TRIzol reagent (Thermo Fisher Scientific, Inc.). The concentration of total RNA was measured by a NanoVue Plus spectrophotometer (GE Healthcare, Inc., Chicago, IN, USA) and total RNA (1 µg) was reverse transcribed into cDNA using a PrimeScript™ RT Reagent kit (Takara Bio, Inc., Otsu, Japan) according to the manufacturer’s protocol. The RT reaction was conducted at 37°C for 15 min and 85°C for 5 sec. Next, the RT-qPCR reaction was performed using SYBR® Premix Ex Taq™ (Takara Bio, Inc.) using a Bio-Rad Real-Time PCR system (Bio-Rad Laboratories, Inc., Hercules, CA, USA). The amplification reaction involved initial denaturation at 95°C for 1 min, and 40 cycles of denaturation at 95°C for 5 sec and 60°C for 30 sec. The primers for the circRNAs were designed in divergent orientation, and the primer sequence pairs are presented in Table I. The relative level of circRNA expression for each target gene is shown as the 2−ΔΔCQ value (18).

Table I

Primers of validated circRNAs.

Name	Forward	Reverse
GAPDH	5′-GTTGTCTCCTGCGACTTCA-3′	5′-GCCCCTCCTGTTATTATGG-3′
mmu_circRNA_010744	5′-GGCTTACACCTCAACTACCATCC-3′	5′-CCCAGCATTCCGAAGAACC-3′
mmu_circRNA_014597	5′-CTTCAATGATTTTCACCTCCAG-3′	5′-AGCGGTCATCCATGCTTATAT-3′
mmu_circRNA_008379	5′-CTGTCCCAACTGTAAAGAAGGTG-3′	5′-CATCGGTTTGGTGCTCCTC-3′
mmu_circRNA_000267	5′-GGCACCCAACCTCATCTTC-3′	5′-TTCGCACTTCTCATAGTCACTCA-3′
mmu_circRNA_015666	5′-ACGCCCAGGATTATGAAGAGG-3′	5′-GGTACTGTGGTGAGTCTGGACGA-3′
mmu_circRNA_013176	5′-CCCCACTGAACAGGAAGAAGA-3′	5′-ACGATACAAGCACCCAGATAGG-3′
mmu_circRNA_016128	5′-CTTCTGGTTCCTGGAGTTGG-3′	5′-AGTTCTGGGATGCTTGCTTTA-3′
mmu_circRNA_018351	5′-GCTCAACCTCATCCAGCACG-3′	5′-GCGAGTCACAGCCTTCCAT-3′
mmu_circRNA_004592	5′-TGTGGTTGGCAAGGATGTC-3′	5′-CGCATGGGAACTGTAGTAAGAT-3′
mmu_circRNA_018701	5′-GAAGATAGTGCCTTTCAGCCATAC-3′	5′-CACTCCCGTTTATCCCACTCA-3′
mmu_circRNA_016636	5′-GGAAAGGAGCAAACACGAAGG-3′	5′-TGTCAGCGGCAAGTCATCG-3′
mmu_circRNA_010120	5′-GCAGTACGCCGTCATCAGTG-3′	5′-GGAAGGCAGCAGACAAGAGC-3′

circRNA, circular RNA.

Bioinformatics analysis

CircRNAs contain miRNA response elements (MREs) and function as miRNA sponges, which regulate miRNA and downstream gene expression through a competing endogenous RNA (ceRNA) mechanism (19). Putative circRNA-miRNA-mRNA networks were generated using various bioinformatics tools. The circRNA-miRNA interaction was predicted using the miRNA target prediction software (Arraystar, Inc.), which is based on the TargetScan (20) and miRanda (21). The top five circRNA-bound miRNAs were identified based on the prediction of miRNA binding sites. To ensure the efficiency and accuracy of the analysis, circRNAs were selected (fold-change ≥2 and P<0.05) and their target genes were predicted using miRWalk 2.0, which integrates miRWalk, RNA22, miRanda and TargetScan (22). The top 10 circRNA target genes for each miRNA were identified and used to construct circRNA-miRNA-mRNA networks using Cytoscape software (version 3.5.1; The Cytoscape Consortium, San Diego, CA, USA). The functional roles of the circRNA-target genes were further analyzed using Gene Ontology (GO; geneontology.org) and were classified into three categories, including biological process (BP), molecular function (MF) and cellular component (CC). In addition, the biological pathways of circRNA-target genes were identified using Kyoto Encyclopedia of Genes and Genomes (KEGG; www.kegg.jp). The GO and KEGG analyses were performed using DIANA-miRPath, version 3.0 (23).

Statistical analysis

Validated data from RT-qPCR are expressed as the mean ± standard error of the mean. Differences between the groups were assessed by SPSS statistical software version 19.0 (IBM Corp., Armonk, NY, USA) using one-way analysis of variance. P<0.05 was considered to indicate a statistically significant difference.

Results

Expression profile of circRNAs

The murine model of hypoxia-induced PH was evaluated by measuring the RVSP and RVH (24). When compared with the control group, the hypoxia-induced PH group had significantly increased levels of RVSP (30.77 vs. 22.96 mmHg, P<0.01) and RVH (0.39 vs. 0.28, P<0.01; Fig. 1A and B). Subsequent to obtaining the aforementioned measurements, the lungs were harvested for circRNA microarray detection. The box plot in Fig. 2A demonstrates that there was no abnormal distribution of circRNA expression in the six samples after normalization. Next, a scatter plot was used to assess the expression variation of circRNAs between the hypoxia-PH and control groups (Fig. 2B). The hierarchical cluster analysis revealed the differential expression profile of circRNAs in these two groups (Fig. 2C). A volcano plot was then constructed using fold-change values and P-values in order to visualize the differential expression of circRNAs between the two groups (Fig. 2D). In total, 64 differentially expressed circRNAs were identified between the hypoxia-PH and control groups, of which 23 circRNAs were upregulated and 41 circRNAs were downregulated. The top 20 upregulated and top 20 downregulated circRNAs are listed in Table II.

Figure 1

Measurement of RVSP and RVH in mice with hypoxia-induced pulmonary hypertension. (A) RVSP was significantly increased in mice with hypoxia-induced PH, compared with that in the control group. (B) RVH was significantly increased in mice with hypoxia-induced PH, compared with that in the control group. Values are presented as the mean ± standard error of the mean (n=10). **P<0.01. RVSP, right ventricular systolic pressure; RVH, right ventricular hypertrophy; PH, pulmonary hypertension; RV/(LV+S), right ventricle/(left ventricle + septum) weight.

Figure 2

Microarray analysis of differential expression of circRNAs in mice with hypoxia-induced PH. (A) Box plot showing that the normalized intensity of the hypoxia-PH and control groups was nearly the same. (B) Scatter plot demonstrating the expression variation of the circRNAs in the two groups. The values of x and y-axes in the scatterplot represent the normalized signal values of the samples (log2 scaled). The red points outside the top green line represent the upregulated circRNAs with ≥1.5-fold change and the blue points outside bottom green line represent the down-regulated circRNAs with ≥1.5-fold change. (C) Hierarchical cluster analysis revealed the expression profile of the dysregulated circRNAs in the two groups. (D) Volcano plot showing the differentially expressed circRNAs in the two groups, with red points representing the significantly different circRNAs. PH, pulmonary hypertension; circRNA, circular RNA.

Table II

Top 20 upregulated and top 20 downregulated circRNAs in mice with hypoxia-induced pulmonary hypertension.

CircRNA	Fold change	P-value	Regulation	CircRNA type	Chromosome	Gene symbol
mmu_circRNA_004592	2.925894	0.001652	Up	Exonic	X	Fam120c
mmu_circRNA_018701	2.717783	0.003849	Up	Exonic	10	Jmjd1c
mmu_circRNA_010744	2.469385	0.03412	Up	Exonic	1	Sntg1
mmu_circRNA_015666	2.211677	0.004139	Up	Exonic	4	Ikbkap
mmu_circRNA_016128	2.168588	0.011483	Up	Exonic	7	Nars2
mmu_circRNA_016636	2.085027	0.006036	Up	Exonic	12	Ttc7b
mmu_circRNA_002378	2.02187	0.000377	Up	Intronic	2	Kif5c
mmu_circRNA_011174	1.834095	0.035723	Up	Intragenic	17	XLOC_010800
mmu_circRNA_004941	1.82896	0.001544	Up	Intragenic	13	Cdc14b
mmu_circRNA_015350	1.759724	0.000263	Up	Exonic	1	Ncoa2
mmu_circRNA_011414	1.756762	0.029613	Up	Exonic	18	Tcf4
mmu_circRNA_008988	1.717354	0.000246	Up	Exonic	4	Unc13b
mmu_circRNA_011332	1.691847	0.003271	Up	Intragenic	8	Large
mmu_circRNA_014444	1.684291	0.025843	Up	Exonic	4	Haus6
mmu_circRNA_006839	1.632387	0.044651	Up	Intronic	17	Crim1
mmu_circRNA_007240	1.632184	0.021443	Up	Exonic	5	Bre
mmu_circRNA_017590	1.627719	0.025193	Up	Exonic	6	Ccdc132
mmu_circRNA_013053	1.59081	0.001276	Up	Exonic	6	Gmcl1
mmu_circRNA_000113	1.566056	0.02078	Up	Antisense	12	Rian
mmu_circRNA_012851	1.561305	0.004609	Up	Exonic	2	Arhgap21
mmu_circRNA_010120	4.557204	0.028091	Down	Exonic	14	Cryl1
mmu_circRNA_018351	3.479684	0.00028	Down	Exonic	8	Zfhx3
mmu_circRNA_013176	2.940519	0.036055	Down	Exonic	4	Frmd3
mmu_circRNA_012978	2.189961	0.025108	Down	Intronic	17	Tsc2
mmu_circRNA_007835	2.181447	0.049532	Down	Intragenic	13	Cdc14b
mmu_circRNA_008379	2.167655	0.010006	Down	Exonic	2	Sp3
mmu_circRNA_001170	2.150695	0.032906	Down	Intragenic	10	Hsp90b1
mmu_circRNA_000267	2.068873	0.007548	Down	Exonic	3	Plekho1
mmu_circRNA_014597	2.014597	0.021489	Down	Exonic	1	Kifap3
mmu_circRNA_003780	1.983159	0.011809	Down	Exonic	13	Cdyl
mmu_circRNA_005086	1.901804	0.004295	Down	Intragenic	10	Tbpl1
mmu_circRNA_010464	1.894022	0.025737	Down	Exonic	6	Ezh2
mmu_circRNA_009705	1.859579	0.020918	Down	Exonic	7	Mettl9
mmu_circRNA_001362	1.809006	0.045429	Down	Exonic	5	Gm15800
mmu_circRNA_018468	1.776092	0.017304	Down	Exonic	18	Arhgap26
mmu_circRNA_006562	1.752095	0.026256	Down	Intragenic	10	Erbb3
mmu_circRNA_002959	1.727222	0.020858	Down	Exonic	8	Phkb
mmu_circRNA_005155	1.710652	0.000376	Down	Exonic	16	Robo1
mmu_circRNA_000270	1.694113	0.026716	Down	Exonic	2	Fam188a
mmu_circRNA_009674	1.687043	0.046010	Down	Exonic	5	Rnf4

circRNA, circular RNA.

Validation of circRNAs expression

In order to validate the results of the circRNA microarray, 12 differentially expressed circRNAs (including six upregulated and six downregulated circRNAs) were selected for RT-qPCR analysis on the basis of the fold-change and P-values. The six upregulated circRNAs were as follows: mmu_circRNA_004592, m mu_ci rcR NA_018701, m mu_ci rcR NA_01074 4, m mu_ci rcR NA_015666, m mu_ci rcR NA_016128 and mmu_circRNA_016636. The six downregulated circRNAs were as follows: mmu_circRNA_010120, m mu_ci rcR NA_018351, m mu_ci rcR NA_013176, mmu_circRNA_008379, mmu_circRNA_000267 and mmu_circRNA_014597. Analysis by RT-qPCR revealed that five of the upregulated circRNAs and two of the down-regulated circRNAs, exhibited changes in their expression levels that were consistent with the results of the circRNA microarray. More specifically, mmu_circRNA_004592, mmu_circRNA_018701, mmu_circRNA_010744, mmu_circRNA_015666 and mmu_circRNA_016636 were found to be significantly upregulated in the hypoxia-PH group compared with their levels in the control group (Fig. 3A–E), while mmu_circRNA_018351 and mmu_circRNA_008379 were significantly downregulated compared with their levels in the control group (Fig. 3F and G).

Figure 3

Reverse transcription-quantitative polymerase chain reaction analysis for validation of differentially expressed circRNAs. A total of 12 differentially expressed circRNAs (6 upregulated and 6 downregulated circRNAs) were selected for validation. (A) mmu_circRNA_004592, (B) mmu_circRNA_018701, (C) mmu_circRNA_010744, (D) mmu_circRNA_015666 and (E) mmu_circRNA_016636 were significantly upregulated, while (F) mmu_circRNA_018351 and (G) mmu_circRNA_008379 were significantly downregulated. *P<0.05 and **P<0.01. PH, pulmonary hypertension; circRNA, circular RNA.

Identification of circRNA-targeting miRNAs

CircRNAs typically have multiple MREs to bind miRNAs, which reduces the expression of miRNAs and indirectly upregulates the expression of miRNA target genes (7). Thus, putative circRNA-miRNA interactions were predicted according to the complementary matching sequence using the TargetScan and miRanda tools. All circRNAs identified in the present study were found to have miRNA binding sites. The top five miRNA binding sites for dysregulated circRNAs were predicted (data not shown). Table III lists the seven validated circRNAs and their top five miRNA binding sites. The top five potential miRNA binding sites for the most significantly upregulated circRNA (namely mmu_circRNA_004592) were mmu-miR-742-3p, mmu-miR-6373, mmu-miR-880-5p, mmu-miR-298-3p and mmu-miR-152-3p. In addition, the top five potential miRNA binding sites for the most significantly downregulated circRNA (namely mmu_circRNA_018351) were mmu-miR-6992-3p, mmu-miR-5133, mmu-miR-6936-5p, mmu-miR-6938-5p and mmu-miR-7015-3p.

Table III

Top five miRNA binding sites for validated circRNAs.

CircRNA	MRE1	MRE2	MRE3	MRE4	MRE5
mmu_circRNA_016636	mmu-miR-7239-3p	mmu-miR-6408	mmu-miR-877-5p	mmu-miR-149-3p	mmu-miR-7047-5p
mmu_circRNA_015666	mmu-miR-1903	mmu-miR-7008-3p	mmu-miR-103-2-5p	mmu-miR-103-1-5p	mmu-miR-679-5p
mmu_circRNA_010744	mmu-miR-8113	mmu-miR-7012-5p	mmu-miR-7674-5p	mmu-miR-432	mmu-miR-6350
mmu_circRNA_018701	mmu-miR-653-5p	mmu-miR-17-3p	mmu-miR-6394	mmu-miR-6948-3p	mmu-miR-148b-3p
mmu_circRNA_004592	mmu-miR-742-3p	mmu-miR-6373	mmu-miR-880-5p	mmu-miR-298-3p	mmu-miR-152-3p
mmu_circRNA_008379	mmu-miR-6997-5p	mmu-miR-378b	mmu-miR-7233-3p	mmu-miR-6359	mmu-miR-1967
mmu_circRNA_018351	mmu-miR-877-3p	mmu-miR-1903	mmu-miR-667-5p	mmu-miR-207	mmu-miR-665-5p

circRNA, circular RNA; MRE, miRNA binding site; miR, microRNA.

Construction of circRNA-miRNA-mRNA networks

The function of circRNA as an miRNA sponge is the main biological role of circRNA in various diseases. Thus, circRNA-miRNA-mRNA network analysis was used to predict the potential function of differentially expressed circRNAs in hypoxia-induced PH. The two most significant circRNAs (mmu_circRNA_004592 and mmu_circRNA_018351) were selected for ceRNA analysis. The circRNA-targeting genes were predicted by integrating four miRNA-target prediction tools, including miRWalk, RNA22, miRanda and TargetScan.

For mmu_circRNA_004592, there was no predicted data for one of the five miRNAs targets (mmu-miR-6373), and therefore this target does not appear in the analysis shown in Fig. 4A. For mmu_circRNA_004592, the circRNA-miRNA-mRNA network was constructed (Fig. 4A). The upregulated mmu_circRNA_004592 was predicted to inhibit the expression level of mmu-miR-742-3p, mmu-miR-880-5p, mmu-miR-298-3p and mmu-miR-152-3p, and further promote the expression of miRNA-targeted genes. The top ten miRNA-targeted genes are presented in Fig. 4A. The sequence analysis of MREs of mmu_circRNA_004592 was presented in Fig. 4B. The 2D structure demonstrated the MRE sequence, the target miRNA seed type and the 3′ pairing sequence for mmu_circRNA_004592. Similarly, the circRNA-miRNA-mRNA network for mmu_circRNA_018351 was constructed (Fig. 5A). The downregulated mmu_circRNA_018351 was predicted to increase the expression level of mmu-miR-6992-3p, mmu-miR-5133, mmu-miR-6936-5p, mmu-miR-6938-5p and mmu-miR-7015-3p, and further reduce the expression of miRNA-targeted genes. The top ten miRNA-targeted genes are presented in Fig. 5A. The sequence analysis of MREs of mmu_circRNA_018351 identified the MRE sequence, the target miRNA seed type and the 3′ pairing sequence for mmu_circRNA_004592 in 2D sequence (Fig. 5B).

Figure 4

CircRNA-miRNA-mRNA network of mmu_circRNA_004592. (A) Putative circRNA-miRNA-mRNA network of mmu_circRNA_004592 (pink circle) was constructed using bioinformatics tools. (B) Predicted miRNA response elements for mmu_circRNA_004592. circRNA, circular RNA; miRNA, microRNAs.

Figure 5

CircRNA-miRNA-mRNA network of mmu_circRNA_018351. (A) Putative circRNA-miRNA-mRNA network of mmu_circRNA_018351 (blue circle) was constructed using bioinformatics tools. (B) Predicted miRNA response elements for mmu_circRNA_018351. circRNA, circular RNA; miRNA, microRNAs.

GO and KEGG analyses for circRNA-targeting genes

CircRNAs with a fold-change of ≥2 and a P-value of <0.05 were selected from the dysregulated circRNAs in the hypoxia-induced PH group, and the role of these circRNAs was further investigated by GO and KEGG analyses. For the target genes of the upregulated circRNAs, the top 10 enriched GO terms in the BP and CC categories, as well as the total eight GO terms in the MF category were identified. The top three enriched GO terms were biological processes, cellular nitrogen compound metabolic processes and anatomical structure development in BP; cellular component, cell and intracellular in CC; molecular function, ion binding and enzyme binding in MF (Fig. 6A). For the target genes of the downregulated circRNAs, the top 10 enriched GO terms in the BP and CC categories, as well as the total five GO terms in the MF category were displayed. The top three enriched were biological processes, anatomical structure development and cellular nitrogen compound metabolic process in BP; cellular component, cell and intracellular in CC; molecular function, ion binding and nucleic acid binding transcription factor activity in MF (Fig. 6B).

Figure 6

GO analysis of the target genes of circRNAs. The target genes of dysregulated circRNAs with 2-fold change and P<0.05 were identified using bioinformatics tools. Functions of the target genes of the (A) upregulated and (B) circRNAs were enriched by GO analysis. The total number of genes is shown next to each bar. GO, Gene Ontology; circRNA, circular RNA.

The top 20 KEGG pathways were also identified for the target genes of the upregulated and downregulated circRNAs, respectively (Fig. 7A–D). The top five enriched pathways for the target genes of the upregulated circRNAs included pathways in cancer, the PI3K-Akt signaling pathway, endocytosis, HTLV-I infection and regulation of actin cytoskeleton (Fig. 7A and B). Furthermore, the top five enriched pathways for the target genes of the downregulated circRNAs involved pathways in cancer, the PI3K-Akt signaling pathway, regulation of actin cytoskeleton, proteoglycans in cancer, and the MAPK signaling pathway (Fig. 7C and D).

Figure 7

KEGG pathway analysis of the target genes of circRNAs. circRNAs with 2-fold change and P<0.05 were selected from the dysregulated circRNAs, and the target genes of these circRNAs were identified using bioinformatics tools. (A) The top 20 pathways of the target genes of upregulated circRNAs were identified using KEGG analysis according to the number of enriched genes. (B) The top 20 pathways of the target genes of upregulated circRNAs were identified using KEGG analysis according to the rich factor, Q value and the number of enriched genes. (C) The top 20 pathways of the target genes of downregulated circRNAs were identified using KEGG analysis according to the number of enriched genes. (D) The top 20 pathways of the target genes of downregulated circRNAs were identified using KEGG analysis according to rich factor, Q value and the number of enriched genes. KEGG, Kyoto Encyclopedia of Genes and Genomes; circRNA, circular RNA.

Discussion

To the best of our knowledge, the present study is the first to perform circRNA microarray analysis to detect the differential expression profile of circRNAs in the lungs of mice with hypoxia-induced PH. In total, 23 significantly upregulated and 41 significantly downregulated circRNAs were identified. Of these, 12 differentially expressed circRNAs were selected for further validation using RT-qPCR. Bioinformatics tools were used to construct the circRNA-miRNA-mRNA networks, and the putative functions and pathways were obtained using GO and KEGG enrichment analyses. The results from circRNA-miRNA-mRNA networks analysis and GO and KEGG analysis suggested that dysregulated circRNAs may serve key roles in the pathogenesis of hypoxia-induced PH.

Pulmonary vascular remodeling is a key feature of hypoxia-induced PH with associated dysfunction of endothelial cells, smooth muscle cells and fibroblasts (4). Recently, novel molecular mechanisms of PH have been identified with particular advances in ‘omics’ approaches (25). Several types of noncoding RNAs, such as long noncoding RNAs and miRNAs, participate in the pathological process of hypoxia-induced PH, and certain of these are considered to be diagnostic or therapeutic targets (26,27). CircRNAs are novel members of the non-coding RNA family and serve an important role in the development of several pulmonary diseases. The biological functions of circRNAs are varied; however, their function as miRNA sponges is likely to be the primary mechanism of action in a number of pulmonary diseases, including non-small cell lung cancer (28,29), pulmonary tuberculosis (30), acute respiratory distress syndrome (31) and chronic thromboembolic PH (32). However, to the best of our knowledge, no previous study to date has investigated the expression of circRNAs in hypoxia-induced PH. Thus, to further elucidate the molecular RNA-based mechanism of PH, the circRNA expression profile in hypoxia-induced PH was determined in the present study, using microarray analysis and bioinformatics to identify several circRNAs that are potentially important in the development of hypoxia-induced PH.

Based on the circRNA target prediction analysis, all of the differentially expressed circRNAs that were identified contained MREs to sponge different miRNAs. Of these, the two best candidate circRNAs were selected, namely mmu_circRNA_004592 and mmu_circRNA_018351, and their circRNA-miRNA-mRNA networks were constructed using bioinformatics tools. The upregulated mmu_circRNA_004592 was predicted to negatively regulate the function of mmu-miR-742-3p, mmu-miR-6373, mmu-miR-880-5p, mmu-miR-298-3p and mmu-miR-152-3p. Although there is no direct evidence to show the role of these miRNAs in hypoxia-induced PH, Wu et al (33) have previously reported that hypoxia significantly inhibited the expression of miR-152 in human umbilical vein endothelial cells and that miR-152 targeted ADAM17 to inhibit cell migration and proliferation. Furthermore, several studies have suggested that miR-152 is a tumor suppressor (34). miR-152 has been reported to inhibit cancer cell proliferation, invasion and migration by targeting Wnt1 (35), PIK3CA (36), DNMT1 (37) and PTEN (38). Additionally, Gu et al (39) identified that miR-152 promoted senescence in stem cells from human dental pulp by targeting SIRT7 expression. Since mmu_circRNA_004592 expression was upregulated in the present study, it can be hypothesized that miR-152 is downregulated, which promotes proliferation, and inhibits senescence and apoptosis of endothelial cells, smooth muscle cells and fibroblasts in the pulmonary artery. These actions may further raise the pressure in the pulmonary vessels.

The current study also demonstrated that mmu_circRNA_018351 was downregulated and had potential binding sites for mmu-miR-877-3p, mmu-miR-1903, mmu-miR-667-5p, mmu-miR-207 and mmu-miR-665-5p. Of these, miR-207 serves a protective role against cell apoptosis following isch-emic stroke (40). In addition, miR-665 regulates insulin-like growth factor 2 to protect from sevoflurane anesthesia-induced cognitive dysfunction through inhibiting apoptosis (41). Based on the bioinformatics analysis of the present study, it is hypothesized that miR-207 and miR-665 are upregulated as a consequence of the downregulated mmu_circRNA_018351 expression to repress the apoptosis of endothelial cells, smooth muscle cells and fibroblasts in hypoxia-induced PH.

Function and pathway analyses for the target genes of the circRNAs were performed to provide further insight into the mechanism of hypoxia-induced PH. Among the enriched GO terms in BP, the cellular nitrogen compound metabolic process was identified. Several studies have demonstrated that the level of nitric oxide (NO) is decreased in patients with PH (42). NO is an important vasodilator for pulmonary vessels, and a decrease in the levels of NO could increase the vascular resistance and vascular pressure (43). Although current therapy targeting the NO pathway improves the clinical outcome of patients with PH, the role of the NO pathway in the pathogenesis of pulmonary vascular remodeling in PH has not yet been fully elucidated (44). Thus, it is proposed that the circRNAs identified in the present study contribute to the pathophysiological effect of NO in PH and are potential clinical biomarkers.

Bioinformatics analysis also identified several hypoxia and vascular remodeling-associated pathways. The PI3K-Akt signaling pathway serves an important role in hypoxia-induced PH through promoting the proliferation of endothelial cells and smooth muscle cells (45,46). Furthermore, the Hippo, Wnt, AMPK and HIF-1 signaling pathways, which were identified in KEGG pathway analysis in the present study, have been reported to have a close association with the development of hypoxia-induced PH (47,48). These findings suggest that circRNAs may be key modulators in the molecular mechanism of hypoxia-induced PH and that, in the future, specific circRNAs could be tested as diagnostic and therapeutic targets of hypoxia-induced PH.

There are a number of limitations in the present study. Firstly, the differentially expressed circRNAs were identified only in a mouse model with hypoxia-induced PH and should be further validated in patients with PH. Secondly, the functions of the differentially expressed circRNAs were predicted using bioinformatics tools, and there were no further studies to demonstrate the role of candidate circRNAs in the pathological process of hypoxia-induced PH. Thirdly, the circRNA-miRNA-mRNA networks were constructed in the current study according to the premise that the majority of circRNAs function as miRNA sponges. However, the molecular functions of >99.9% of identified circRNAs remain unknown, and numerous other functions, such as translation and regulation of transcription and alternative splicing, have been reported (49). Thus, in the future, functional experiments should be conducted to elucidate the roles of these circRNAs in hypoxia-induced PH.

In conclusion, the present study is the first to profile circRNA expression in hypoxia-induced PH. In total, 23 significantly upregulated and 41 significantly downregulated circRNAs were identified. The circRNA-miRNA interactions were also predicted and circRNA-miRNA-mRNA networks were constructed for two leading candidate circRNAs. Furthermore, GO and pathway analyses revealed the potential role of the target genes of circRNAs in hypoxia-induced PH. The findings suggest that circRNAs serve a key role in the pathological process of hypoxia-induced PH. Thus, candidate circRNAs are potential targets for treatment and potential biomarkers for diagnosis of patients with hypoxia-induced PH.