Screening of differentially expressed circRNAs
We employed the Arraystar Human circRNA Array Analysis to identify TOF-related circRNAs. The box plot showed the distribution of circRNA expression profiling. After the normalization, the distributions of log2 ratios across different samples were displayed in Figure. 1A. We then used the scatter-plot to identify differentially expressed circRNAs between TOF group and control group (Figure. 1B). Unsupervised hierarchical clustering analysis shows the dysregulated circRNA expression pattern in clinical samples (Figure. 1C).
The general data, including the functional classification and chromosome location of these differentially expressed circRNAs were summarized. Overall, 276 circRNAs were found to be significantly differentially expressed (|logFC|-value >2.0, P-value <0.05) (Supplementary Table 1). Compared to control group, 214 circRNAs, including 183 exonic, 13 intronic, 16 sense overlapping, 1 antisense and 1 intergenic were up-regulated and 62 circRNAs, including 53 exonic, 3 intronic, 5 sense overlapping and 1 antisense were down-regulated in CHD group (Figure. 2A & 2B). Figure. 2C and 2D demonstrated that the up‑ and down-regulated circRNAs were located in human chromosomes.
CircRNA-miRNA-mRNA network in TOF
A total of 793 miRNA response elements (MREs) of the differentially expressed circRNAs were identified by utilizing the miRNA target prediction tool derived from Arraystar (Figure.3A). Next, to consolidate the identified MREs that are associated with TOF and screen the more significant miRNA, we used GSE35490 from the NCBI/GEO database to screen 123 differentially expressed miRNAs in TOF group compared with health groups. (|logFC|-value >1.5, P-value <0.05) (Figure.3B). Using venn diagram analysis, 25 shared miRNAs (key miRNA signature) were obtained and presented in Venn diagram (Figure.3C).
In order to predict the key target mRNAs, we futher uploaded the abovementioned 25 key miRNAs to the miRDB, miRTarBase and TargetScan databases for searching. 34 key mRNAs that interacted with 8 of the 25 key miRNAs in all 3 datasets were selected. After removing the remaining 17 key miRNA and the corresponding circRNAs, 19 circRNAs, 8 miRNAs and 34 mRNAs were finally obtained to construct the ceRNA network related to TOF (Figure 4; Supplementary Table 2).
Functional GO terms and pathway enrichment analyses
To characterize the functional consequences of our identified genes, we performed an enrichment analysis of the ChIP-X (ChEA) database with Enrichr for all of the 34 key genes involved in the ceRNA network. The results revealed that GO terms with the most statistical significance was heart development (P-value=0.002538) (Figure.5A).
In addition, we characterized the putative pathogenic gene of TOF based on the CTD database. A total of 34 genes had curated or inferred association with “Congenital Heart Defects” ranked by inference score. Among these genes, HIF1A was the most significant marker gene (Figure.5B), confirmed by three studies [11, 26, 27]. To make our results more reliable, we futher assessed the association of 34 genes with TOF by using public databases GSE35776 [28]. Unsupervised clustering resulted in two distinct subgroups according to the genes expression model (Figure.5C).
Up-regulation of hsa_circ_0007798 expression in TOF heart tissues
To validate the previous study results, we expanded the sample size. The expression levels of hsa_circ_0007798 in 10 TOF heart tissues and 10 health heart tissues were measured by qPCR method, which showed that hsa_circ_0007798 expression was significantly up-regulated in TOF group (P < 0.001, Figure.6A). By means of inquiring from circBase (http://www.circbase.org) , we knew that hsa_circ_0007798 was located at chromosome 6 and was composed of 6 exons ontaining the full 805-bp sequence (Figure.6B). CircRNAs mostly act by competitively sponging miRNAs to regulate the expression of their target genes according to ceRNA theory. Therefore, basd on our screening strategy, we firstly consider hsa_circ_0007798/miR-199b-5p/HIF1A (Figure.6C). In addition, the results of the correlation analysis of qPCR verified the positvie linear relationship between hsa_circ_0007798 and HIF1A (pearson r coefficient = 0.6291) (Figure.6D). These results could make it possible for us to understand how predicted and screened key circRNAs-miRNA-mRNA are related to TOF progression.