Identification of miRNA–mRNA–TFs Regulatory Network and Crucial Pathways Involved in Tetralogy of Fallot

Tetralogy of Fallot (TOF) is the most common cyanotic congenital heart disease. However, its pathogenesis remains unknown. To explore key regulatory connections and crucial pathways underlying the TOF, gene or microRNA expression profile datasets of human TOF were obtained from the National Center for Biotechnology Information (NCBI) Gene Expression Omnibus (GEO) database. The differentially expressed mRNAs (DEmRNAs) and microRNAs (DEmiRs) between TOF and healthy groups were identified after data preprocessing, followed by Gene Ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis. Then, we further constructed protein–protein interaction (PPI) network and subnetwork of modules. Ultimately, to investigate the regulatory network underlying TOF, a global triple network including miRNAs, mRNAs, and transcription factors (TFs) was constructed based on the integrated data. In the present study, a total of 529 DEmRNAs, including 115 downregulated and 414 upregulated DEmRNAs, and 7 significantly upregulated DemiRs, including miR-499, miR-23b, miR-222, miR-1275, miR-93, miR-155, and miR-187, were found between TOF and control groups. Furthermore, 22 hub genes ranked by top 5% genes with high connectivity and six TFs, including SRF, CNOT4, SIX6, SRRM3, NELFA, and ONECUT3, were identified and might play crucial roles in the molecular pathogenesis of TOF. Additionally, an miRNA–mRNA–TF co-regulatory network was established and indicated ubiquitin-mediated proteolysis, energy metabolism associated pathways, neurodevelopmental disorder associated pathways, and ribosomes might be involved in the pathogenesis of TOF. The current research provides a comprehensive perspective of regulatory mechanism networks underlying TOF and also identifies potential molecule targets of genetic counseling and prenatal diagnosis for TOF.


INTRODUCTION
Cardiac development is an extremely complex process involving the integration of multiple cell lineages into the threedimensional organ and its seamless connection to the vascular system (Liu and Olson, 2010). Proper embryological heart development and function are precisely controlled by an evolutionarily conserved network of transcription factors (TFs), which are triggered by upstream signaling systems involved in cardiac myogenesis, morphogenesis, and contractility (Olson, 2006;Liu and Olson, 2010). Abnormalities in the networks during the cardiac formation process might lead to the failure of proper cardiac differentiation, migration, and apoptosis, which finally results in congenital heart disease (CHD).
Tetralogy of Fallot (TOF) is the most common cyanotic CHD, occurring in about 3 of every 10,000 live births and accounting for 7-10% of CHD (Bailliard and Anderson, 2009). TOF is caused by the non-uniform separation between bulbos arteriosus and truncus arteriosus during the early stage of embryonic development, resulting in malformations such as obstruction of the right ventricular outflow tract, ventricular septal defect (VSD), override of the ventricular septum by the aortic root, and right ventricular hypertrophy (Bailliard and Anderson, 2009;Wang et al., 2018). The clinical, pathological, phenotypical, and biologic heterogeneity of TOF suggests the complex pathogenesis. The precise molecular pathogenesis of TOF has not been completely identified. However, the latest research is advancing the understanding of the molecular mechanisms of TOF. Today, most investigations into the pathogenesis of TOF have focused on mutations in specific protein-coding genes, such as GATA4, NKX-2.5, JAG1, FOXC2, TBX5, TBX1, etc. (Lalani and Belmont, 2014). In addition, microdeletions of chromosome 22q11.2 and several other copy number variations, including15q11.2, 1q21.1, 2p13.3, 3p25.1, and 16p13.11, are found implicated in the pathological mechanisms of TOF (Lalani and Belmont, 2014).
Recently, several studies have revealed that microRNAs (miRNAs) play crucial roles in cardiac signaling and transcriptional pathway modulating multiple cardiac development, function, and disease (Liu and Olson, 2010). Cardiac development depends on the proper spatiotemporal expression of particular miRNAs. miR-1 was the first reported miRNA to be implicated in cardiac development, which targets DLL1 and HAND2 and further promotes the differentiation of embryonic stem cells into cardiac lineage (O'Brien et al., 2012). MiR-133a-1 and miR-133a-2 double-knockout mice have abnormal expression of SRF and CCND2, leading to late embryonic and neonatal lethality due to chamber dilatation and VSD (Liu et al., 2008). Deletion of either miR-106a-363 or miR-106b-25 combined with the miR-17∼92 null allele can result in embryonic lethality, accompanied by severe VSD, atrial septal defects, and thin-walled myocardium (Ventura et al., 2008). Taken together, these studies revealed that miRNAs play vital roles in vertebrate heart development. However, a comprehensive profiling of miRNA regulatory networks in TOF has not been reported to date.
The molecular mechanism and critical regulators involved in TOF are poorly deciphered. In the present study, to investigate the regulatory networks and key pathways underlying TOF, we constructed the protein-protein interaction (PPI) and an integrative miRNA-mRNA-TF regulatory network to further elucidate the key mechanisms underlying TOF.

Collection of Datasets
All of the TOF datasets were obtained from the Gene Expression Omnibus (GEO) (Clough and Barrett, 2016) with the following criteria: (1) mRNA or miRNA expression profiles by microarray or RNA-seq could be available; and (2) TOF patients and healthy controls were estimated. Consequently, a total of four microarray datasets [GSE35776 (O'Brien et al., 2012), GSE26125 (Bittel et al., 2011), GSE40128 (Zhang et al., 2013), and GSE35490 (Bittel et al., 2014)] and one RNA-seq dataset [GSE36761 (Grunert et al., 2014)] were included in this study, in which GSE40128 and GSE35490 were used for its miRNA expression profile. A total of 107 samples (32 controls and 75 cases) were analyzed in this study; details are provided in Table 1.

Identification of Differentially Expressed mRNAs and miRNAs
For microarray data, background correction, quantile normalization, and determination of expression levels were conducted using R/Affy (Gautier et al., 2004) or Oligo (Carvalho and Irizarry, 2010) packages. Differential expression analysis of mRNAs or miRNAs between the TOF and control groups was conducted using the R/Limma package (Diboun et al., 2006).
For RNA sequencing dataset (GSE36761), adapter sequences of raw reads were trimmed by Trimmomatic (Bolger et al., 2014). STAR was utilized to align the short reads to the reference human genome (hg38) (Dobin et al., 2013). Accepted mapped reads were then summarized at gene level counts using featureCounts  where its size is proportional to the number of input genes falling into that term, and its color represents its cluster identity (i.e., nodes of the same color belong to the same cluster). Terms with a similarity score >0.3 are linked by an edge (the thickness of the edge represents the similarity score). (B) The top 20 KEGG enrichment analysis for DEmRNAs. The numbers on the x-axis are the names of pathways. The numbers on the y-axis are gene counts.
FIGURE 2 | The protein-protein interaction (PPI) network of differentially expressed genes (DEmRNAs). Larger circles with red colors represent hub genes, and smaller hexagons with green colors represent non-hub genes.
software (Liao et al., 2014). Differential gene expression analysis between the TOF and control groups was conducted using R package DESeq2 (Love et al., 2014). The adjusted P < 0.05 and |log Fold Change (FC)| > 2 were set as the cut-off criteria for differentially expressed mRNAs (DEmRNAs), as well as differentially expressed miRNAs (DEmiRs) with thresholds of |log FC| > 1.5 and adjusted P < 0.05.

miRNA-Target Regulatory Network Analysis
miRWalk software is a comprehensive database that provides predicted and validated miRNA binding sites of known genes of human, rat, cow, mouse, and dog (Dweep and Gretz, 2015). Experimentally verified miRNA-gene regulatory pairs were considered as targets of DEmiRs. Then, we extracted the genes as significant DEmRNAs by intersecting the target genes of the DEmiRs and the overlapping genes of the DEmRNAs. The miRNA-target gene regulatory network was constructed using the Cytoscape software (Shannon et al., 2003).

Gene Ontology and Kyoto Encyclopedia of Genes and Genomes Pathway Enrichment Analysis
To further understand the gene function in TOF, Gene Ontology (GO) analysis (Ashburner et al., 2000) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis (Kanehisa and Goto, 2000) of the DEmRNAs were performed using Metascape , with P < 0.05 as the cutoff criterion.

PPI Network Construction and Identification of Hub Genes
The STRING (Search Tool for the Retrieval of Interacting Genes/Proteins) database provides a critical assessment and integration of PPI, including direct (physical) as well as indirect (functional) associations (Szklarczyk et al., 2015). The PPI networks of DEmRNAs were identified using the STRING database, and PPI score threshold was set as 0.4 to obtain interactions with higher confidence (>0.4). As nodes with a higher degree of connectivity contributing more to the stability of the network, the top 5% DEmRNAs ranked by the degree connectivity were defined as hub genes by using the plugin CentiScaPe (Scardoni et al., 2009). The Molecular Complex Detection (MCODE) (Bader and Hogue, 2003) plugin was applied to identify significant modules with node score cutoff = 0.2, K-Core = 4, max depth = 100 and degree cut-off = 4 set as the cut-off criteria. The PPI network was visualized by Cytoscape software (Shannon et al., 2003).

miRNAs-mRNAs-TFs Regulatory Network Analysis
MicroRNAs can regulate gene expression at the posttranscriptional stage (Slezak-Prochazka et al., 2010), whereas TFs can activate or repress transcription at a pretranscriptional stage (Zaret and Carroll, 2011). The miRNA-target regulatory network was visualized to depict interactions between miRNAs and their potential targets using Cytoscape 3.7.1.
The plugin iRegulon (Janky et al., 2014) in Cytoscape software, which included multiple human databases of the TF-target pairs such as Transfac and Encode, was used to identify TFs regulating DEmRNAs in the miRNA-target regulatory network. In this study, the parameters for TF motif enrichment analysis were set as the minimum identity between orthologous genes ≥0.05, maximum false discovery rate (FDR) on motif similarity ≤0.001, and TF motifs with normalized enrichment score (NES) >3 were considered as the threshold value for the selection of potential associations. TransmiR is a TF-miRNA regulation database, through which regulatory relations between TFs and miRNAs can be identified (Tong et al., 2019). Ultimately, the miRNA-mRNA--TF regulatory network was constructed utilizing Cytoscape 3.7.1.

Identification of DEmRNAs and DE miRNAs (DEmiRs) for TOF
Compared with the normal samples, following quality control and data normalization, a total of 1735 DEmRNAs including 1315 upregulated and 420 downregulated DEmRNAs were identified  Figure S1. We obtained 6261 targeted genes of the 7 significant DEmiRs from the miRWalk database. Four hundred and fourteen genes were found as the intersection genes between these target genes and the overall upregulated DEmRNAs. Therefore, a total of 529 DEmRNAs, including the 414 upregulated and 115 downregulated genes, were identified as the final sets of significant DEmRNAs (Supplementary Table S1).

GO and KEGG Pathway Enrichment Analysis of DEmRNAs
Functional and pathway enrichment analyses were applied to explore the pathogenesis in TOF. The DEmRNAs enriched according to the P-value are visualized in Figure 1A (Supplementary Table S2). The significantly top five GO Frontiers in Genetics | www.frontiersin.org terms in DEmRNAs were striated muscle cell differentiation, mitochondrial matrix, cofactor metabolic process, myofibril assembly, and generation of precursor metabolites and energy. The significant top five PPI pathways of DEmRNAs were ubiquitin-mediated proteolysis, non-alcoholic fatty liver disease, Parkinson's disease, Huntington's disease, and citrate cycle ( Figure 1B and Supplementary Table S3).

PPI Network Construction and Module Analysis
To investigate the relationship of the DEmRNAs, a PPI network was constructed utilizing Cytoscape software based on the STRING database. In total, 448 nodes and 1693 PPI relationships were obtained (Figure 2). Twenty-two genes ranked by top 5% genes with high connectivity were defined as hub genes for TOF, all of which were upregulated ( Table 2 and Supplementary  Table S4). The 22 hub genes interact directly with 337 edges. GO analysis results showed that hub genes were significantly enriched in ubiquitin-protein transferase activity and ubiquitin protein ligase activity in the molecular function ontology, cytosolic part and ubiquitin ligase complex in cellular components ontology, and protein polyubiquitination and cellular protein catabolic process in biological processes ontology. The results of KEGG analysis showed that hub genes were significantly enriched in ubiquitin-mediated proteolysis and protein processing in endoplasmic reticulum (ER) and ribosomes (Figure 3 and Supplementary Table S5). These hub genes might play a crucial role in TOF and should be studied further.
A total of three modules from the PPI network satisfied the criteria of an MCODE computed K-score = 4 and degree cut-off = 4 (Figure 4). Fifty-nine MCODE-related genes were identified and which contained the 17 hub genes. Module 1 (Figure 4A) included 21 nodes and 210 PPI pairs; KEGG analysis showed that these genes were significantly enriched in ubiquitinmediated proteolysis; and GO analysis showed that these genes were mainly associated with protein polyubiquitination, proteasomal protein catabolic process, contractile fiber part, sarcomere, and ubiquitin-protein transferase activity ( Figure 4D). Module 2 included 27 nodes and 84 PPI pairs (Figure 4B), and KEGG analysis showed that these genes were enriched in energy metabolism associated pathways (oxidative phosphorylation, citrate cycle, and carbon metabolism), and neurodevelopmental disorders (Parkinson's disease, Huntington's disease, and Alzheimer's disease), and GO terms were primarily in cellular respiration, organelle localization, oxidoreductase complex, mitochondrial protein complex, oxidoreductase activity, and ATPase activity ( Figure 4E). Module 3 included 11 nodes and 32 PPI pairs ( Figure 4C); KEGG analysis showed that these genes were significantly enriched in ribosomes, and GO analysis revealed that module 3 related genes were mainly associated with establishment of protein localization to organelles, translation initiation, cytosolic part chaperone complex, structural constituents of ribosomes, and unfolded protein binding ( Figure 4F). The detailed GO and KEGG pathway analysis information of three modules is shown in Supplementary Table S6.

miRNA-TF-Target Regulatory Network Analysis
To further understand the regulatory relationship between DEmiRs and DEmRNAs, the miRNA-target regulation network was constructed as shown in Figure 5. miR-222 had the highest connectivity with 529 target DEmRNAs. Moreover, several miRNAs were predicted to have common targets, for example, FIGURE 5 | miRNA-target regulatory networks of differentially expressed genes (DEmRNAs). The pink V-shapes represent miRNAs, red circles represent upregulated DEmRNAs, and green circles represent downregulated DEmRNAs. miR222 and miR93 have 139 common target genes, and miR222 and miR23b have 90 common target genes (Figure 5).

DISCUSSION
Heart development is an extremely complex process involving networks of various signaling pathways, the downstream TFs, protein-coding genes, and posttranscriptional regulation for cardiac myogenesis, morphogenesis, and contractility. Owing to the lack of specific molecule mechanisms, genetic counseling and prenatal diagnosis of TOF have been a major challenge. Therefore, it is imperative to explore the specific molecular mechanism for TOF, which will be valuable for risk assessment and even prognosis prediction. Recently, certain dysregulated mRNAs, miRNAs, and TFs have been reported in TOF individually. Notably, transcriptional and posttranscriptional regulation plays important roles in TOF development. In our study, an integrative analysis of microarray and RNAsequencing datasets was conducted to identify the regulatory network underlying TOF. A total of 529 DEmRNAs, 7 DEmiRs, and 6 TFs were identified as significant elements for the TOF. Based on 529 DEmRNAs, we obtained 22 hub genes and three significant modules from PPI networks using MCODE analysis. Eleven hub genes were significantly enriched in protein polyubiquitination, including SKP1, FBXO32, UBE3A, PSMC2, RPS27A, SIAH2, UBE2D3, UBE2N, FBXW11, DLD, and ZBTB16. The ubiquitin-proteasome system (UPS) is important for the regulation of cellular protein degradation and synthesis (Al-Hassnan et al., 2016). Impaired UPS has been involved in the etiology of human cardiovascular diseases, such as cardiomyopathy (CMP), heart failure, and cardiac ischemia. Skp1/Cul1/F-box (SCF) ubiquitin E3 ligase complex plays a crucial role in ubiquitination of cardiac proteins (Jang et al., 2011). FBXO32, which encodes a muscle-specific ubiquitin ligase, has been implicated in the pathogenesis of CMP through the UPS (Al-Hassnan et al., 2016). UBE3A is one of the important members of UPS; Ube3a mutant mice had deficits in Ca 2+ /calmodulin-dependent kinase II (CaMKII), while CaMKII played a vital role in cardiac hypertrophy (Cheng et al., 2019). ZBTB16, a Kruppel-like zinc finger protein, highly expressed in the heart, is an important TF that regulates cardiac hypertrophy through the Ang II receptor 2 in response to angiotensin II (Wang et al., 2012). GO function annotation also indicates that many hub genes (HSPA8, PSMC2, RPL7A, RPL8, RPS2, RPS24, RPS27A, HSPA9, TCP1, CCT5) are significantly enriched in regulation of transcription, posttranscriptional regulation and translational regulation, including mRNA catabolic process, RNA catabolic process, ribosomes, translational initiation, nucleartranscribed mRNA catabolic process, regulation of mRNA stability, DNA biosynthetic process, and rRNA processing. Elevated rates of protein synthesis have been demonstrated to contribute to cardiomyocyte hypertrophy. Accelerated rates of protein synthesis can be achieved by altering the translational efficiency of existing ribosomes and/or by accelerating the rate of synthesis of new ribosomes, indicating a crucial role for ribosome biogenesis in the regulation of cardiac growth (Brandenburger et al., 2003). Abnormal expression of ribosomes can affect the development and function of the cardiomyocyte and may lead to right ventricular hypertrophy of TOF. Members of the heat shock protein 70 family (HSPA8 and HSPA9) have also been proven to be essential to maintain cardiomyocyte protein homeostasis (Fang et al., 2017). Underexpression of HSPA8 and HSPA9 will disrupt cardiomyocyte homeostasis and interfere with critical cellular functions by small aggregates of proteins, which finally contributes to heart failure and hereditary heart diseases. However, the biochemical property and the molecular function of other hub genes have not been fully elucidated, and further functional studies are crucial to better understand the mechanism underlying the pathogenesis of TOF.
Gene Ontology analysis of DEmRNAs showed that DEmRNAs are mainly involved in striated muscle cell differentiation, mitochondrial matrix, cofactor metabolic process, myofibril assembly, and generation of precursor metabolites and energy, and these biological processes mainly contribute to cardiac cell division and differentiation. KEGG pathway analysis demonstrated abnormal molecular expression of several pathways may contribute to the pathogenesis of TOF, such as ubiquitin-mediated proteolysis pathways, energy metabolism associated pathways (citrate cycle, cysteine and methionine metabolism, glycolysis/gluconeogenesis, and PPAR signaling pathway), and pathways associated with neurodevelopmental disorders (Parkinson's disease, Huntington's disease, and Alzheimer's disease). Ubiquitin-mediated proteolysis is essential to maintain the integrity of cellular proteins that make up FIGURE 7 | miRNA-mRNA-TF regulatory network of differentially expressed genes (DEmRNAs). The pink diamonds represent miRNAs, orange triangles represent TFs, red circles represent upregulated genes, and green circles represent downregulated genes. TFs, transcription factors. sarcomeres, mitochondria, ER, and membranes and therefore is critical to cardiomyocyte survival and functioning (Pagan et al., 2013;Li et al., 2018). Insufficient ubiquitin-mediated proteolysis might be a pathogenic factor leading to cardiac pathogenesis (Li et al., 2018). Myocardial energy metabolism is an important determinant of cardiac structure and function (Finck, 2007). Energy metabolism associated pathways affect the normal cellular biological processes, including cell survival/death, differentiation/proliferation, and DNA repair (Dhar-Chowdhury et al., 2007;Hanahan and Weinberg, 2011;Tran and Wang, 2019). Pathological alterations of energy metabolism pathways could result in impaired signaling transduction, perturbed ion and redox homeostasis, and contractile dysfunction. Future studies in the specific role of energy metabolism will shed new light on the pathogenesis of TOF. Fifteen genes, including SLC25A4, ATP5F1A, ATP5F1B, NDUFA10, NDUFB5, POLR2E, REST, SDHD, SOD2, UQCRC2, LPL, VDAC1, VDAC2, UQCR11, and GRIN2A, were also enriched in neurodevelopmental disorders associated pathways. New research shows neurodevelopmental disorders occur in 10% of all children with CHD and in 50% of children with severe CHD (Marino et al., 2012), and de novo mutations account for at least 20% of CHD patients with associated extracardiac or neurodevelopmental abnormalities (Homsy et al., 2015), indicating a striking shared genetic etiology between TOF and neurodevelopmental disorders.
Among the six TFs predicted in the present study, SRF has been previously reported to contribute to differentiation, maturation, and homeostasis of cardiomyocytes, Srf knockout mice show severe defects vascular smooth muscle cells lacking organizing actin/intermediate filament bundles and die at the gastrulation stage (Schlesinger et al., 2011). CNOT4 functions as a ubiquitin-protein ligase (E3) (Deshaies and Joazeiro, 2009), interacting with a subset of ubiquitin−conjugating enzymes (E2s). SIX6 can mediate the repression of p27/Kip1, an inhibitor of the cell cycle G1/S transition, thereby regulating cellular organ−specific proliferation (Conte et al., 2005). SRRM3 is a alternative splicing regulator required for motor coordination (Nakano et al., 2019), potentially related to myocardial systolic motion. NELFA, a member of the NELF, cause promoter proximal pausing on the HSP70 gene (Wu et al., 2003), while HSPA8 and HSPA9 are identified hub genes of higher degrees. ONECUT3 is essential in multiple tissues for proper specification and development of unique cell types (Kropp and Gannon, 2016). However, whether these TFs are implicated in the pathogenesis of TOF can be investigated in future studies.
MicroRNAs, one of the most important endogenous epigenetic factors, inhibit posttranscriptional gene expression of target genes. Recently, several studies have revealed miRNAs as governing gene expression during cardiac development (Small and Olson, 2011). Seven upregulated DEmiRs were identified in the present study, including miR-499, miR-155, miR-23b, miR-222, miR-1275, miR-93, and miR-187. Among of them, miR-499, miR-155, miR-23b, and miR-222 have been demonstrated as participating in cardiac development and function. miR-499 regulates the mitochondrial fission machinery in cardiomyocytes by suppressing calcineurinmediated dephosphorylation of DRP1 . Similarly, miR-1275 was found to correlate with decreased mitochondrial membrane potential (Nymark et al., 2015). miR-155 is a critical microRNA that represses MEF2A expression (Seok et al., 2011), and Mef2a knockout mice exhibit dilated right ventricle, myofibrillar fragmentation, mitochondrial disorganization, and cardiac sudden death (Naya et al., 2002). miR-222 is upregulated in ventricular outflow tract tissues from infants with non-syndromic TOF, and miR-222 is involved in proliferation and differentiation of cardiomyocytes in in vitro experiments (Zhang et al., 2013). miR-23b was found upregulated in cyanotic CHD patients, and the upregulated miR-23b increases cardiomyocyte apoptosis and reduces cell growth under hypoxic conditions (Li et al., 2017). miR-93 was identified as highly expressed in the fetal mouse heart at four key time points (embryonic day E12.5, E14.5, E16.5, and E18.5) (Cao et al., 2012). miR-187 and its target genes Itpkc, Tbl1xr1, and Lrrc59 were reported to play regulatory roles in cardiomyocyte apoptosis and cardiac inflammation (Ektesabi, 2019). Therefore, we speculate that these miRNAs mainly contribute to abnormal cardiac cell division, differentiation, and apoptosis for TOF pathogenesis.

CONCLUSION
In summary, we performed an integrated analysis, based on four microarray datasets and one RNA-seq dataset, to identify specific TFs, miRNAs, and mRNAs that may play a pivotal role in TOF. Moreover, we established a TF-miRNA-mRNA network in TOF. However, there are certain limitations in our present study. The sample number of TOF specimens from patients was relatively small, which may reduce the credibility of the miRNA-mRNA-TF coregulatory network analysis. In addition, the present study lacked further experiments to verify the analysis results as a solid foundation. All in all, the current research provides a comprehensive perspective of regulatory mechanism networks underlying TOF and also identifies potential molecule targets of genetic counseling and prenatal diagnosis.

DATA AVAILABILITY STATEMENT
The datasets generated for this study can be found in the Gene Expression Omnibus (GEO): GSE35776, GSE36761, GSE26125, GSE40128, GSE35490.

AUTHOR CONTRIBUTIONS
FL conceived and designed the experiments. GY, BZ, BW, and QF analyzed and interpreted the data. GY and BZ wrote the manuscript. All the authors reviewed the manuscript.

SUPPLEMENTARY MATERIAL
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene. 2020.00552/full#supplementary-material FIGURE S1 | Volcano plots of differentially expressed genes or microRNAs. The x-axis represents the fold change (log-scaled), and the y-axis represents the adjusted P-value (log-scaled). Each symbol represents a different gene. The red color of the symbols means upregulated, whereas the green color of the symbols means downregulated.