Comprehensive Identification of Bridge Genes to Explain the Progression from Chronic Hepatitis B Virus Infection to Hepatocellular Carcinoma

Introduction

Hepatocellular carcinoma (HCC) is the second leading cause of cancer-related death worldwide.^1–3 Surgical treatments have been developed to complement standard treatments, such as hepatectomy, liver transplantation and palliative treatment. These surgical procedures have improved prognosis of patients with HCC, but the overall survival rate is still poor due to the high mortality rate.³ HCC onset typically occurs after prolonged liver disease and is associated with numerous risk factors, which include chronic infection with hepatitis B virus (HBV) or hepatitis C virus (HCV), alcohol abuse, and nonalcoholic fatty liver disease.⁴

A World Health Organization survey in 2015 showed that approximately 257 million people (3.5% of the global population) suffer from chronic HBV infection.⁵ Globally, HBV infection is a major driver of HCC development.⁶ This is particularly true in Asia, where HBV infection is endemic and 70–80% of Asians having chronic HBV infection suffers from HCC.⁷ In Africa, HBV-associated HCC remains a leading cause of death among young people.⁵ However, the molecular link between HBV and HCC is unknown.

The HBV genome can integrate into the host DNA, modifying host genes near the gene insertion point, which destabilizes the host genome and produces an oncogenic fusion protein.⁸ In addition to altering cell survival, chronic HBV infection changes the liver environment.⁹ HBV infection induces an immune response, leading to chronic liver inflammation and fibrosis. Not only does this deplete the local immune response, but it also creates hypertonic conditions in the liver that lead to shrinkage and dysfunction of hepatocytes. This microenvironment can promote hepatocyte mutations that trigger apoptosis and malignant transformation.^10–13 However, the networks of genes that serve as a bridge from HBV-mediated liver inflammation to HCC development remain unclear.

To explore these networks, we analyzed gene expression profiles from patients with chronic HBV infection, HCC, or both. Differentially expressed genes (DEGs) were assigned to different modules based on predicted protein-protein interactions, and modules were identified that likely contribute to HBV-mediated HCC. Potential key bridging factors were confirmed by pathway enrichment and correlation with patient prognosis. Furthermore, potential transcriptional and post-transcriptional regulators of suspected bridge genes were identified. The results of these analyses may provide insights into the pathogenesis of HCC and help identify prognostic markers for earlier diagnosis.

Materials and Methods

Data Resources

To systematically analyze the molecular linkages between chronic HBV infection and HCC, we collected expression profiles of liver mRNA from the Gene Expression Omnibus. GSE83148¹³ contains 122 HBV-positive hepatitis patients and 6 normal control liver tissue samples, while GSE121248¹⁴ includes 70 HBV-positive HCC and 37 HBV-positive paracancerous liver tissue samples. All data were assessed using the Affymetrix Human Genome U133 Plus 2.0 Array. Next, high throughput RNA sequencing (RNA-seq) data from 341 HCC and 50 healthy samples were obtained from the Cancer Genome Atlas (TCGA) database. These HCC samples came from stage I (n = 171), stage II (87) or stages III–IV (83). The quality of the data samples was evaluated using principal component analysis (PCA).¹⁵

Identification of Differentially Expressed Genes (DEGs)

Genes differentially expressed (DEGs) among patients with HBV or HCC or both were calculated using the limma package¹⁶ in R. First the data underwent background correction. Then the “normalizebetweenarrays” function filtered the data based on levels of control and low-expression probes. We defined statistically significant DEGs as those with p < 0.05 and |log (fold change [FC]) | > 0.5.

For the RNA-seq gene expression profile data from TCGA, the 0 genes were filtered out by 50%, and the difference analysis was performed using degseq2.¹⁷ Each HCC sample, regardless of stage, was analyzed for differences from the corresponding control samples. DEGs that met a statistical threshold of p < 0.05 and |log₂FC| > 0.5 were analyzed further.

Gene Set Enrichment Analysis

Gene set enrichment analysis (GSEA) software computed enrichment scores for DEGs.¹⁸ Briefly, genes were ranked by significance after a signal-to-noise test was performed. Overlap of genes across datasets increases the level of significance. The magnitude of either the increase or decrease in the statistical significance of each step corresponds to the degree of change in gene expression under each condition. Next, the statistical significance of the enrichment score was estimated using a phenotype-based displacement test procedure. Then, the enrichment score of each sample was normalized based on sample number size, resulting in a standardized enrichment score (NES). We also calculated the false discovery rate (FDR). The pathways were analyzed in detail if their NES was associated with p < 1.0 and nominal (NOM) p < 0.05.

Identification of Interactions Between HBV- and HCC-Related Proteins

All human protein-protein interactions were identified based on the STRING database,¹⁹ which indicated 405 916 interaction pairs involving 10,514 proteins. The union of inflammation and cancer DEGs was mapped onto the protein-protein interaction network and maximal centrality clique (MCC). The MCODE algorithm²⁰ identified functional modules related to inflammation and cancer by applying default parameters and screening interaction pairs with a score >900. Cytoscape software²¹ visualized the interactions between HBV- and HCC-related proteins. The ClusterONE algorithm²² selected modules with degree > 30, and gene sets with p < 0.05 were considered important genes within the modules.

Analysis of Protein Crosstalk Between Chronic HBV Infection and HCC

We constructed a comprehensive protein pathway crosstalk network with 1000 random protein-protein interaction networks, where network size and the degree of each node remained constant. Significance of crosstalk for each pair of modules shared between HBV and HCC was calculated by p = n/N. The number of inter-module interaction pairs in N random networks (N = 1000) is estimated to be higher than the number of inter-module interaction pairs in the real system (denoted as n). When p < 0.01, then crosstalk module pairs are more significant than random. We found that the number of interaction pairs between modules was greater than the number of interaction pairs in a random background. Cytoscape software created complex relationship maps illustrating the crosstalk between the significant module genes.

Signaling Pathway Enrichment Analysis

The Clusterprofiler package²³ in R used protein-protein interaction modules to enrich for pathways of the Kyoto Encyclopedia of Genes and Genomes (KEGG) (p cut-off = 0.01, q cutoff = 0.01). Subsequently, we extracted HBV- and HCC-associated signaling pathways and considered these to be molecular bridges between the two diseases at the functional and pathway levels. To further evaluate the potential diagnostic value of bridge genes for HBV-related HCC, receiver operating characteristic curve (ROC) analysis was performed using the pROC routine²⁴ in R.

Pivot Analysis Predicts Transcriptional Regulators of Module Genes

Pivot nodes, which are nodes that interact with at least two modules, may be transcription factors or non-coding RNAs that regulate the expression of genes in modules related to HBV or HCC. The Python program identified pivot nodes, and significance (p < 0.01) was calculated using the hypergeometric test. This analysis allowed the prediction of regulators of bridge genes involved in HBV-associated HCC.

Molecular Docking Analysis

Docking is the process of bringing one molecule close to another. The PDB files of DOCK2 and its targets were downloaded from the Protein Data Bank (www.rcsb.org).²⁵ The docking of DOCK2 and its targets was performed using hex 8.0.0.0 software,²⁶ and the results were visualized using PyMOL software.²⁷

Results

Clustering of Genes Involved in Viral Hepatitis Type B and HCC

We analyzed the quality of HBV, HCC, and HBV-associated HCC gene expression profiles and signaling pathways corresponding to groups (modules) of genes (Figure S1 and 2). HBV-associated genes were significantly involved in the ErbB and Wnt signaling pathways (Figure S2A). The different stages of HCC (I–IV) were significantly associated with genes involved in the following pathways: stage I, genes associated with the cell cycle; stage II, RNA polymerase; stage III, pyrimidine metabolism; and stage IV, Notch signaling (Figure S2B, Table S1). The HBV-associated HCC group contained genes significantly involved in ubiquitin-mediated proteolysis and basal transcription factors (Figure S2C).

Examination of significant DEGs identified 1555 in the HBV group (Figure 1A). The HBV-associated HCC group included 2386 significant DEGs (Figure 1B). Within the HCC group, stage I contained 9122 significant DEGS; stage II, 10,249; stage III, 10,347; and stage IV, 8723 (Figure 1C). There were 5823 genes shared among the 4 stages, and 901 of these genes were progressive up- or down-regulation moving from stage I to stage IV (Figure 1D). We found 65 genes shared among the three groups of patients with HBV, HCC or both, and these genes showed different expression levels across the three groups (Figure 1E). We hypothesized that these genes act as “bridges” to progress from chronic HBV infection to HCC.

Figure 1 Genes differentially expressed in patients with HBV, HCC or both. Levels of mRNA were analyzed to identify differentially expressed genes, and results are shown as volcano plots comparing (A) HBV and healthy samples, (B) HBV-related HCC and healthy samples, or (C) stage-specific HCC and healthy samples. Red represents up-regulated genes, and blue represents down-regulated genes. (D) Venn diagram showing 5823 differentially expressed genes among stages I–IV of HCC. (E) Venn diagram showing 901 genes differentially regulated between disease and healthy states, and the subset of 65 genes shared among patients with HBV or HCC or both. (F) Differentially expressed genes clustered into 10 modules (M) based on protein-protein interactions. The middle cluster contains genes overlapping among all three groups (HBV, HCC, HBV-associated HCC). Red nodes contain differentially expressed genes present in all three groups; orange nodes, genes in two groups; blue nodes, genes only in HBV; yellow nodes, genes only in HBV-associated HCC; purple nodes, genes only in HCC.

Interactions Between Proteins Encoded by Overlapping DEGs

We searched for genes that were differentially expressed in HBV or HBV-associated HCC or that were persistently dysregulated in HCC. These genes were then screened to determine whether the encoded proteins were potentially involved in progression from chronic HBV infection to HCC. A protein-protein interaction network based on 4042 genes showed that 352 genes clustered into 10 modules associated with dysfunction (Figure 1F). The DEGs associated with HBV and HCC clustered tightly within the same module, and all modules contained genes associated with HBC, HCC or HBV-associated HCC. This suggests that gene modules mediate a link between chronic HBV infection and HCC.

Significant Crosstalk and Shared Pathways Between Liver Inflammation and HCC

We examined connections among the 10 modules and found 7 with significant crosstalk (Figure 2A). The module M1 had the highest connectivity, suggesting that its genes may reveal critical insights into how liver inflammation may lead to HCC. KEGG pathway analysis revealed that two gene modules were involved in HCC (M7 and M5), progesterone-mediated oocyte maturation (M1 and M4), p53 signaling pathway (M1 and M3), Staphylococcus aureus infection (M3 and M6), and RNA transport (M1 and M2). We found that most of the signaling pathways associated with the module genes were related to microenvironment and infection immunity (Figure 2B), including glycolysis, TGF-beta signaling, and Kaposi sarcoma-associated herpesvirus infection.

Figure 2 Protein crosstalk and signaling pathways associated with genes bridging chronic HBV infection and HCC. (A) Protein-protein interactions were assessed by measuring crosstalk between modules shared by HBV and HCC. Red nodes, module M1; blue nodes, M3; pink nodes, M5; green nodes, M6; yellow nodes, M7; orange nodes, M9; and purple nodes, M10. (B) KEGG analysis identified signalling pathways involving module genes from all groups. The most common pathways are related to cancer (orange arrow), infection (purple triangle), metabolism (green square), and immune system (red diamond). Yellow circles are labeled with their corresponding module number.

Genes Bridging Chronic HBV Infection and HCC

We identified 65 DEGs whose expression was higher in HCC and HBV-associated HCC than in HBV (Figure 1E). Of these, we identified 32 showing sustained up- or down-regulation suggesting that they are critical for HBV-associated HCC (Figure 3A).

Figure 3 Identification of bridge genes for HBV-induced HCC. (A) Heatmap depicting 32 dysregulated genes in the progression from HBV to HCC. Red represents up-regulated and blue represents down-regulated. (B) Line graph depicting survival probabilities over time for HCC patients showing high or low expression of the respective gene. (C) Receiver operating characteristic curves showing prognostic potential of bridge genes in HBV-related HCC. (D) Heat map depicting correlations of bridge genes. Red stands for positive correlation and blue stands for negative correlation.

Next, we analyzed network connectivity of the genes clustered into modules. Genes with higher connectivity are thought to have greater regulatory potential in the network. CDC20 had the maximum network connectivity. Other genes, including CCNB1, NUP133, KIF2C, and NUP107, also showed high connectivity. Moreover, we found that the drug Lomustine may target the protein product of the pivot gene CCNB1, based on the DRUGBANK database (https://go.drugbank.com/).

Genes with high network connectivity were considered to be potentially associated with patient prognosis. CPEB3, RAB26, SLCO1B1, ST3GAL6, XK, CDC20 and NUP107 were significantly linked to low survival rates in patients with HCC (Figure 3B). Cases with high expression of RAB26 showed better survival. The prognostic value of bridge genes in HBV-associated HCC was tested directly using ROC curves. CDC20 (area under ROC = 0.919) was a superior biomarker than other bridge genes (Figure 3C). Further analyses showed that bridge genes were significantly associated with phenotypic markers present in HBV and HCC liver tissue, including markers corresponding to different HCC stages (Figure 3D).

Non-Coding RNA and Transcriptional Factors Regulate Bridge Genes in HBV-Associated HCC

To explore the potential transcriptional regulation of bridge genes, we applied prediction analysis to the dysfunction modules based on transcriptional and post-transcriptional regulatory relationships. We found that a total of 202 microRNAs (miRNAs), 7 long non-coding RNAs (lncRNAs) and 44 transcription factors were associated with our modules. These equated to 2086 miRNA-module interaction pairs, 169 lncRNA-module interaction pairs, and 260 transcription factor-module interaction pairs.

In-depth analysis of these predicted interactions revealed that significant regulatory effects were exerted on the bridge gene CCNB1 by 28 miRNAs, 2 lncRNAs and 6 transcription factors (Figure 4). Three miRNAs and one lncRNA substantially regulated NUP107. NUP133 was significantly impacted by 10 miRNAs and 1 lncRNA. These transcription regulators may mediate the occurrence and development of HBV-associated HCC. Furthermore, to verify if there were targeted bindings among molecules shown in Figure 5, molecular docking was performed. Negative binding energy, which suggests possible binding, was identified for CCNB1 and YBX1 (energy = −700.656), CCNB1 and NFKB1 (energy = −727.3), CCNB1 and E2F4 (energy = −572.42), CCNB1 and 2EBT (energy = −641.16), CCNB1 and KLF5 (energy = −754.62), as well as CCNB1 and RELA (energy = −366.21) (Figure 5).

Figure 4 Pivot regulators of bridge genes within the module. Regulators interacting with multiple modules at least twice were identified as potential regulators of bridge genes. Green squares represent bridge genes in modules; blue hexagons, microRNAs; yellow triangles, long non-coding RNAs; pink arrows, transcription factors; and red circles, signaling pathways involved in bridge genes.

Figure 5 Molecular docking between pivot regulators and the bridge genes in the module. Docking of CCNB1 with (A) YBX1, (B) NFKB1, (C) E2F4, (D) KLF5, (E) TP53, or (F) RELA.

Discussion

Approximately half of all HCC cases are associated with HBV, making it the most common cause of cancer worldwide.²⁸ Therefore, the present study integrated resources from multiple databases to determine which genes may be responsible for the transition from HBV-induced liver inflammation to HCC onset. We assessed changes in gene transcription levels in healthy control patients and patients with HBV, HCC or both. Protein-protein interaction networks, crosstalk analysis, transcriptional regulation, and post-transcriptional regulation were combined to identify genes bridging chronic HBV infection and HBV-associated HCC. These findings provide new insights into the gene and protein interactions common to these diseases. The predicted upstream non-coding RNA and transcription factors provide numerous leads and hypotheses for subsequent research.

Our analyses found unique and overlapping DEGs in liver samples from patients with HBV, HCC or both. Then the corresponding proteins were examined for enrichment in molecular pathways and involvement in protein-protein interactions. Based on this, the DEGs were clustered into modules. Some genes were identified as key players in bridging HBV inflammation to HCC, including CDC20, CCNB1, and KIF. These genes have well established roles in HCC development. For example, CDC20 expression is elevated in HCC tissues compared to paracancerous tissue, and we and others showed that its expression is altered in HBV-positive HCC.^29–31 Expression of both CCNB1 and KIF2 is associated with HCC stage.^32,33 These genes and their encoded proteins are linked to poor prognosis and may be meaningful biomarkers for HCC diagnosis. In addition, molecular docking identified several potential binding partners of CCNB1: YBX1, NFKB1, E2F4, 2EBT, KLF5, and RELA.

We narrowed down the list of bridge factors by examining protein crosstalk, ie proteins shared or interacting across molecular pathways. We identified four pathways involving significant crosstalk: cancer-related signaling, inflammation-related signaling, immune-related signaling, and microenvironment-related signaling. N-acetyltransferase10 (NAT10), ubiquitin-binding enzyme E2T (UBE2T) and Bcl2-like protein 12 (Bcl2L12) were included in these pathways, and these genes have been implicated in HCC.^34,35 The protein p53 was also identified, reflecting its well-established multi-faceted roles in cancer development.³⁴ RNA transport plays a key role in the process of HBV-induced HCC and leads to further deterioration in HCC disease.³⁶ Kinesin superfamily proteins (KIFs) largely serve as molecular motors on the microtubule system and transport various cellular proteins, potentially linking HBV infection to HCC.³⁷ Identifying key genes and proteins responsible for bridging HBV to HCC development indicates how these diseases are related, which may improve our understanding of the pathogenesis of HBV-related HCC.³⁸

Beyond understanding HBV-associated HCC pathogenesis, our analyses may also help identify prognostic markers to improve therapy and survival. We found 65 genes that were dysregulated among patients with HBV, HCC or HBV-related HCC, among which 32 were dysregulated during HBV-induced HCC onset. Five of these genes with large areas under the ROC curve have significant effects on the prognosis of HCC: CPEB3, RAB26, SLCO1B1, ST3GAL6, and XK. We hypothesize that these genes play an important role in the development of HBV-induced HCC. We found that CPEB3 can be targeted by miR-107, and the CPEB3-EGFR axis can be targeted by miR-452-3p. Both interactions promote HCC tumor growth,^38,39 suggesting that CPEB3 is a novel tumor regulatory factor. RAB26 mediates targeted autophagy to play an important regulatory role in the tumor microenvironment.⁴⁰ SLCO1B1 affects drug absorption, distribution, metabolism, and excretion in patients with HCC.⁴¹ ST3GAL6 appears to influence HCC cell growth, migration, and invasion.⁴² Finally, the Xkr 8 protein mediates the display of phosphatidylserine on the cell surface in response to apoptotic stimuli, which represents an “eat me” signal for phagocytes, giving Xkr 8 a significant regulatory role in the development of HCC.⁴³ In summary, these genes affect prognosis of patients with HBV-related HCC via interactions with other genes.

Transcriptional and post-transcriptional regulation has been recognized as key factors in the development and progression of many cancers. To elucidate the transcriptional regulatory factors associated with the molecular links between the two diseases, we performed a pivotal regulator analysis. We found that the following regulators had significant regulatory effects on bridge factors in liver tissue of HCC patients: lncRNAs MALAT1, TUG1, and CRNDE; miRNAs miR-590-3p, miR-340-5p, miR-139-5p, and miR-93-5p; and transcription factors E2F4, NFKB1, and RELA. In addition, molecular docking identified YBX1, NFKB1, E2F4, 2EBT, KLF5, and RELA as potential binding partners of CCNB1, further supporting the pivot analysis. Interestingly, previous work has described a potential role for all these mediators in HBV-associated HCC development.^44–46 For example, lncRNA-MALAT1 targets CCNB1 to block cancer stem cell development induced by the hepatitis B virus protein HBx through the PI3K/Akt signaling pathway.^47,48 Hepatitis B virus promotes cancer cell migration by down-regulating miR-340-5p and up-regulating STAT3.⁴⁹ Polymorphism in the NFKB gene and its promoter is significantly associated with the occurrence of liver cancer in patients chronically infected with HBV of genotype C.^50,51 We propose that these factors are key regulators of the development of HBV-associated HCC. In-depth study of how these pivotal regulators modulate the expression of bridge genes will help us fully understand the pathogenesis of the disease.

There are some limitations in this study. Firstly, the data set used in this study is not very large and based on the same platform. In order to increase the robustness of the conclusions, we hope to verify the results of this study in a different batch and a larger data set in the future. Secondly, existing studies have shown that Hbx protein is related to the progress of HBV or HCC, but we have not studied it in this study. We hope to make up for this shortcoming in the next study. Thirdly, the results of this study are based on biological informatics analysis. Due to our lack of funding, we are currently unable to carry out further research. But this research provides a preliminary basis for us to carry out molecular experiments in the future. We also hope that molecular experiments can further confirm the results of this study in the future. Fourthly, the data that we collected did not report indicators about infection, so we considered viral infection to be absent. We hope that we can conduct the study based on a dataset including HCC patients without HBV or HCV infection.

Conclusions

We identified candidate bridge genes and proteins as well as their candidate regulators that may explain the progression from chronic HBV infection to HCC. These results provide testable hypotheses to guide research into HBV-associated HCC as well as efforts to identify diagnostic and prognostic markers.