Integrative module analysis of HCC gene expression landscapes.

Despite hepatocellular carcinoma (HCC) being a common cancer globally, its initiation and progression are not well understood. The present study was designed to investigate the hub genes and biological processes of HCC, which change substantially during its progression. Three gene expression profiles of 480 patients with HCC were obtained from the Gene Expression Omnibus database. Subsequent to performing functional annotations and constructing protein-protein interaction (PPI) networks, 657 differentially expressed genes were identified, which were subsequently used to screen candidate hub genes. PPI networks were modularized using the weighted gene correlation network analysis algorithm, the topological overlapping matrix and the hierarchical cluster tree, which were utilized via STRING. Clinical data obtained from The Cancer Genome Atlas were then analyzed to validate the experiments performed using six hub genes. Additionally, a transcription factor and microRNA-mRNA network were constructed to determine the potential regulatory mechanisms of six hub genes. The results revealed that the oxidation-reduction process and cell cycle associated processes were markedly involved in HCC progression. Six highly expressed genes, including cyclin B2, cell division cycle 20, mitotic arrest deficient 2 like 1, minichromosome maintenance complex component 2, centromere protein F and BUB mitotic checkpoint serine/threonine kinase B, were confirmed as hub genes and validated via experiments associated with cell division. These hub genes are necessary for confirmatory experiments and may be used in clinical gene therapy as biomarkers or drug targets. Copyright: © Li et al.

Introduction

Hepatocellular carcinoma (HCC) is the third leading cause of cancer-associated mortality globally and ranks as the 5th and 7th most common cancer type occurring in men and women, respectively (1,2). It is mainly attributed to hepatitis B virus (HBV) or hepatitis C virus (HCV) infections and patients with cirrhosis are more likely to develop HCC when compared with those without (3,4). The exact molecular mechanism underlying HCC has not yet been fully elucidated. However, liver transplantation remains the most effective treatment due to the poor prognosis and ineffective therapy for the disease (5,6). Recurrence rates following liver transplantation remain high and result in patient mortality in a majority of cases, despite a united therapeutic approach being utilized (7,8). The research and development of effective drugs for HCC therapy remains excessive and time-intensive. Therefore, novel treatments of HCC are urgently required.

Although a number of studies have investigated the molecular mechanisms of HCC (9–11), the effective molecular targets of drugs and biomarkers, in addition to the biological processes and signaling pathways of the disease, remain obscure. Therefore, the present study selected 480 chip datasets of HBV- and HCV-associated HCC tumor tissues and liver cirrhosis non-tumor tissues (samples without cirrhosis) from the Gene Expression Omnibus (GEO) database for bioinformatics analysis.

The limma package was used to filter the differentially expressed genes (DEGs), which were assessed via Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses for significantly expressed biological processes and signaling pathways. Furthermore, Gene Set Enrichment Analysis (GSEA) was applied to evaluate the statistically significant results at the transcriptome level. DEGs were also used to construct a protein-protein interaction (PPI) network via the STRING database to identify hub genes. The weighted gene correlation network analysis (WGCNA) algorithm was utilized to adjust the PPI network to model the dynamics of proteome changes. Finally, the clinical outcomes of hub genes that are involved in HCC progression were assessed via The Cancer Genome Atlas (TCGA) survival analysis and experiments. Overall, the present study was designed to develop a new method to find the biomarkers of HCC, which serve important functions in cancer detection and treatment.

Materials and methods

Microarray data source and pre-processing

The gene expression profiles of HCC were constructed by collecting the following three data sets, which were based on the Affymetrix HT HG-U133A and HG-U133A 2.0 Array (National Center for Biotechnology Information GEO database): GSE14323 (12), GSE14520 (13) and GSE17967 (14,15). A total of 480 biochips obtained from patients with resected HCC were analyzed and included 237 HCC tumor samples and 243 non-tumor samples. The raw data of three sets were downloaded from GEO and read via Simpleaffy, which is used as an R package for Affymetrix quality control and data analysis (16). Annotations were made using gene symbols and their platform annotations, following which a united gene expression matrix, including all 480 samples, was compiled. The mean value of gene expression was used in multiple probe sets with one gene symbol allocated. Normalization and batch rectification were performed prior to analysis.

Differential expression and functional analysis

DEGs were screened out using the limma package with the Bayesian method in R (17), and Log fold change (FC)>1 and an adjusted P-value of P<10−5 were considered to indicate a statistically significant difference. GO (18,19) and KEGG (20) analyses were performed using the Database for Annotation, Visualization and Integrated Discovery (DAVID) to investigate the significant biological processes and signaling pathways associated with DEGs (21,22). The enriched results of DAVID were presented using GOPlot, a package that visualizes the functional analysis of omics data (23).

GSEA

GSEA calculates whether a priori defined gene set is statistically significant, and determines concordant differences among biological processes (24). In the present study, GSEA was used to evaluate the differences between HCC tissues and adjacent non-tumor tissues in the whole transcriptome, beyond the DEGs, to avoid individual bias. Permutations of 105 were serviced in the progress. Leading edge analysis was performed for the enriched core genes of HCC.

Co-expression and PPI network analysis

To cluster the functional genes of DEGs, a PPI network was constructed from the STRING database (25), and WGCNA (26) was also performed using the dynamic tree cut package (27) at a minimum height of 0.2 for each module. A topological overlapping matrix (28) was utilized to screen networks. Pearson's correlation coefficient was used to collect the eigengene and interactors of DEGs. Finally, the annotations of each module were completed using clusterProfiler (29) and visualized in Cytoscape V_3.6.0 software (National Institute of General Medical Sciences, Bethesda, MD, USA) (30).

TCGA survival analysis

It is important to validate biomarkers for clinical outcomes in cancer prognosis. To perform survival analysis and risk assessment, the clinical data and expression profiles of HCC were downloaded from TCGA database via TCGAbiolinks (31).

Prediction of candidate gene transcription factors (TFs) and microRNAs (miRNAs)

The transcription factors of six candidate genes were determined using ENCODE Chip-seq data using human transcription factor information included in NetworkAnalyst (32,33).

The miRNA of five target genes were predicted using miRWalk 3.0 (http://mirwalk.umm.uni-heidelberg.de/) (34). The filter criterion for miRNA were confirmed using >3 of the following databases: miRwalk, miRanda (35), RNA22 (36) and TargetScan (37). The TF-gene and miRNA-gene network were visualized in Cytoscape software separately.

Patients and sample collection

A total of 47 HCC samples and their paired adjacent normal tissues were collected between February 2014 and May 2017 at the People's Hospital of Liaoning Province (Liaoning, China). Following surgical resection, the tissues were stored at −80°C. All specimens were subsequently evaluated by two independent pathologists via hematoxylin and eosin staining. The present study was approved by the Ethics Committee of the People's Hospital of Liaoning Province, and all patients provided written informed consent.

Reverse transcription-quantitative PCR (RT-qPCR)

Total RNA was extracted from 47 HCC and adjacent samples using TRIzol® (Invitrogen; Thermo Fisher Scientific, Inc.). RT was performed at 50°C for 3 min to synthesize cDNA. The six hub genes identified via bioinformatics were analyzed via RT-qPCR using a One-Step qPCR kit (cat. no. 11746100; Invitrogen; Thermo Fisher Scientific, Inc.) and a CFX Connect™ Real-Time System (Bio-Rad Laboratories, Inc.). The following thermocycling conditions were used for qPCR: Initial denaturation at 95°C for 15 min; 40 cycles of denaturation at 95°C for 10 sec, annealing at 60°C for 30 sec, and extension at 72°C for 20 sec; with 2X SYBR®-Green reaction mix (Invitrogen; Thermo Fisher Scientific, Inc.) RT-qPCR data were analyzed using the 2−ΔΔCq method (38), with β-actin used as a reference gene. The following primers were used in the current study: CCNB2 forward, 5′-CCGACGGTGTCCAGTGATTT-3′ and reverse, 5′-TGTTGTTTTGGTGGGTTGAACT-3′; CDC20 forward, 5′-GCACAGTTCGCGTTCGAGA-3′ and reverse, 5′-CTGGATTTGCCAGGAGTTCGG-3′; MAD2L1 forward, 5′-GTTCTTCTCATTCGGCATCAACA-3′ and reverse, 5′-GAGTCCGTATTTCTGCACTCG-3′; MCM2 forward, 5′-ATGGCGGAATCATCGGAATCC-3′ and reverse, 5′-GGTGAGGGCATCAGTACGC-3′; CENPF forward, 5′-CTCTCCCGTCAACAGCGTTC-3′ and reverse, 5′-GTTGTGCATATTCTTGGCTTGC-3′; BUB1B forward, 5′-AAATGACCCTCTGGATGTTTGG-3′ and reverse, 5′-GCATAAACGCCCTAATTTAAGCC-3′); and β-actin forward, 5′-CATGTACGTTGCTATCCAGGC-3′ and reverse, 5′-CTCCTTAATGTCACGCACGAT-3′.

Western blotting

HCC tissue samples were lysed using Tissue Extraction Reagent I (Invitrogen; Thermo Fisher Scientific, Inc.) with protease and phosphatase inhibitors. Protein concentrations were determined using a bicinchoninic acid assay kit (Thermo Fisher Scientific, Inc.). A total of 20 µg lysate proteins of every sample were separated using 6% (for CENPF) and 12% (for others) SDS-PAGE and transferred onto polyvinylidene difluoride membranes (Invitrogen; Thermo Fisher Scientific, Inc.). The membranes were blocked by 5% nonfat dry milk (Thermo Fisher Scientific, Inc.) for 30 min at room temperature. All antibodies were diluted at 1:1,000 in 5% nonfat dry milk (Thermo Fisher Scientific, Inc.), 1X TBS and 0.1% Tween-20 (Thermo Fisher Scientific, Inc.) at 4°C with gentle shaking, overnight. The following antibodies were used: CCNB2 (cat. no. ab10839; Abcam), CDC20 (cat. no. 4823; Cell Signaling Technology, Inc.), MAD2L1 (cat. no. 4636; Cell Signaling Technology, Inc.), MCM2 (cat. no. 4007; Cell Signaling Technology, Inc.), CENPF (cat. no. 58982; Cell Signaling Technology, Inc.), BUB1B (cat. no. 4116; Cell Signaling Technology, Inc.) and β-actin (cat. no. 4970; Cell Signaling Technology, Inc.). The membranes were subsequently incubated with horseradish peroxidase-conjugated secondary antibodies (cat. no. 7074S; anti-rabbit IgG; species cross-reactivity: Human, mouse, rat; and cat. no. 7076S; anti-mouse IgG; species cross-reactivity: Human, mouse, rat; Cell Signaling Technology, Inc.) at room temperature for 1 h. The immunoblots were visualized with an enhanced chemiluminescence detection reagent (EMD Millipore). Band intensities were quantified using the Image Lab 2.0 software (Bio-Rad Laboratories, Inc.).

Statistical analysis

An paired Student's t-test was utilized for comparisons between two groups (tumor tissue and its adjacent tissue). Data are presented as the mean ± standard deviation, except when indicated otherwise. The Kaplan-Meier method was used to perform survival analysis, and two subgroups were compared with Breslow. P<0.05 was considered to indicate a statistically significant difference. All analyses in the present study were performed with R software (version 3.5.0, http://www.R-project.org).

Results

The raw data of three data sets were downloaded from GEO and read into R using the Simpleaffy package for quality control and normalization (16). The combat method was used to eliminate the batch effects (39). A total of 480 gene expression profiles, obtained from HCC and non-tumor tissues, were analyzed in the present study.

HCC DEGs and functional analysis

A heatmap and volcano plot were constructed to present the variation in DEGs between HCC and non-tumor tissues (Fig. 1A). In total, 657 genes, including 386 upregulated genes and 271 downregulated genes, were differentially expressed in HCC, with a log2FC>2 and an adjusted P-value of P<10−5.

Figure 1.

Heat map and functional enrichment of differentially expressed genes in hepatocellular carcinoma. (A) Plots demonstrating the expression values of all top 800 changed genes detected by a microarray. (B) Gene Ontology enrichment of first seven functions. Colors represent log fold change. (C) A total of 24 statistically significant pathways are listed and their colors represent P-values. KEGG, Kyoto Encyclopedia of Genes and Genomes.

The molecular mechanism underlying HCC progression is yet to be fully elucidated. Therefore, the present study aimed to determine the functions of hub genes in HCC by analyzing associated DEGs. GO and KEGG pathway analyses were also performed to understand this process. The results of GO revealed that the most enriched GO targets were those involved in oxidation-reduction processes, the epoxygenase P450 pathway and cell cycle associated processes (Fig. 1B). Additionally, it was determined that a number of the genes involved in oxidation-reduction processes were downregulated. KEGG analysis demonstrated that enrichments were associated with the biosynthesis of antibiotics, fatty acid degradation, the cell cycle and DNA replication (Fig. 1C).

To confirm the functional associated genes in the whole transcriptome, as opposed to DEGs alone, GSEA was performed using the expression matrix between HCC tissue and their adjacent non-tumor types. HCC was singularly associated with downregulated genes associated with the oxidation-reduction process and fatty acid derivative metabolic processes. Furthermore, numerous upregulated genes were enriched in the cell cycle process, the regulation of cell division and the G2M checkpoint, which is concordant with the data obtained from DEG GO and KEGG enrichment analysis (Fig. 2A). GSEA leading edge analysis also determined that CDC20, nucleolar and spindle associated protein 1, NIMA related kinase 2, kinesin family member 11, baculoviral IAP repeat containing 5, TTK protein kinase and CENPF were exhibited in five gene sets, with a further six genes appearing in four gene sets (Fig. 2B and C).

Figure 2.

Gene Set Enrichment Analysis of whole transcriptome expression profiling of hepatocellular carcinoma. (A) GO enrichment plot-associated cell cycle gene sets. The normalized enrichment score and the false discovery rates are presented. Each member of gene set was represented by a bar at the bottom of plot. (B) Leading edge analysis among gene sets. The overlap of each subset was demonstrated by color intensity; the darker the color, the greater the overlap between the subsets. (C) The bar graph represents each gene and the number of subsets in which it appears. GO, Gene Ontology.

Integrative analysis elucidates the advanced functional modules of HCC

A novel method was utilized to simulate the dynamics of proteome alterations during the cancerous progression of HCC, which included a hierarchical cluster tree, topological overlapping matrix, weighted gene correlation network analysis and PPI network analysis. The present study simplified networks by modulating the correlative proteins to functional modules, which were involved in similar biological processes. As a result, 121 modules were established with members ranging from 20 to 2 (Fig. 3A). The majority of modules were extremely interconnected through their core nodes, which were considered to be candidate hub genes (Fig. 3C). Modules were annotated using clusterProfiler with GO terms and KEGG pathways. The results revealed that numerous modules were markedly enriched in cell cycle-associated progression, including module 2, 4, 5, 37 and 54 (Fig. 3B). In addition, module 17 was involved in the epoxygenase P450 pathway, module 27 was involved in fatty acid beta-oxidation using acyl-CoA dehydrogenase and module 37 was involved in the long-chain fatty acid biosynthetic process (Fig. 3D). In summary, the progression of HCC may occur via the rebalanced regulation and extensive reprogramming of mutually connected functional modules.

Figure 3.

Expression profiling of proteome reveals co-expression clusters and functional modules in HCC. (A) Distribution of modules size. Modules were identified by the superimposition of proteins in HCC onto the PPI network. The numbers of members from each module are exhibited in the figure. (B) Five interconnected modules of cell cycle derived from HCC, presenting the protein names and representative functional terms. (C) Allocation of 121 modules. Each node represents the individual module and their interactions by the module size. Edges connect modules that share PPIs. Boxed modules are further enlarged in (B) and (D). (D) Modules of DNA replication, epoxygenase P450 pathway, fatty acid beta-oxidation using acyl-CoA dehydrogenase and long-chain fatty acid biosynthetic process derived from HCC. HCC, hepatocellular carcinoma; PPI, protein-protein interaction.

To validate the hub genes of HCC, 360 HCC clinical and expression data were downloaded from TCGA database. Six hub genes were notable in the survival analysis from 36 candidate genes, which were significantly associated with patient prognosis. The aforementioned hub genes included: CCNB2, CDC20, MAD2L1, MCM2, CENPF and BUB1B (Fig. 4).

Figure 4.

Survival curves of 424 patients with hepatocellular carcinoma. (A) CCNB2, cyclin B2; (B) CDC20, cell division cycle 20; (C) MAD2L1, mitotic arrest deficient 2 like 1; (D) MCM2 minichromosome maintenance complex component 2; (Ε) CENPF, centromere protein F; (F) BUB1B, BUB mitotic checkpoint serine/threonine kinase B. P-value, hazard ratio and 95% CI are presented. CI, confidence interval.

TF and miRNA-gene network construction

Based on the former analyses, the present study aimed to determine the TFs and miRNAs associated with the identified hub genes. The ENCODE Chip-seq database was used to identify the TFs of the six hub genes. The results revealed structural maintenance of chromosomes 3 and tripartite motif containing 28 targeting CDC20, histone deacetylase 1 targeting CCNB2, regulatory factor X associated ankyrin containing protein targeting MAD2L, ZFP64 zinc finger protein targeting MCM2 and cAMP responsive element binding protein 3 targeting BUB1B were significantly upregulated in HCC (Fig. 5A). Furthermore, NetworkAnalyst and three other miRNA databases identified four hub miRNAs in HCC. Among them, hsa-mir-215-5p and hsa-mir-192-5p interacted with MAD2L1, CDC20, CENPF and BUB1B. Additionally, hsa-mir-1-3p interacted with MAD2L1, MCM2 and CENPF. Hsa-let-7c-5p interacted with CDC20, MCM2 and CCNB2 (Fig. 5B).

Figure 5.

TF and miRNA-hub gene networks in hepatocellular carcinoma. (A) TF-hub gene network. (B) miRNA-mRNA hub gene network. miRNA which intersect with >3 target genes were defined as hub miRNA. TF, transcription factor; miRNA, microRNA.

Validation of hub genes via RT-qPCR and western blotting

To confirm the hub genes identified using the aforementioned analyses and to determine the variations of expression at the mRNA and protein level, RT-qPCR and western blotting were performed (Fig. 6). The results revealed that the six hub genes were significantly upregulated in HCC when compared with the control group (P<0.001), which was concordant with the results obtained from the bioinformatics analysis.

Figure 6.

Validation of six hub genes through experiments. (A) Validation of hub genes by reverse transcription-quantitative PCR. Boxplots indicate the medians and dispersions of 47 HCC and their adjacent tissue samples. P-values were determined using a Student' t-test. ***P<0.001 and ****P<0.0001 with comparisons shown by lines. (B) Western blotting detection of hub genes. Lysates from 3 pairs of HCC and adjacent tissues were subjected to western blotting with antibodies against CCNB2, CDC20, MAD2L1, MCM2, CENPF and BUB1B. β-actin was used as the reference gene. HCC, hepatocellular carcinoma; CCNB2, cyclin B2; CDC20, cell division cycle 20; MAD2L1, mitotic arrest deficient 2 like 1; MCM2 minichromosome maintenance complex component 2; CENPF, centromere protein F; BUB1B, BUB mitotic checkpoint serine/threonine kinase B.

Discussion

It is well known that HCC is a prevalent primary cancer of the liver that occurs in numerous countries (40). However, there is a lack of consensus about its therapeutic practice, which is primarily due to the insufficient elucidation of the molecular mechanism underlying HCC (41). For clinical treatment, the detection of HCC in its early stage is important and providing treatment based on the molecular mechanism of HCC would be ideal. The present study comprehensively analyzed the gene expression profiles of HCC to determine the hub genes that may serve an indispensable function in HCC gene therapy.

In the present study, three data sets from the GEO database, GSE14323, GSE14520 and GSE17967, were used due to their large sample size to satisfy this analysis, and which all use the Affymetrix HG-U133A series annotation set to avoid losing numerous genes when merging the data. In addition, approximately three-quarters of HCCs are attributed to chronic HBV and HCV infections (42). The present study considered this factor when selecting the data, and thus these three data sets contain patients with HCC infected with these two viruses. By performing routine analysis, the present study identified 657 DEGs in HCC tissues when compared with their adjacent non-tumor tissues. Among them, 386 genes were upregulated and 271 genes were downregulated. The results of GO and KEGG analyses revealed that biological processes and pathways were enriched in cell cycle-associated processes, which is congruent to the results obtained in previous studies (9,43). The present study also demonstrated that the oxidation-reduction process was the most significant in HCC. A previous study determined that the reduction of fatty acid oxidation is associated with HCC progression (44). However, the traditional analysis of gene expression microarray data tends to identify individual genes that exhibit different expressions between two groups (40). This may be insufficient, as certain biological processes may be influenced by a slight change of individual genes or an entire gene network (43). GSEA represents a good tool to combat this problem. Therefore, to validate the results DEG analysis, the present study used GSEA to interpret the gene expression matrix at the whole transcriptome level. The results of GSEA were similar to those obtained via GO and KEGG DEG enrichment analysis. To further analyze co-expression and PPI networks, the WGCNA algorithm was superimposed onto the PPI STRING database. Cluster modules were screened out using the hierarchical cluster tree and topological overlapping matrix, and annotated via GO and KEGG. The resulting intricate network was predigested into modules, which were used to analyze core connections or hub genes. Survival analyses in 36 candidate hub genes were performed in the present study. A total of six genes were determined to be significantly associated with clinical prognosis. In addition, the transcription factors and miRNAs of six hub genes were predicted. Each are involved in the transcriptional and post-transcriptional regulation of gene expression. They may serve important regulatory functions in numerous biological processes, including differentiation, metabolism, development and cellular signaling (45). Thus, the identification of gene targets is important for the functional characterization of transcription factors and miRNAs and gives novel insights into the biological processes that may produce biomarkers and predictors of drug response for the disease.

The results of the present study demonstrated that CCNB2, CDC20, MAD2L1, MCM2, CENPF and BUB1B were important to HCC and may directly or indirectly regulate its progression. These results were partially congruent with a previous study, which indicated that CENPF and BUB1B were good predictors of HCC therapy (9). In addition, CCNB2 is key protein of the cyclin family and regulates the progression of the G2/M transition during the cell cycle (46). Furthermore, CCNB2, MAD2L1 and MCM2 are mitotic checkpoint associated proteins that are overexpressed in numerous types of human cancer (47–49). The dysregulation of mitotic checkpoints may also result in aneuploidy and the promotion of tumorigenesis (50). Next, the two genes CCNB2 and CDC20 are discussed. It was reported that CCNB2 was upregulated and was associated with degree of differentiation, tumor size, lymph node metastasis, distant metastasis and clinical stage, which represented a poor prognosis in patients with non-small cell lung carcinoma by survival analysis (51). The invasion and metastasis of bladder cancer was also inhibited by decreasing CCNB2, which prolonged the survival time of mice (52,53). On the contrary, CCNB2 overexpression promoted the cell proliferation and tumor growth of gastric cancer (46). Shubbar et al (54) reported that CCNB2 overexpression is an independent prognostic marker for breast cancer disease-specific survival time, as the c-index of CCNB2 alone is 0.662 and the prediction accuracy is improved with the passage of time. In the BEL-7404 HCC cell line, the downregulation of CCNB2 promotes apoptosis and may explain why the overexpression of CCNB2 results in the malignant potential of HCC (55). Another previous study also confirmed that CCNB2 knockdown inhibits tumor metastasis and prolongs the survival time of tumor-bearing nude mice (52). Based on these results, it was concluded that CCNB2 may be a key oncogenic target and is associated with HCC prognosis. CDC20 is one of the anaphase promoting complex (APC) activators and performs its functions via the ubiquitination and degradation of downstream substrates (56). Mounting evidence has determined that CDC20 is an accelerator of tumorigenesis and is overexpressed in numerous types of human cancer (57,58). For example, CDC20 expression is higher in pancreatic carcinoma compared with normal pancreatic tissue or chronic pancreatitis tissue (59). The depletion of CDC20 has also been demonstrated to contain cell growth and promote G2/M arrest (60,61). The expression of CDC20 has also been positively correlated with Tumor-Node-Metastasis stage and HCC differentiation (61). Considering the crucial function of CDC20 in tumorigenesis, an inhibitor of CDC20 may afford a medicinal window in a number of different human malignancies. To this end, the discovery and development of small molecule inhibitors that specifically target CDC20 oncoproteins may be a novel strategy for the treatment of a variety of human cancer types. For example, Zeng et al (62) proved that a small molecule, known as tosyl-L-arginine methyl ester, may weaken the interaction between APC and CDC20 and subsequent inhibit APC E3 ligase activity. Withaferin A is extracted from Withania somnifera, which has been identified to possess anti-tumor properties (56). It was reported that Withaferin A may gave rise to G2/M phase arrest and apoptosis in colorectal cancer (63). Furthermore, withaferin A may result in mitotic delay by degrading CDC20 and MAD2L1, indicating that inhibiting CDC20 may be a mechanism underlying the anti-cancer character of withaferin A (63). Other small molecules, including N-alkylated amino acid-derived sulfonamide hydroxamate (64), Ganodermanontriol (65), CFM-4 and BCHHD (66) were also proven to be anti-tumor drugs by targeting CDC20. Given the evident function of CDC20 in tumorigenesis, CDC20 may serve as a biomarker or drug target of HCC gene therapy.

The present study reported the biomarkers of HCC, which serve important functions in cancer detection and treatment. Cancer biomarkers are designed from tumor tissues, serum, DNA, enzymes, transcription factors and other proteins that may be measured, estimated and utilized as indicators of important biological process, pathways or pharmacological responses (67). Altogether, the integrated analysis of microarray data revealed six hub genes involved in cell cycle associated processes, which may be good indicators for HCC detection or therapy. Despite the oxidation-reduction process being notably involved in HCC, the present study failed to screen the hub genes as biomarkers for clinical prognosis. Future confirmatory experiments are therefore required to validate the results of the present study.