Construction of a potential microRNA, transcription factor and mRNA regulatory network in hepatocellular carcinoma.

BACKGROUND: Hepatocellular carcinoma (HCC) is one of the most common malignancies worldwide and the third leading cause of cancer-related death. MicroRNAs and transcription factors (TFs) cooperate to regulate the same target gene, thus affecting the progression of HCC.
METHODS: Differentially expressed miRNAs and mRNAs were screened. Functional enrichment analysis of these HCC-related mRNAs was performed, and a protein-protein interaction network was constructed. TFs that regulate these miRNAs and hub genes were also screened.
RESULTS: Ten differentially upregulated miRNAs and 5 differentially downregulated miRNAs were screened. Additionally, 183 downregulated mRNAs and 303 upregulated mRNAs that are potentially bound to these differentially expressed miRNAs were identified. The Kyoto Encyclopedia of Genes and Genomes (KEGG) results showed that the differentially expressed mRNAs were significantly enriched in pathways in cancer, the Wnt signaling pathway, and the Rap1 signaling pathway. Then, 220 TFs were identified for 5 candidate genes of the downregulated mRNAs, and 258 TFs were identified for 9 candidate genes of the upregulated mRNAs. Finally, the 9 upregulated hub genes were related to higher overall survival (OS) in the low-expression group, and 4/5 downregulated hub genes were related to higher OS in the high-expression group.
CONCLUSIONS: This study constructed a potential regulatory network between candidate molecules and that need to be further verified. These regulatory relationships are expected to clarify the new molecular mechanisms of the occurrence and development of HCC. 2020 Translational Cancer Research. All rights reserved.

Introduction

Hepatocellular carcinoma (HCC) is the fifth most common malignant cancer worldwide, and its incidence has been increasing (1-3). Due to the absence of obvious symptoms at the initial stage of HCC, patients are often already at an advanced stage at the time of diagnosis, so radical resection cannot be performed, and the prognosis is very poor. In recent decades, great progress has been made in the early diagnosis and treatment of HCC. For patients who cannot undergo surgical resection, chemoembolization and systemic symptomatic treatment is adopted (2). The specific mechanism of HCC has not been clarified, and no effective therapeutic target has been found, so the long-term prognosis is still not ideal. As a molecular targeted drug, the multikinase inhibitor sorafenib has benefited some patients with liver cancer, revealing the importance of molecular therapy in HCC, although the benefit has been limited so far (4-6). Therefore, the exploration of markers with high sensitivity and specificity in the process of HCC is of great significance for early diagnosis and treatment.

MicroRNAs (miRNAs) are a class of noncoding single-stranded RNA molecules mainly involved in the posttranscriptional regulation of target genes (7). Specifically, miRNAs directly degrade mRNAs or inhibit translation through complementary base pairing (complete or incomplete) binding to the 3‘ untranslated region (3'UTR) of mRNA (7,8). Previous studies have indicated that miRNAs are involved in the occurrence and development of tumors by affecting the expression level of target genes and activating or inhibiting relevant signaling pathways (7,9). Because the transcription level of miRNAs is correlated with the prognosis of patients, miRNAs have been used as tumor biomarkers (10). More encouragingly, miRNAs that inhibit tumor progression (such as miR-34 and miR-122) have been used for targeted therapy, and clinical trials have been conducted (8,10). Besides, Wu et al. combined deep sequencing data of 15 pairs liver cancer tissue samples and TCGA data to screen out 9 pairs of miRNA-mRNAs, and explored their potential ability to predict the progression and prognosis of HCC (11). Interestingly, Tu et al. found different miRNA-mRNA pairs using 5 datasets from the Gene Expression Omnibus (GEO) database (12), suggesting that we might be able to mine novel potentially valuable molecules by combining different datasets. Transcription factors (TFs) bind to specific nucleotide sequences upstream of genes and regulate transcription by regulating the binding of RNA polymerases to DNA templates. TFs and miRNAs cooperate to regulate the same target gene: TFs bind to the promoter region where gene transcription is initiated, and miRNAs perform modifications at the posttranscriptional level (13-15). TFs play important roles in regulating the progression of HCC, and inducing the differentiation of key cell TFs may be a potential tumor therapy strategy (16,17). It has been confirmed in studies that miRNAs and TFs can regulate each other and thus affect the occurrence and development of cancers (18-22). Therefore, the exploration of target genes regulated by miRNAs and upstream TFs and the establishment of a mutual regulatory network are expected to clarify the molecular mechanisms of HCC progression and provide new therapeutic strategies.

In the current study, differentially upregulated and downregulated miRNAs were screened after integrating The Cancer Genome Atlas (TCGA) database and the GSE36915 dataset. In addition, upregulated and downregulated mRNAs in the TCGA database that are potentially bound to these differentially expressed miRNAs were identified. Then, candidate hub genes of the differentially downregulated and upregulated genes were identified. TFs that regulate these miRNAs and hub genes were also screened, and further biological analysis of these candidate molecules was performed to construct a network of mutually regulated relationships between them. Finally, 5 downregulated hub genes (ESR1, FOXO1, IGF1, SERPINE1, CAT) and 9 upregulated genes (CDK1, CCNE2, CHEK1, CCNE1, CDC25B, CDC6, CDC25A, CCNF, BIRC5) were found to serve as tumor biomarkers, the correlation between the differential expression of the hub genes and the prognosis of HCC patients was also demonstrated.

Methods

HCC-related mRNA and miRNAs

The RNA-Seq and miRNA-Seq data from normal samples and HCC tumor samples were downloaded from the TCGA database (https://portal.gdc.cancer.gov/). Then, the limma package and edgeR package were used to identify differentially expressed mRNAs and miRNAs in R (21,22). The mRNAs and miRNAs with a P value <0.05 and log fold change (logFC) >1 were chosen.

In addition, the miRNA expression profiles of HCC and adjacent tissues were also obtained from the GSE36915 dataset in the National Center of Biotechnology Information GEO database (https://www.ncbi.nlm.nih.gov/geo/). The web tool GEO2R was used to identify HCC-related miRNAs (https://www.ncbi.nlm.nih.gov/geo/geo2r/). The miRNAs with a P value <0.05 and logFC >1 were chosen, and the intersecting miRNAs in the TCGA and GEO databases were selected for further analysis. MiRNAs common to both datasets were defined as HCC-related miRNAs in this study. Venn diagrams were performed by a web tool (23).

Construction of a protein-protein interaction (PPI) network

The target genes of the candidate HCC-related miRNAs were predicted by using the Cytoscape tool, which links to miRTarBase and TargetScan for the predictions. Next, the miRNA-mRNA pairs were identified by taking the intersection with potential target genes and the differentially expressed genes (DEGs) of the TCGA database. The selected mRNAs were uploaded to the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING) database (http://string-db.com/), an online database for the retrieval of interacting genes and the construction of PPI networks, and then analyzed by using the multiple protein online tool. The interactions with a combined score greater than 0.4 were identified and visualized by Cytoscape, and the hub genes in the PPI network were screened using Molecular Complex Detection (MCODE).

Gene Ontology (GO) annotation and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of the HCC-related genes

The selected mRNAs were subjected to GO and KEGG analysis using the Database for Annotation, Visualization, and Integrated Discovery (DAVID) (https://david-d.ncifcrf.gov/). A P value <0.05 was identified as statistically significant, and the top 10 pathways in the gene count were visualized.

Construction of a TF-related network

To predict the upstream TFs, the hub genes were uploaded to NetworkAnalyst (https://www.networkanalyst.ca/) for target screening. In this online platform, three TF-gene interaction databases (ENCODE, JASPAR and ChEA) were included. Only a peak intensity signal <500 and a score <1 are used in the ENCODE database. JASPAR is a TF binding site profile database. ChEA is a database that contains information inferred from integrating literature-curated Chip-X data. Additionally, FunRich software (24) was used to predict the potential TFs of the HCC-related miRNAs, with a P value <0.05 as the threshold to select TFs for further analysis.

The hub gene-TF pairs and miRNA-TF pairs were visualized. Then, TFs that can regulate both the hub gene and miRNA were identified. Finally, a potential regulatory network was constructed according to the miRNA-hub gene pairs, miRNA-TF pairs and hub gene-TF pairs.

The Kaplan-Meier plotter analysis

The Kaplan-Meier plotter (http://kmplot.com/analysis/), with the data from the TCGA project and GEO, is capable of assessing the effect of a variety of genes on survival in 21 cancer types. In this study, to understand the prognostic role of the hub genes in HCC patients, Kaplan-Meier survival curves were constructed, and the cutoff values for the most statistically significant differences were used to divide HCC patients into a high-expression group and a low-expression group. A log-rank P value <0.05 was considered statistically significant.

Results

HCC-related mRNAs and miRNAs

We screened HCC-related differentially expressed mRNAs and miRNAs to identify the pair of interactions, and the results were shown in . In the TCGA dataset, 24 downregulated and 232 upregulated miRNAs were found in 374 tumor samples and 50 normal samples (). A total of 53 downregulated and 14 upregulated miRNAs were found in the GSE36915 dataset (), which contained 68 tumor samples and 21 nontumor samples. Further analysis showed that 5 downregulated miRNAs (hsa-miR-424, hsa-miR-1258, hsa-miR-511, hsa-miR-195, and hsa-miR-214) and 10 upregulated miRNAs (hsa-miR-10b, hsa-miR-224, hsa-miR-183, hsa-miR-452, hsa-miR-96, hsa-miR-182, hsa-miR-1180, hsa-miR-877, hsa-miR-190b, and hsa-miR-551b) overlapped in TCGA dataset and GSE36915 dataset (). These 15 miRNAs were defined as HCC-related miRNAs and used for the prediction of target genes ().

Figure 1

Hepatocellular carcinoma-related mRNAs and miRNAs. (A) The volcano plot of the differentially expressed miRNAs from The Cancer Genome Atlas (TCGA). (B) The volcano plot of the differentially expressed miRNAs from GSE36915. (C) The volcano plot of the differentially expressed mRNAs from TCGA. (D) The intersection of the downregulated miRNAs in TCGA and GSE36915. (E) The intersection of the upregulated miRNAs in TCGA and GSE36915. (F) The intersection of the upregulated mRNAs in TCGA and potential targets of the downregulated miRNAs. (G) The intersection of the downregulated mRNAs in TCGA and potential targets of the upregulated miRNAs. (H) The potential targets of the downregulated miRNAs. (I) The potential targets of the upregulated miRNAs.

Figure S1

The Kaplan-Meier survival curve of up-regulated and down-regulated miRNA in hepatocellular carcinoma patients.

In this study, potential target genes of the HCC-related miRNAs were screened from miRTarBase and TargetScan by using Cytoscape software. For the 5 downregulated miRNAs, 2,056 candidate mRNAs were identified and are shown in . Moreover, 3,592 genes were predicted to potentially target the 10 upregulated miRNAs. In addition, 1,053 downregulated mRNAs and 3,794 upregulated mRNAs were identified in the TCGA dataset (). Further analysis showed that 183 of 1,053 downregulated mRNAs were targeted by the 10 upregulated miRNAs (), and 303 of 3,794 upregulated mRNAs were targeted by the 5 downregulated miRNAs (). These 486 mRNAs were defined as HCC-related mRNAs for further analysis.

GO annotation and KEGG pathway enrichment analysis of the HCC-related genes

GO enrichment analysis for the 183 downregulated mRNAs and 303 upregulated mRNAs was performed to identify the significant biological processes (BPs), molecular functions (MFs), and cellular components (CCs). A P value <0.05 was identified as statistically significant, and the top 10 pathways in the gene count are shown in . In detail, for BP analysis, the downregulated mRNAs were mostly enriched in the positive regulation of transcription from the RNA polymerase II promoter, which was the same for the upregulated mRNAs (). CC analysis showed that these downregulated and upregulated mRNAs were mostly enriched in the plasma membrane (). MF analysis demonstrated that the downregulated mRNAs were significantly enriched in TF activity and sequence-specific DNA binding (), and the upregulated mRNAs were significantly enriched in metal ion binding ().

Figure 2

Gene Ontology annotation and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analysis. (A,B) The top 10 enriched biological processes items of the downregulated and upregulated candidate genes. (C,D) The top 10 enriched cellular components items of the downregulated and upregulated candidate genes. (E,F) The top 10 enriched molecular functions items of the downregulated and upregulated candidate genes. (G,H) The top 10 enriched Kyoto Encyclopedia of Genes and Genomes pathways of the downregulated and upregulated candidate genes.

Additionally, KEGG pathway enrichment analysis was performed, and the results showed that the downregulated mRNAs were significantly enriched in pathways in cancer, proteoglycans in cancer, and the Wnt signaling pathway (), while the upregulated mRNAs were significantly enriched in proteins in cancer, the Rap1 signaling pathway, and valine, leucine and isoleucine degradation ().

Construction of the PPI network

PPI networks were further constructed for the candidate target genes of the downregulated and upregulated mRNAs using STRING (). The results demonstrated that 118 of 183 downregulated mRNAs and 188 of 303 upregulated mRNAs were found to interact with each other (). For further module analysis, these interacting pairs were uploaded to Cytoscape. The module with the top score according to MCODE was identified and is shown in . The module of the downregulated PPI network contains 5 hub genes (ESR1, FOXO1, IGF1, SERPINE1, and CAT), and the module of the upregulated PPI network contains 9 hub genes (CDK1, CCNE2, CHEK1, CCNE1, CDC25B, CDC6, CDC25A, CCNF, and BIRC5).

Figure 3

The protein-protein interaction network. (A) The protein-protein interaction network of the upregulated mRNAs. (B) The hub genes of the upregulated mRNAs. (C) The protein-protein interaction network of the downregulated mRNAs. (D) The hub genes of the downregulated mRNAs.

Construction of the hub gene-TF-miRNA-related network

To further understand the roles of the hub genes in HCC, the upstream TFs of these 14 hub genes were predicted by using NetworkAnalyst. A total of 220 TFs were found for the 5 candidate genes of the downregulated mRNAs, and 258 TFs were found for the 9 candidate genes of the upregulated mRNAs. HNF4A and PPARG regulate all five downregulated hub genes, and MYC has the largest number of targets of upregulated genes. The results were visualized by Cytoscape and are shown in , .

Figure 4

Construction of the hub gene-transcription factor-miRNA-related network. (A) The network of the upregulated hub genes and potential target transcription factors. (B) The network of the downregulated hub genes and potential target transcription factors. (C) The network of the downregulated miRNAs and potential target transcription factors. (D) The network of the upregulated miRNAs and potential target transcription factors. (E) The network of the upregulated hub genes, downregulated miRNAs, and transcription factors. (F) The network of the downregulated hub genes, upregulated miRNAs, and transcription factors.

Table 1

Potential target transcription factors of the upregulated hub genes

Transcription factors	Numbers	Hub genes
MYC	9	BIRC5, CDC6, CDC25A, CDK1, CCNF, CCNE2, CCNE1, CHEK1, CDC25B
E2F1	8	CDC6, CDC25A, CDK1, CCNE2, CCNE1, CHEK1, CDC25B, CCNF
KDM5B	8	BIRC5, CDC6, CDC25A, CDK1, CCNF, CCNE2, CCNE1, CHEK1
MYCN	7	BIRC5, CDC6, CDC25A, CCNF, CCNE2, CCNE1, CHEK1
SOX2	7	BIRC5, CDC25A, CDK1, CCNE2, CCNE1, CDC25B, CHEK1
KLF4	7	BIRC5, CDK1, CCNF, CCNE2, CCNE1, CHEK1, CDC25B
CREM	7	BIRC5, CDC6, CDC25A, CCNF, CCNE2, CCNE1, CHEK1
CREB1	7	BIRC5, CDC6, CDC25A, CCNF, CCNE1, CHEK1, CDC25B
FOXC1	6	CDK1, CDC25A, CCNE2, CHEK1, CDC6, CCNE1
SAP30	6	BIRC5, CCNE1, CCNE2, CDC25A, CDC6, CHEK1
PHF8	6	CDC25A, CCNE2, BIRC5CCNE1, CDC6, CHEK1
CEBPB	6	CDC6, CDC25A, CCNF, CDC25B, CHEK1, CCNE2
PPARG	6	CDC6, CDK1, BIRC5, CCNE2CDC25B, CHEK1

Table 2

Potential target transcription factors of the downregulated hub genes

Transcription factors	Numbers	Hub genes
HNF4A	5	ESR1, FOXO1, IGF1, SERPINE1, CAT
PPARG	5	ESR1, FOXO1, IGF1, SERPINE1, CAT
NFE2L2	4	ESR1, IGF1, SERPINE1, CAT
SMAD4	4	ESR1, FOXO1, IGF1, SERPINE1
EGR1	4	ESR1, FOXO1, SERPINE1, CAT
STAT3	4	ESR1, FOXO1, IGF1, SERPINE1
TP63	4	FOXO1, IGF1, SERPINE1, CAT
E2F1	4	FOXO1, IGF1, SERPINE1, CAT
TCF4	3	ESR1, IGF1, CAT
RUNX1	3	ESR1, SERPINE1, CAT
FOXA2	3	ESR1, IGF1, SERPINE1
BACH1	3	ESR1, FOXO1, IGF1
CREB1	3	ESR1, FOXO1, IGF1
FLI1	3	FOXO1, IGF1, SERPINE1
SPI1	3	FOXO1, SERPINE1, CAT
SOX2	3	FOXO1, IGF1, SERPINE1
TFAP2A	3	ESR1, FOXO1, IGF1
GATA2	3	ESR1, SERPINE1, CAT
RAD21	3	FOXO1, IGF1, CAT
ZNF354C	3	ESR1, IGF1, SERPINE1
FOXC1	3	ESR1, IGF1, CAT

To screen the upstream TFs of the HCC-related miRNAs, the 5 downregulated miRNAs and 10 upregulated miRNAs were uploaded to FunRich. Two of the 5 downregulated miRNAs were found to have 19 significant potential TFs (). Seven of the 10 upregulated miRNAs were found to have 74 significant potential TFs (). Six miRNAs (hsa_miR_214, hsa_miR_224, hsa_miR_1258, hsa_miR_511, hsa_miR_1180 and hsa_miR_877) were not found to have a potential regulatory TF. However, we find that hsa_miR_214 and EIRC5, hsa_miR_224 and SERPINE1 show potential regulatory relationships in NetworkAnalyst databases. Thus, although no TF was found that might regulate hsa_miR_214 and hsa_miR_224, they were also shown in the .

In this study, the TFs that could modulate hub genes and miRNAs were selected, and networks were constructed. As shown in , the downregulated mRNA network consists of 5 downregulated hub genes, 7 HCC-related miRNAs and 16 TFs, and the upregulated mRNA network consists of 9 upregulated hub genes, 3 HCC-related miRNAs and 6 TFs ().

Figure S2

The Kaplan-Meier survival curve of the TFs that interacted with both up-regulated and down-regulated miRNA in hepatocellular carcinoma patients. TF, transcription factor.

Correlation between the hub gene transcripts and the prognosis of patients

As seen from the Kaplan-Meier survival curve in , the cutoff values for the most statistically significant differences were used to divide HCC patients into the high-expression group and the low-expression group. The overall survival of the 9 differentially upregulated hub genes was higher in the low-expression group than in the high-expression group [hazard ratio (HR): 1.28–3.44, P<0.05]. The gene with the highest HR value was CDC6 [HR: 2.44 (1.7–3.5), P=5.6e-07]. Conversely, 4 of the 5 differentially downregulated hub genes had a higher overall survival (OS) in the high-expression group than in the low-expression group (HR: 0.28–0.82, P<0.05).

Figure 5

The Kaplan-Meier survival curve of the transcription of hub genes.

Discussion

According to previous research, the posttranscriptional regulation of target genes by miRNAs directly degrades mRNA or inhibits translation through complementary base pairing (complete or incomplete), binding to the 3'UTR of the mRNA (7,8). So far, the miRNA-gene regulatory network in HCC has not been fully elucidated. Some of the differentially expressed miRNAs selected in our study were confirmed by previous studies. The downregulation of miR-195 leads to the accelerated progression of HCC due to the direct binding of CCNE1 and CDC25A, which affect the cell cycle (23). The downregulation of miR-195 can promote the angiogenesis and invasion of HCC, and the upregulation of miR-195 may potentially improve prognosis (24,25). The downregulation of miR-214 promoted tumor angiogenesis (26) and was associated with tumor cell invasion and HCC recurrence (27). The downregulation of miR-424 inhibits the migration and invasion of HCC, suggesting a better prognosis (28,29). MiR-224 is upregulated in HCC and affects key BPs such as apoptosis and proliferation (30,31). Furthermore, research has shown that the high expression of miR-224 in plasma can be used as a noninvasive biomarker to predict the early development of HCC (32). MiR-10b is upregulated in tumor tissues and promotes the progression and metastasis of HCC through a variety of pathways, such as CAMD1 (33), RhoC, uPAR and MMPs (34). MiR-183/96/182 were upregulated in cancer tissues and directly bound to FOXO1, promoting tumor progression and serving as a biomarker for poor prognosis (35).

It can be seen that −log10(P value) of hsa-miR-424 in TCGA data is greater than that in GSE36915 dataset in this study. There was also more miRNA significantly up-regulated in TCGA dataset than GSE36915. The possible reason is that TCGA has a larger number of samples than GSE36915, because the TCGA dataset contains sequencing data of 371 tumor samples and 50 normal samples, while there are only 72 tumor tissues and 21 normal tissues in GSE36915. In addition, the difference of detection technology (HCC data in TCGA from high-throughput sequencing, and data from miRNA arrays in GEO) is also an important potential factor. Besides, hsa-miR-424 has been demonstrated to be down-regulated in liver cancer (36), non-small cell lung cancer (37), skin hemangioma (38), osteosarcoma (39), breast cancer (40), etc., playing an suppressor role.

As seen from the above findings, both upregulated and downregulated miRNAs are directly involved in the occurrence and progression of HCC and can be used as biomarkers for prognosis. On the one hand, a miRNA can bind multiple mRNAs and regulate different molecular signaling pathways, affecting the progression of HCC. The innovation of our study was the construction of a network and discovery of new potential miRNA-mRNA pairs. Moreover, some pairs have been demonstrated in other cancers. For example, miR-424 can bind CDK1 and reduce translation, activate the Hippo pathway and extracellular signaling-mediated kinase pathway, and inhibit the proliferation of breast cancer cells (41). Xie et al. (41) and Tu et al. (12). identified potentially valuable miRNAs using different datasets. This suggests that combining different databases to screen new potentially valuable factors is a feasible approach. We combined the data from GSE36915 and TCGA and screened out novel miRNA-mRNA pairs. Therefore, the role of miRNA in HCC needs further exploration and verification.

TFs play important roles in regulating the progression of HCC, and inducing the differentiation of key cell TFs may be a potential tumor therapy strategy (16,17). It has been confirmed in studies that miRNAs and TFs can regulate each other and thus affect the occurrence and development of cancers (18-22). In our study, FOXA1, MEF2A and NKX6-1 were mutually regulated with differentially expressed miRNAs. Previous studies have shown that miR-182 and miR-183 can negatively regulate EGR1 in breast cancer (42). In pituitary adenoma cells, EGR1 binds to the promoter of miR-183, reducing its transcription (43).

HNF4A and PPARG regulate all five downregulated hub genes in the TF-related regulatory network and have been shown to be involved in the development and progression of HCC. HNF4A was found to induce HCC cells to differentiate into mature cells, thus inhibiting the occurrence and development of HCC (16,17). For TFs in combination with upregulated hub genes, MYC had the largest number of targets. MYC is a widely studied oncogenic gene. At the same time, as a TF, it can also show anticancer activity by inducing apoptosis, cell transformation and cell aging (44,45).

Our study screened 183 downregulated mRNAs and 303 upregulated mRNAs. KEGG enrichment analysis showed that pathways in cancer, proteoglycans in cancer, the Wnt signaling pathway, the cell cycle and the Hippo signaling pathway are the primary significant pathways for the upregulated genes. The downregulated genes are mainly enriched in proteoglycans in cancer and the Rap1 signaling pathway. Previous studies have shown that Wnts and their downstream effectors regulate various processes that are important for tumor initiation and growth, cell aging and death, differentiation, and metastasis (46,47). The Hippo signaling pathway is associated with liver cell division and differentiation as well as tumorigenesis (48). Rap1 and Ras have high sequence similarity; in other words, they have overlapping binding sites. Rap1 has been shown to both oppose and mimic Ras-driven cancer progression by synergistically initiating and sustaining ERK signaling (49,50). These studies support that the pathways of the differentially expressed hub genes we screened are closely related to the occurrence and development of HCC. In addition, the abnormal expression of these genes may predict poor prognosis and may serve as a relevant target for tumor therapy.

We screened out 5 hub genes with low expression in tumor tissues and high expression in adjacent tissues. The Kaplan-Meier survival curve showed that HCC patients with a high expression of the hub genes (except SERPINE1) had a better prognosis. This finding suggests that these genes may serve as biomarkers for the prognosis of HCC. Therefore, it is reasonable to speculate that these downregulated genes may be tumor-suppressor genes whose decreased expression in tumors leads to tumor progression and a worse prognosis in patients. As reported by previous studies, ESR1 may mediate the susceptibility of hepatitis B virus (HBV) carriers to HCC (51). The high expression of miR-182 and low expression of FOXO1 activated the PI3K/AKT signaling pathway to promote the proliferation of hepatic fibrosis cells (52). SERPINE1 behaves differently, and further research may be needed to clarify its function. The HR values of the 9 hub genes (CDC25A, CCNE1, CDK1, CCNE2, CHEK1, CDC25B, CDC6, CCNF, and BIRC5) with high expression in tumor tissues and low expression in para-cancerous tissues ranged from 1.83 to 2.40 (P<0.05), indicating that their high expression could be used as biomarkers for poor prognosis. Therefore, it can be speculated that the upregulation of these hub genes may be the main cause of tumor progression. It has been proven that upregulated CDC25A can be used as an independent prognostic marker for poor prognosis in HCC (53). CCNE1 overexpression is associated with hepatocyte chromosomes and genomes, increasing the risk of developing HCC in HBV patients (54).

By integrating all the results obtained, a new potential mutual regulation network was constructed by miRNA-TF pairs, hub gene-TF pairs and miRNA-TF pairs. Most miRNAs, mRNAs, and TFs have been shown to affect HCC progression, but their mutual regulatory networks have not been fully elucidated. The innovation of this study was the exploration of the potential regulatory relationships between these candidate molecules related to HCC based on the constantly enriched big data platform. Moreover, only a few of these screened regulatory relationships have been verified in HCC, and the remaining regulatory relationships need to be further verified. These regulatory relationships are expected to clarify the newly identified molecular mechanisms of the occurrence and development of HCC and develop new therapeutic targets.

This study also has some shortcomings. Our studies were all based on sequencing results and chip results, without experimental verification with clinical samples. In addition, the direct regulatory relationships among the screened miRNAs, mRNAs and TFs have not been verified. Finally, the functional verification of related factors needs further experimental study. Therefore, our future work involves designing experiments to verify the regulatory relations and functions of candidate molecules.

Conclusions

This study constructed a potential regulatory network between candidate molecules and that need to be further verified. These regulatory relationships are expected to clarify the new molecular mechanisms of the occurrence and development of HCC.