The Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/) was searched to obtain the gene expression datasets by using the search terms ((cardiac hypertrophy OR cardiac remodeling OR cardiac fibrosis) AND angiotensin II). The obtained datasets were filtered according to the following inclusion criteria: (i) the organism was filtered by “mus musculus”; (ii) the study type was “Expression profiling by array”; (iii) mRNA expression profiles derived from the samples from both the Ang II-infused mice (cardiac remodeling) and the saline-infused mice (control); (iv) the raw data of the series matrix file(s) should be provided for reanalysis. Datasets that did not meet the aforementioned criteria were excluded. After careful review and comparison, the gene expression dataset GSE47420 was downloaded from the GEO database for subsequent analysis. GSE47420 was processed based on the GPL6887 platforms and contained heart samples of 3 Balb/C mice with Ang II-infused cardiac remodeling and 3 age-matched control mice [20]. Before analysis, these data were normalized using R’s limma package (R version 4.2.2 and limma version 3.54.0, Bell Labs, Lucent Technologies, Murray Hill, NJ, USA).

We first used GEO2R (https://www.ncbi.nlm.nih.gov/geo/geo2r/) to identify the differentially expressed genes (DEGs) between the Ang II-infused mice (cardiac remodeling) and the saline-infused mice (control) in GSE47420. GEO2R (version 1, National Center for Biotechnology Information, Bethesda, MD, USA) is an interactive web tool that allows for comparison between two groups of Samples from normalized GEO datasets. The Benjamini and Hochberg false discovery rate (FDR) method was used in GEO2R to adjust the p-value to reduce false positivity. Then the significant DEGs between the Ang II and control groups were obtained based on the following inclusion criteria: a threshold of adjusted p 0.05 and the absolute value of log2 fold change (log2FC) $>$0.585. Volcano maps and heatmaps of DEGs were visualized by TBtools (https://github.com/CJ-Chen/TBtools/releases) [21].

The pyroptosis-related genes were downloaded from the GeneCards database (https://www.genecards.org/). The GeneCards database is a searchable, integrative database that can be used to integrate gene-centric data from approximately 150 web sources, including genomic, transcriptomic, proteomic, genetic, clinical, and functional information. Pyroptosis-related genes were intersected with the total DEGs in GSE47420 to obtain the pyroptosis-related DEGs. Venn diagram and heatmaps of pyroptosis-related DEGs were visualized by TBtools.

The online biological tool Metascape (https://metascape.org/gp/index.html#/main/step1) was used to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of pyroptosis-related DEGs in GSE47420 [22]. GO enrichment analysis included biological process (BP), cellular component (CC), and molecular function (MF). The adjusted p value 0.05 was considered statistically significant. Bubble plots of GO and KEGG pathway enrichment analyses of pyroptosis-related DEGs were drawn using the online bioinformatics analysis platform SangerBox (http://www.sangerbox.com/) [23].

In order to investigate the mechanism of interaction among the pyroptosis-related DEGs in Ang II-induced cardiac remodeling, a PPI network was constructed by the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING) online database (http://stringdb.org/) and then visualized using Cytoscape v3.8.2 software (https://cytoscape.org/). The Cytohubba plug-in of Cytoscape software was used to screen out the top 15 key genes with high connectivity in this PPI network using the Maximal Clique Centrality (MCC) method. Pyroptosis-related hub genes were identified by using multiple algorithms including the MCC, Maximum Neighborhood Component (MNC), and Degree, and those genes ranked in the top 15 in at least three types of algorithms were considered to be hub genes for Ang II-induced cardiac remodeling.

To investigate the interaction between the upstream miRNAs and the pyroptosis-related hub genes, the upstream miRNAs that regulate the pyroptosis-related hub genes were predicted using the online StarBase v2.0 dataset based on one/multiple target-predicting programs (TargetScan, PicTar, RNA22, PITA, miRanda, and mirSVR) [24]. Subsequently, the predicted miRNAs were intersected with the differentially expressed miRNA (DE-miRNAs) in the dataset GSE69559 to obtain the pyroptosis-related DE-miRNAs. In GSE69559, the DE-miRNAs derived from the heart samples from both the Ang II-infused mice (cardiac remodeling) and the saline-infused mice (control), and there were three miRNAs upregulated and nineteen miRNAs downregulated [25], as shown in Supplementary Table 1. The obtained overlapping DE-miRNAs were filtered according to the following inclusion criteria: the upregulated DE-miRNAs in GSE69559 target the downregulated pyroptosis-related DEGs in GSE47420, whilst the downregulated DE-miRNAs target the upregulated pyroptosis-related DEGs. Subsequently, the miRNA-target gene pairs were constructed based on the potential interaction between the DE-miRNAs and their corresponding target genes. A miRNA-target gene interaction sub-network was visualized using Cytoscape software.

To investigate the interaction between the upstream transcription factors (TFs) and the pyroptosis-related DE-miRNAs, the upstream TFs that regulate the pyroptosis-related DE-miRNAs were predicted using the online TransmiR v2.0 dataset [26]. Currently, TransmiR v2.0 contains 3730 literature-curated TF-miRNA regulations, which cover ~623 TFs, ~785 miRNAs, 19 organisms, and 1349 publications. Subsequently, the predicted upstream TFs were intersected with the total DEGs in GSE47420 to obtain the differentially expressed TFs (DE-TFs). The TF-miRNA pairs were constructed based on the potential interaction between the DE-TFs and their corresponding target miRNAs. The predicted interaction between the DE-TFs and miRNAs by TransmiR2.0 includes four action types: regulation, activation, repression, and feedback. The obtained TF-miRNA pairs were filtered according to the following inclusion criteria: the TF-miRNA pair whose action type is activation or repression must be differentially expressed in terms of transcript expression in GSE47420 and GSE69559 in the same or opposite directions. A TF-miRNA interaction sub-network was visualized using Cytoscape software.

A TF-miRNA-target gene Interaction Network in Ang II-induced cardiac remodeling was constructed based on the TF-miRNA pairs and the miRNA-target gene pairs. The integrated TF-miRNA-target gene regulatory network was visualized by Cytoscape software. The integrated TF-miRNA-target gene regulatory network was filtered according to the following inclusion criteria: the miRNA-target gene pairs that had no predicted regulatory TFs will be discarded.

The correlation between the target genes and cardiac remodeling was evaluated by the online platform CTD (https://ctdbase.org/), scored, and visualized as bar graphs.

To screen potential drugs targeting cardiac remodeling, the identified target genes in our regulatory network were introduced into the Drug-Gene Interaction Database (DGIdb, https://dgidb.genome.wustl.edu/) as drug targets to search for potential drugs that target them. Those drugs that show well-known interaction types with the identified target genes were considered potential drugs for the treatment of cardiac remodeling.

In this study, an Ang II-induced cardiac remodeling in mice was used to analyze and compare the expression of key hub DEGs, DE-miRNAs, and DE-TFs. The animal experiments were approved by the Institutional Animal Use and Care Committee at the Hubei University of Science and Technology and followed the National Institutes of Health Guide for the Care and Use of Laboratory Animals (8th Ed., National Research Council, National Academy Press, Washington, DC, 2011). 12 healthy C57BL/6 male mice (8-week-old, 20–25 g), purchased from the Center for Disease Control of Hubei Province (Wuhan, China), were randomly assigned to two groups (n = 6 per group): the control group and the Ang II group. All animals were subcutaneously implanted with osmotic pumps (ALZET Model 2004, Durect Corporation, Cupertino, CA, USA.) to continuously infuse Ang II (1.5 mg/kg per day for 2 weeks, Ang II group) or equal volumes of PBS (vehicle, pH 7.2, control group) as described previously [27].

Mouse body weights were obtained before the animals were euthanized. At the end of the study, all the animals were killed by cervical dislocation, and the hearts were quickly dissected and weighed to calculate the ratio of heart weight to body weight (HW/BW, mg/kg). Heart tissues were quickly harvested, frozen, and stored at –80 °C until assays of gene expression analyses, or embedded in 4% paraformaldehyde for wheat germ agglutinin (WGA) staining.

Average CSA was measured by WGA staining. Left ventricular tissues were fixed with 4% paraformaldehyde, embedded in paraffin, and sectioned into 5-µm slices. After being deparaffinized and rehydrated according to the usual procedure, tissue sections were stained with an Alexa Fluor 488–conjugated WGA (20 µg/mL, Sigma-Aldrich) at 4 °C overnight to visualize cell membranes and DAPI for 15 min to stain cell nuclei. The average CSA from approximately 200 cells per slide was measured from the image captured with the fluorescence microscope and analyzed by Image-Pro Plus 5.0 software (Media Cybernetics, Inc., Silver Spring, MD, USA).

The total RNA was extracted from left ventricular tissues using TRIzol Reagent (cat. no. 15596018, Invitrogen, Carlsbad, CA, USA) according to the manufacturer’s instructions. For mRNA detection of DEAD-Box Helicase 3 X-Linked (DDX3X), Tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein zeta (YWHAZ, also named 14-3-3$\zeta$), CCAAT/enhancer-binding protein beta (CEBPB), and BCL11B, SuperScript™ III reverse transcriptase (cat. no. 18080–051, Invitrogen, Carlsbad, CA, USA) was utilized to synthesize cDNA. The synthesized cDNA was amplified via real-time polymerase chain reaction (ABI Q6, Applied Biosystems Inc., Waltham, MA, USA). The amplification reaction was carried out as follows: initial denaturation at 94 °C for 3 min followed by 35 cycles of denaturation at 94 °C for 30 s, annealing at 58 °C for 30 s, and extension at 72 °C for 45 s. Primer sequences used for DDX3X, YWHAZ, CEBPB, BCL11B, atrial natriuretic peptide (ANP), $\beta$-myosin heavy chain ($\beta$-MHC), and GAPDH were listed in Table 1.

For quantitative real-time PCR of miRNAs, total RNA was reversed-transcribed to cDNA by a TaqMan miRNA reverse transcription kit (Applied Biosystems, CA, USA) using specific stem-loop primers. After single-stranded cDNA was synthesized, quantitative real-time PCR (qRT-PCR) amplifications were performed using an iCycleriQ Real-Time Detection System (Bio-Rad Laboratories, USA). Primer sequences used for miR-16-5p, miR-22-3p, and U6 were listed in Table 1. PCR amplification protocol consisted of initial denaturation at 94 °C for 5 min followed by 35 cycles of denaturation at 94 °C for 10 s, annealing at 55 °C for 30 s, and extension at 72 °C for 30 s. Data were calculated by the 2${}^{-\mathrm{\Delta }\mathrm{\Delta }\text{CT}}$ method [28] and normalized to GAPDH before being compared with the control. All qRT-PCRs were carried out in triplicate.

All data were presented as mean $\pm$ S.E.M. Statistical analysis for animal experiments was performed using GraphPad Prism 5.0 (GraphPad Software, San Diego, CA, USA). The Student’s t-test was used to compare gene expression between the control group and the Ang II group. Differences were considered significant when p 0.05. Pearson correlation coefficients were estimated for correlation analysis. When the GEO2R database was used to identify the differentially expressed mRNAs and miRNAs, a significant level of the statistical analysis was set at adjusted p 0.05. When the Metascape database was used to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis, a significant level of the statistical analysis was set at p 0.01.

The dataset GSE47420 contained the expression data of mRNA transcripts from the Ang II-infused mice and the control mice. The online GEO2R tool was used to identify the differentially expressed genes (DEGs) between the two groups. Based on the inclusion criteria (adjusted p 0.05 and the absolute value of log2 fold change (log2FC) $>$0.585), a total of 1289 DEGs (796 upregulated and 493 downregulated) were identified. The DEGs were visualized in the volcano plots (Fig. 2A) and heatmap (Fig. 2B) using the bioinformatics analysis tool TBtools. Subsequently, 372 pyroptosis-related genes were downloaded from the GeneCards database (Supplementary Table 2) and intersected with the total DEGs in GSE47420 to obtain the pyroptosis-related DEGs. A total of 32 pyroptosis-related DEGs (20 upregulated and 12 downregulated) were identified based on the aforementioned intersection (Fig. 2C, Table 2). A heatmap was generated to reveal the expression profiles of the 32 pyroptosis-related DEGs using TBtools (Fig. 2D).

Fig. 2.

Identification of differentially expressed genes (DEGs) in Ang II-induced cardiac remodeling. (A) Volcano plot of total DEGs in GSE47420 associated with Ang II-induced cardiac remodeling. (B) Heatmap of total DEGs in GSE47420. (C) Venn diagram demonstrating the intersection between total DEGs in GSE47420 and pyroptosis-related genes from GeneCards. (D) Heatmap of 32 pyroptosis-related DEGs in GSE47420. The red arrow indicated pyroptosis-related DEGs, blue arrows indicated downregulated pyroptosis-related DEGs.

Gene ontology (GO) annotation and KEGG pathway enrichment analysis were used to deduce the function of the 32 pyroptosis-related DEGs. The results of GO and KEGG enrichment analyses were shown in Fig. 3. For the BP category of GO enrichment, the pyroptosis-related DEGs were mainly enriched in regulation of defense response, negative regulation of cellular catabolic process, negative regulation of catabolic process, negative regulation of response to external stimulus, regulation of response to biotic stimulus, etc. (Fig. 3A). The pyroptosis-related DEGs were enriched in the CC category of GO analysis such as transcription regulator complex, perinuclear region of cytoplasm, chromatin, extracellular matrix, and external encapsulating structure (Fig. 3B), and MF category of GO analysis such as transcription factor binding, ubiquitin-like protein ligase binding, ubiquitin protein ligase binding, DNA binding transcription factor binding, and protein kinase binding (Fig. 3C). KEGG pathway enrichment analysis revealed that the pyroptosis-related DEGs were enriched in NOD-like receptor signaling pathway, salmonella infection, PI${}_3$K-Akt signaling pathway, pathway in cancer, IL-17 signaling pathway, necroptosis, etc. (Fig. 3D).

Fig. 3.

GO (Gene Ontology) and KEGG pathway enrichment analysis of 32 pyroptosis-related DEGs, and PPI network. (A) Top 10 GO terms for biological process (BP). (B) Top 10 GO terms for cellular component (CC). (C) Top 10 GO terms for molecular function (MF). (D) Top 10 KEGG pathways. (E) PPI networks visualized using Cytoscape software and hub genes. Red round rectangles indicated representative GO enrichment and KEGG pathway enrichment terms. GO, gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; DEGs, differentially expressed genes; PPI, protein-protein interaction.

To further investigate the protein-protein interaction (PPI), the 32 pyroptosis-related DEGs were introduced into the String database to construct a PPI network with the cut-off criterion of confidence score $>$0.4. As shown in Fig. 3E, 27 of 32 pyroptosis-related DEGs were eventually identified to generate a PPI network with 27 nodes and 65 edges. Unrelated five genes (HTRA1, FAT1, FNDC5, TCEA3, and TUBB6) in the PPI network were not included in further analysis. Next, the PPI network was visualized using Cytoscape software. Based on at least three types of algorithms (MCC, MNC, Degree), the top-ranked 13 genes (STAT3, HSP90AB1, HSP90AA1, MAPK14, PARP1, MDM2, YWHAZ, EED, ANXA2, CEBPB, ETS1, DDX3X, and ELAVL1) were screened out using Cytoscape software as the pyroptosis-related hub DEGs (Fig. 3E, Table 3).

To generate a miRNA-target gene network, the 13 pyroptosis-related hub DEGs were introduced into the online StarBase v2.0 dataset to predict the upstream miRNAs that target these hub DEGs. A total of 300 miRNAs were predicted using the StarBase v2.0 dataset based on one/multiple target-predicting programs (TargetScan, PicTar, RNA22, PITA, miRanda, and mirSVR). Subsequently, the predicted 300 miRNAs were intersected with the 22 differentially expressed miRNAs (DE-miRNAs) in GSE69559 to obtain 10 overlapping pyroptosis-related DE-miRNAs (Fig. 4A, Table 4). In these overlapping DE-miRNAs, miR-208b-3p was upregulated in GSE69559, whilst the other nine miRNAs were downregulated. The upregulated miR-208b-3p was predicted to target two pyroptosis-related hub DEGs, PARP1 and ETS1. The downregulated miR-191-5p was predicted to target upregulated CEBPB. The other eight upregulated miRNAs (let-7g-5p, let-7i-5p, miR-133a-3p, miR-146a-5p, miR-16-5p, miR-22-3p, miR-30b-5p, miR-378b) were predicted to target five pyroptosis-related hub DEGs (DDX3X, ELAVL1, YWHAZ, STAT3, EED). These detailed miRNA-target gene regulatory pairs were outlined in Table 4. Subsequently, the miRNA-target gene pairs were integrated to generate a miRNA-target gene sub-network and visualized using Cytoscape software (Fig. 4B). Three unrelated miRNA-target gene pairs (miR-208b-3p-PARP1, miR-208b-3p-ETS1, miR-191-5p-CEBPB) were eliminated in further analysis. Thus, the miRNA-target gene sub-network eventually contained eight DE-miRNAs, five target DEGs related to pyroptosis, and fourteen interactions between the DE-miRNAs and the target DEGs.

Fig. 4.

Construction of TF-miRNA-target gene regulatory network. (A) Venn diagram demonstrating the intersection between total DE-miRNAs in GSE69559 and the predicted miRNAs by StarBase2.0. (B) miRNA-target gene sub-network. (C) Venn diagram demonstrating the intersection between total DEGs in GSE47420 and the predicted TFs by TransmiR2.0. (D) TF-miRNA sub-network. (E) TF-miRNA-target gene regulatory network. Blue nodes indicate the target genes that were targeted by miRNAs, green nodes indicate the miRNAs and pink nodes indicate the TFs. Red-capped lines indicated overlapping differentially expressed miRNAs (DE-miRNAs) and TFs (DE-TFs). TF, transcription factor; miRNA, microRNA; DE-miRNA, differentially expressed microRNA; DEG, differentially expressed gene.

To generate a TF-miRNA network, the eight DE-miRNAs that targeted five pyroptosis-related hub DEGs in the miRNA-target gene sub-network were introduced into the online TransmiR v2.0 dataset to predict the upstream TFs that interact with these DE-miRNAs. A total of 265 TFs were predicted to interact with the eight DE-miRNAs. Subsequently, the predicted 265 TFs were intersected with the total DEGs in GSE47420 to obtain 16 overlapping differentially expressed TFs (DE-TFs) (Fig. 4C, Table 5). Subsequently, the TF-miRNA pairs were constructed based on the potential interaction between the overlapping DE-TFs and their corresponding target DE-miRNAs. The TF-miRNA pairs were integrated to generate a TF-miRNA sub-network and visualized using Cytoscape software. The TF-miRNA sub-network eventually contained sixteen DE-TFs, six DE-miRNAs, and thirty-five interactions between the DE-TFs and their target DE-miRNAs (Fig. 4D).

Finally, the miRNA-target gene sub-network and the TF-miRNA sub-network were merged using the network merge operation of Cytoscape to build a TF-miRNA-target gene regulatory network (Fig. 4E). The TF-miRNA-target gene regulatory network was composed of 26 nodes, 49 edges, and two types of regulatory relationships, including the TF-miRNA and miRNA-target gene actions.

In addition, the topological properties of the regulatory network were analyzed with the NetworkAnalyzer plugin in the Cytoscape. The top-ranked four key target genes, miRNAs, and TFs with high degrees were listed in Table 6. To ensure the rationality of the finding, the expression of the top-ranked CEBPB-miR-16-5p-DDX3X and BCL11B-miR-22-3p-YWHAZ regulatory axes were further verified in the following animal experiments using qRT-PCR.

Left ventricular hypertrophy (LVH) is a characteristic component of left ventricular remodeling. To confirm the functional roles of the target genes in cardiac remodeling, we employed the CTD database to confirm their correlation with cardiovascular diseases via evaluating the LVH score. The genes with an LVH score $>$15 are identified as significantly correlated genes and are visualized as red bars in Fig. 5. Among the above-mentioned five target genes, DDX3X, ELAVL1, YWHAZ, STAT3, and EED scored 31.09, 16.22, 36.13, 81.7, and 16.11, respectively, and were considered to be crucial cardiac remodeling-related genes.

Fig. 5.

Left ventricular hypertrophy (LVH) score determined by the CTD database. The correlation between the genesis of LVH and DDX3X (A), ELAVL1 (B), YWHAZ (C), STAT3 (D), and EED (E) was evaluated by the online platform Comparative Toxicogenomics Database (CTD), scored, and shown as bar graphs. Red rectangles indicated the LVH scores.

To explore the interrelationship between the central genetic aspects of cardiac remodeling and the available drugs, the drug-target gene interaction network was predicted using the DGIdb database. As shown in Table 7 (Ref.[29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40]), 26 drugs targeting 3 hub genes (ELAVL1, STAT3, and EED) were identified. Except for 5 known types, most drug-gene interaction types were unknown. The 5 known interaction types were all inhibitory for ELAVL1 and STAT3. The 5 drugs (dehydromutactin, celecoxib, Azd-1480, ochromycinone, and cosmomycin C) with inhibitory effects on ELAVL1 and STAT3 were considered as potential drugs for the treatment of cardiac remodeling. STAT3 was screened for a wide range of relationships with multiple drugs, which may provide potential targets for the treatment and prognosis of cardiac remodeling.

To validate the key nodes and regulatory relationship in cardiac remodeling, qRT-PCR was used to analyze the expression of the core CEBPB-miR-16-5p-DDX3X and BCL11B-miR-22-3p-YWHAZ regulatory axes in an Ang II-induced cardiac remodeling in mice. In this study, the 8-week-old C57BL/6 male mice were subcutaneously implanted with osmotic pumps to continuously infuse Ang II (1.5 mg/kg per day for 2 weeks) to establish an Ang II-induced cardiac remodeling, as evidenced by the elevated ratios of heart weight to body weight (HW/BW), increased average CSA, and enhanced mRNA expression of hypertrophic genes (ANP, $\beta$-MHC) in the Ang II-infused mice (cardiac remodeling) compared with the saline-infused mice (control, Fig. 6A–E). qRT-PCR analysis of the genes in the left ventricular tissues revealed that the mRNA expressions of the pyroptosis-related hub genes DDX3X (p = 0.034) and YWHAZ (p = 0.018) and transcription factor CEBPB (p = 0.021) were significantly upregulated, while the transcript expressions of BCL11B (p = 0.020), miR-16-5p (p = 0.016), and miR-22-3p (p = 0.023) were significantly downregulated in the Ang II-infused group than in the control group (Fig. 6F–K). The results were consistent with those obtained via bioinformatic analysis. Additionally, Pearson correlation analysis revealed that miR-16-5p showed a significant negative correlation with DDX3X (r = –0.79, p = 0.008) and CEBPB (r = –0.86, p 0.001), while miR-22-3p exhibited a significant negative correlation with YWHAZ (r = –0.74, p = 0.009) and positive correlation with BCL11B (r = 0.85, p 0.001) (Fig. 6L).

Fig. 6.

Verification of key differentially expressed genes related to Ang II-induced cardiac remodeling and pyroptosis. (A) Representative WGA staining for cross-sectional area (CSA) measurement. Bars = 20 µm. (B) Quantification of cardiomyocyte CSA. Averaged data results from 20 fields per heart. (C) Quantification of the ratio of heart weight to body weight (HW/BW). (D,E) Expression of ANP and $\beta$-MHC mRNA. (F,I) Expression of target genes DDX3X and YWHAZ mRNA. (G,J) Transcript Expression of miR16-5p and miR-22-3p. (H,K) Expression of TFs CEBPB and BCL11B mRNA. (L) Pearson correlation analysis of mRNAs and miRNAs. Red circles indicated a correlation between miR-16-5p and DDX3X/CEBPB. Blue circles indicated a correlation between miR-22-3p and YWHAZ/BCL11B. *p 0.05; **p 0.01; ***p 0.001; ****p 0.0001.

Since the changing trend of most predicted key nodes in our animal experiments was similar to that in the microarray analysis in GSE47420 and GSE69559, our results revealed that the relationships between the predicted miRNAs and their target genes or the corresponding TFs in the TF-miRNA-target gene network were reliable.

Although we for the first time investigated the potential TF-miRNA–target gene regulatory network in Ang II-induced cardiac remodeling by integrating multiple microarray datasets from the GEO platform, there were several limitations existing in the present study. Firstly, the sample size of the GEO dataset included in this study is not big enough, only containing 6 tissue samples. Secondly, construction of these miRNA-target gene pairs and TF-miRNA pairs in the regulatory network was only based on predictions from public databases and lacking experimental validation of these pairs in this study. Thirdly, we did not perform in vitro and in vivo function experiments on these miRNA-mRNA regulatory pathways in cardiac remodeling.

Future investigations with more tissue samples and corresponding experimental validations in vivo and in vitro are still needed to validate our in silico analysis. Additionally, the clinical relevance of this study is not clear and requires further studies.