The Tohoku Journal of Experimental Medicine
Online ISSN : 1349-3329
Print ISSN : 0040-8727
ISSN-L : 0040-8727
Regular Contribution
Identification and Analysis of Subtypes of Liver Cancer Based on Genes Related to E3 Ubiquitin Ligases and Deubiquitinating Enzymes
Yunming HuPeipei HuangFeizhao Jiang
Author information
JOURNAL OPEN ACCESS FULL-TEXT HTML
Supplementary material

2024 Volume 262 Issue 2 Pages 75-84

Details
Abstract

Recent studies have reported a correlation between ubiquitination or deubiquitination and cancer development. But mechanisms underlying the roles of genes associated with E3 ubiquitin ligases and deubiquitinating enzymes (DUB) in liver cancer remain to be explored. We analyzed and screened differentially expressed genes related to E3 ubiquitin ligases and DUB in liver cancer on the basis of public databases. Cluster analysis was utilized to classify liver cancer samples into different subtypes. Survival analysis, immune analysis, and pathway enrichment analysis were performed on the subtypes. We constructed a protein-protein interaction network using STRING to screen hub genes. Finally, we used the Connectivity Map (CMap) database to predict targeted small molecules. The results show that a total of 139 differentially expressed E3/DUB genes in liver cancer were screened. Then, liver cancer was classified into two subtypes, cluster 1 and cluster 2, based on E3-related and DUB-related genes. Patients in cluster 1 had higher survival rates and immune levels than those in cluster 2. Four hub genes (RPSA, RPS5, RPL30, and RPL8) significantly affecting the survival of the two subtypes of liver cancer patients were identified based on cluster 1 and cluster 2. Finally, the CMap database predicted that small-molecule drugs including probenecid, dexamethasone, and etomidate may improve the prognosis of liver cancer patients. These findings may offer a reference for risk stratification studies and drug development in liver cancer.

Introduction

Liver cancer is a common and lethal malignancy that poses a serious threat to human health around the world (Siegel et al. 2020). Global Cancer Statistics (GLOBOCAN) estimated 906,000 new cases of liver cancer and 830,000 deaths worldwide in 2020, with 632,000 new cases occurring in men and 273,000 in women, and 578,000 deaths occurring in men and 253,000 in women (Sung et al. 2021). With the advancement of technology and the accumulation of medical experience in recent years, the treatment of liver cancer has shifted from primarily surgical treatment to the comprehensive use of multiple methods (Anwanwan et al. 2020). However, due to the various risk factors associated with liver cancer, the 5-year survival rate of liver cancer patients is only 12.5%, and the prognosis of liver cancer patients remains poor (Luo et al. 2021). Therefore, studying the regulatory mechanisms underlying liver cancer development and discovering new drug targets are of crucial importance for the diagnosis and treatment of liver cancer.

Ubiquitination is a post-translational modification in which E3 ubiquitin ligase catalyzes the final step of ubiquitination, which is also a key step in ubiquitin conjugation, directly determining the substrate pathway mediated by ubiquitin (Buetow and Huang 2016; Rennie et al. 2020). Deubiquitinating enzymes (DUBs) can specifically hydrolyze ubiquitin molecules from ubiquitinated protein substrates to maintain free ubiquitin levels in cells (Hussain et al. 2010; Nakamura 2018). Previous studies have shown that abnormal expression of E3 ligases and DUBs is implicated in adverse outcomes in a variety of cancers. For example, Li et al. (2021) presented that high expression of DUB ubiquitin specific peptidase 39 (USP39) and E3 ligase tripartite motif containing 26 (TRIM26) in hepatocellular carcinoma (HCC) tissues drives cancer cell proliferation and migration via manipulation of ZEB1 protein level. In addition, several studies have shown that the ubiquitin system involving E3 and DUB is critical in tumor prognosis. Xu et al. (2021b) found associations between the abnormal expression of multiple ubiquitination-related genes and the prognosis of lung adenocarcinoma patients, which can predict the overall survival of patients. Wang et al. (2023) found that USP39, a DUB family member, can reduce the degradation of β-catenin through its deubiquitination function and drive cell proliferation and migration through repression of splicing of E3 ligase TRIM26 in cancer. Moreover, ubiquitin ligases and DUBs can regulate patients’ immune response by interacting with other factors in vivo, thereby affecting cancer progression (Damgaard et al. 2012; Fiil et al. 2013; Hrdinka et al. 2016). Therefore, understanding the ubiquitination status in liver cancer patients has important value in the diagnosis, monitoring, and exploration of tumor progression and changes in the immune microenvironment of cancer.

This study collected genes related to E3 and DUB in liver cancer from public databases. Subsequently, a series of bioinformatics analyses (differential analysis, subgroup analysis, survival analysis, immune analysis, mutation analysis, and drug sensitivity analysis) were systematically conducted to elucidate the mechanism of ubiquitination-related genes in liver cancer. Some potential drug targets for lung cancer were identified. This study used novel approaches to investigate the molecular mechanism of ubiquitination in liver cancer and targeted therapy.

Materials and Methods

Data collection

RNA sequencing transcriptomic data and clinical data related to liver cancer were downloaded from The Cancer Genome Atlas (TCGA) (https://portal.gdc.cancer.gov/). The dataset included 374 liver cancer samples and 50 normal samples. E3/DUB-related genes were identified from Integrated annotations for Ubiquitin and Ubiquitin-like Conjugation Database (IUUCD; https://iuucd.biocuckoo.org/) (Xu et al. 2021b).

Screening and functional analysis of E3/DUB-related differentially expressed genes (DEGs)

Based on the liver cancer transcriptomic data (normal and tumor groups), the R package edgeR (Robinson et al. 2010) was utilized for differential analysis [|log fold-change (FC) | > 1, False Discovery Rate (FDR) < 0.05] to identify DEGs in the tumor group. Subsequently, these DEGs were overlapped with E3/DUB-related genes to obtain E3/DUB-related DEGs in liver cancer. The protein-protein interaction (PPI) network, which reflected the intrinsic relationship between E3/DUB genes, was integrated and constructed using STRING (https://cn.string-db.org/cgi/input.pl). Then, following previous research methods (Bian et al. 2021), Chi-square tests were conducted on the copy number variations (CNV) of E3/DUB-related DEGs in liver cancer based on CNV data from TCGA. The R package RCircos (Zhang et al. 2013) was used to plot the specific location of these genes on chromosomes. The R package clusterprofiler (Yu et al. 2012) was applied for Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses of these genes.

Identification and immune analysis of E3/DUB gene-based subtypes

Based on expression data of E3/DUB DEGs, we performed cluster analysis on tumor samples using the R package ConsensusClusterPlus (Hu et al. 2022). Following that, we employed R package survival for survival analysis on samples from each cluster. The ESTIMATE method was utilized to compute the ESTIMATE score, immune score, stromal score, and tumor purity of the samples from each cluster. Additionally, according to gene expression levels in 29 immune-related gene sets, we performed single-sample gene set enrichment analysis (ssGSEA) using the R package GSEAbase (Xu et al. 2021a). We also calculated expression levels of immune checkpoint genes in the cluster samples and plotted box plots.

Gene set enrichment analysis (GSEA)

Following previous research methods (Yao et al. 2021), we used GSEA software (http://software.broadinstitute.org/gsea/index.jsp) and the c2.cp.kegg. V 2022.1. Hs.symbols.gmt gene set database to perform GSEA on the mRNAs related to different groups using Java under 1,000 random sample permutations.

Mutation analysis

The corresponding somatic mutation information for liver cancer was downloaded from TCGA dataset. Tumor mutation burden (TMB) was defined as the number of somatic, coding, base substitution, and insertion/deletion mutations per million bases in the genome detected using nonsynonymous and frame-shifting insertions and deletions (indels) at a 5% detection threshold. We used the R package maftools (Xu et al. 2021a) to compute the number of somatic nonsynonymous point mutations in each sample. Subsequently, we performed Wilcoxon tests on the TMB values of different subtype groups and plotted violin plots. R package GenVisR (Skidmore et al. 2016) was applied to plot waterfall plots of the top 30 mutated genes in different subtype groups to present the mutation landscape.

Hub gene identification

To explore the differences between gene subgroups, we employed the edgeR package to conduct differential analysis on data from cluster 1 and cluster 2 (|logFC| > 1.0, FDR < 0.05). With cluster 1 as the control group, DEGs in cluster 2 were screened out. We performed interaction analysis on DEGs in cluster 2 on the STRING website and generated a PPI network. Subsequently, we used the CytoHubba plugin in Cytoscape to screen hub genes in the PPI network (using the MCC algorithm, node = 10). Based on hub gene levels, we used timeROC package (Blanche et al. 2013) to draw corresponding receiver operation characteristic (ROC) curves to predict whether patients have the disease based on gene expression. Furthermore, we tested gene expression with an area under the curve (AUC) value greater than 0.7 in both tumor samples and normal samples, and plotted box plots to further verify the predictive results of the ROC.

Drug sensitivity analysis

To identify new potential targets and more effective therapeutic drugs, we uploaded the top 150 upregulated genes from cluster 2 to the CMap database (https://clue.io/) to identify putative small molecule drugs for the treatment of liver cancer patients. Drugs with negative scores indicate that they can reverse the upregulated expression of the input genes (Zhu et al. 2022).

Results

E3/DUB-related DEGs in liver cancer

By performing differential analysis on data from normal and tumor groups of liver cancer in TCGA, 4,932 DEGs were identified (3,892 upregulated and 1,040 downregulated) (Fig. 1A). Taking the intersection of DEGs and E3/DUB-related genes (1,024 genes), the upset plot illustrated that 139 E3/DUB DEGs were acquired in liver cancer (Fig. 1B). The PPI network diagram presented 39 nodes and 76 interactions, indicating complex interactions between these E3/DUB DEGs (Fig. 1C). To explore the functional features and biological effects of E3/DUB DEGs, we performed GO and KEGG enrichment analyses. The GO enrichment results showed substantial enrichment of E3/DUB DEGs in biological functions such as protein autoubiquitination, protein deubiquitination, and positive regulation of cellular catabolic processes (Fig. 1D). The KEGG enrichment analysis revealed the enrichment of genes in ubiquitin-mediated proteolysis, JAK-STAT signaling, TNF signaling, and Notch signaling pathways (Fig. 1E). The chromosome distribution map revealed that the CNVs of E3/DUB DEGs were distributed in multiple locations of the chromosome (except for chromosomes 11 and 13), and DEGs had copy number gains or losses (Fig. 1F). Taken together, these analyses suggested that these E3/DUB DEGs were implicated in liver cancer progression.

Fig. 1.

Screening and functional analysis of E3/DUB-related differentially expressed genes (DEGs).

(A) Volcano plot of DEGs in liver cancer. (B) Upset plot of DEGs and E3/DUB genes. (C) The protein-protein interaction (PPI) network of E3/DUB DEGs. (D) Gene Ontology (GO) enrichment analysis of E3/DUB DEGs. (E) Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis of E3/DUB DEGs. (F) Chromosomal copy number variations (CNV) distribution of E3/DUB DEGs, where the outer circle represents the 24 chromosomes and the inner circle shows the distribution of CNV.

Identification of tumor subtypes based on E3/DUB ubiquitination-related genes

Following cluster results of the differential expression matrix of E3/DUB-related genes, samples were divided into cluster 1 and cluster 2 (Fig. 2A). Kaplan-Meier (K-M) survival analysis was subsequently done on samples from cluster 1 and cluster 2. As shown in Fig. 2B, cluster 1 had a higher survival rate than cluster 2.

Fig. 2.

Tumor subtype identification based on E3/DUB-related genes.

(A) Gene typing cluster analysis of liver cancer samples. (B) Survival analysis of E3/DUB-related subtypes.

Immune differences between E3/DUB ubiquitination subtypes

The ESTIMATE algorithm can estimate stromal and immune cells in malignant tumors, inferring tumor purity from the unique properties of transcriptional profiles. ssGSEA heatmap depicted that immune level of cluster 2 was lower than that of cluster 1 (Fig. 3A). The violin plot showed that immune, stromal, and ESTIMATE scores of cluster 2 were noticeably lower, and the tumor purity score was significantly higher as compared with cluster 1 (Fig. 3B). Immune checkpoint box plot presented that the expression levels of most immune checkpoints (TNFSF13B, CD160, CD200, CD200R1, CD27, CD274, CD28, CD40, CD40LG, CD80, IDO2, PDCD1LG2, TMIGD2, TNFRSF25, TNFRSF8, TNFSF18) in cluster 2 were significantly lower than in cluster 1, but ADORA2A and TNFRSF18 levels in cluster 2 were notably higher than in cluster1 (Fig. 3C). Taken together, these analyses suggested that cluster 1 had a higher immune level than cluster 2. The higher immune, stromal, and ESTIMATE scores, as well as immune checkpoint expression levels may be related to a higher survival rate in the cluster 1 subtype.

Fig. 3.

Immune differences between E3/DUB ubiquitination subtypes.

(A) Single-sample gene set enrichment analysis (ssGSEA) heatmap. (B) Violin plot of the immune score, stromal score, ESTIMATE score, and tumor purity score. (C) Box plot of immune checkpoints. *P < 0.05, **P < 0.01, ***P < 0.001.

Pathway enrichment differences between different subtypes

GSEA unveiled the main enrichment of cluster 2 in pathways such as Pyrimidine Metabolism, Oxidative Phosphorylation, Glutathione Metabolism, and Other Glycan Degradation (Fig. 4). As a result, we speculated that differences in prognosis between the two groups may be regulated by life activities such as pyrimidine metabolism, oxidative phosphorylation, glutathione metabolism, and other glycan degradation.

Fig. 4.

Gene set enrichment analysis (GSEA) results.

Mutational differences in different subtypes

To explore the tumor mutation situation in liver cancer patients, we analyzed the TMB of liver cancer samples grouped by different subtypes. TMB of cluster 2 was significantly higher than that of cluster 1 (Fig. 5A). Subsequently, we analyzed the mutation frequency of the top 30 genes with single nucleotide variants (SNVs) in the dataset and found that the mutation frequency of some genes varied between the two groups. TP53 had the highest mutation frequency in cluster 1, while CTNNB1 had the highest mutation frequency in cluster 2 (Fig. 5B, C).

Fig. 5.

Tumor mutation analysis.

(A) Tumor mutation burden (TMB). (B) TMB waterfall plot of the top 30 mutated genes in cluster 1. (C) TMB waterfall plot of the top 30 mutated genes in cluster 2.

Differences in key hub genes between different subtypes

To explore the differences between gene subgroups, we conducted differential analysis on the data from cluster 1 and cluster 2. The results showed that cluster 2 contained 1,829 DEGs, including 609 upregulated genes and 1,220 downregulated genes (Supplementary Table S1). The PPI network showed complex interactions between DEGs. Using an interaction score threshold of > 0.9, we obtained 553 nodes and 1,798 associated interaction connections (Fig. 6A). Furthermore, we counted the top 10 hub genes in the PPI network ( RPSA, RPS19, RPL32, RPS5,RPL30, RPS20, RPL14, RPL13, RPL8, andRPL7), and drew a new PPI network (Fig. 6B). Further analysis showed that the AUC values of these 10 hub genes were all greater than 0.519, and that the AUC values of RPSA,RPS5, RPL30, and RPL8 were all greater than 0.7 (Fig. 6C), indicating that these genes could be good predictors of liver cancer. Similarly, the box plot revealed that levels of these four hub genes with AUC values greater than 0.7 were notably higher in tumor samples than in normal samples (Fig. 6D). Based on these analyses, we speculated that the high expression of four hub genes in tumor samples was implicated in the development of liver cancer.

Fig. 6.

Subtype grouping differential analysis and gene interactions.

(A) The protein-protein interaction (PPI) network of differentially expressed genes (DEGs) between different subtypes. (B) PPI network of hub genes. (C) Receiver operating characteristic (ROC) curve of hub genes. (D) Box plots showing the expression levels of the 4 hub genes with area under the curve (AUC) values greater than 0.7 in tumor and normal samples.

Drug sensitivity analysis

To explore potential therapeutic targets and more effective drugs for liver cancer, we used the CMap database to predict the top 150 upregulated genes in cluster 2. The results showed that there were six targeted small-molecule drugs with potential therapeutic effects on liver cancer, namely probenecid, sitagliptin, sibutramine, dopamine, dexamethasone, and etomidate (Table 1). Moreover, the scores of these six candidate small-molecule drugs were negative, indicating that these drugs may act on multiple target genes in liver cancer. Combined with previous studies (Zhu et al. 2022), drugs with negative scores can reverse their upregulation of genes. Therefore, we speculated that these six candidate small-molecule drugs could improve survival in liver cancer patients.

Table 1.

CMap prediction results.

Discussion

E3 ubiquitin ligases and DUBs are pivotal participants of the ubiquitin-proteasome system, where they play critical roles in cellular protein homeostasis (Satija et al. 2013). Abnormal expression of E3 ligases and DUBs has been implicated in poor outcomes of various cancers, including glioblastoma (Bellail et al. 2012), breast cancer (Lee et al. 2017), and HCC (Li et al. 2021). Feng et al. (2020) reported that the E3 ubiquitin ligase A20 can suppress HCC cell proliferation and migration by downregulating glycolytic enzyme phosphofructokinase-1 (PFK1). Chen et al. (2023) demonstrated that overexpression of DUB ubiquitin specific peptidase 15 (USP15) facilitates growth and migration and is implicated in unfavorable prognosis of non-small cell lung cancer patients. However, research on ubiquitination- and deubiquitination-related genes is still limited. Therefore, this study stratified liver cancer patients based on E3-related and DUB-related genes and identified effective biomarkers for prognosis prediction and treatment.

In this study, a comprehensive analysis of E3-related and DUB-related genes in liver cancer was conducted using public databases, and 139 E3/DUB DEGs were identified, which exhibited complex interactions among each other. Enrichment analysis revealed the main enrichment of these genes in life activities like protein ubiquitination, deubiquitination, and cellular catabolism. In addition, over the past few decades, there has been increasing research on the tumor microenvironment. Immune cells and stromal cells are two major types of tumor components that are of paramount significance in the diagnostic and prognostic assessment of tumors. Therefore, as early as 2013, Yoshihara et al. (2013) developed an algorithm (ESTIMATE) to estimate the immune cell components in tumor tissues. In this study, it was found that cluster 1 had a higher survival rate, immune score, and stromal score than cluster 2. As a result, it was speculated that patients with higher immune components had a better prognosis, which was congruous with previous research (Lei et al. 2021).

To further analyze the reasons for the differences in survival among gene subgroups, we conducted differential expression analysis on the genes of each subgroup and constructed a PPI network. From this network, we identified four hub genes (RPSA,RPS5, RPL30, and RPL8) with AUC values greater than 0.7 that were significantly overexpressed in tumors. Aberrant expression of these genes is implicated in tumor occurrence and development. Membrane receptor ribosomal protein RPSA can promote cancer progression by stimulating cancer cell migration and invasion (Lefebvre et al. 2020). Wu et al. (2019) found that high expression of RPSA can affect the migration and invasion of pancreatic cancer cell lines via manipulation of PI3K and MAPK signaling pathways, and patients with high RPSA expression often have unfavorable prognoses. However, inhibiting RPSA significantly weakens invasion and migration and leads to a decrease of p-AKT and p-ERK1/2. RPS5 is also a ribosomal protein gene and a haploinsufficient tumor repressor for various cancer types (Fancello et al. 2017). Pan et al. (2022) identified RPS5 as a key gene in circulating tumor cells of colorectal cancer through bioinformatics analysis. RPL30 and RPL8 are also important ribosomal proteins. Previous studies have found that the ferroptosis marker RPL8 is overexpressed in multiple cancers. For example, Fan et al. (2023) found through comprehensive multi-omics analysis that higher RPL8 expression is associated with shorter overall survival in HCC patients and that RPL8 may be a fresh therapeutic target for the treatment and prognosis evaluation. Hence, the four hub genes identified in this study can affect the prognosis of multiple cancers and have the potential to become new biological markers that influence the prognosis of liver cancer.

To identify new potential targets and more effective therapeutic drugs, we screened the CMap database and identified six potential small molecule drugs: probenecid, sitagliptin, sibutramine, dopamine, dexamethasone, and etomidate. Previous studies have shown that probenecid can sensitize neuroblastoma cancer stem cells to cisplatin, suggesting that combining probenecid with cisplatin could be a new cancer treatment strategy (Campos-Arroyo et al. 2016). Additionally, Salah et al. (2021) found that sitagliptin exerts its anti-tumor activity by inducing apoptosis and inhibiting oxidative stress and inflammation. Similarly, dexamethasone and other glucocorticoids have been widely used in the treatment of cancer. For example, Motafeghi et al. (2023) found that dexamethasone alone or in combination with etoposide can induce oxidative stress and exert anti-cancer effects by disrupting DNA. In addition, coenzyme Q10 exerts antioxidant effects on dexamethasone-induced liver toxicity through the improvement of mitochondrial function and reduction of cysteine dioxygenase-3 activity. Based on these findings, we speculated that these small molecule drugs may be potential anti-cancer drugs for liver cancer treatment.

In summary, this study analyzed E3-related and DUB-related genes in liver cancer and used their expression to cluster liver cancer patients into subtypes. We then performed survival analysis, immune analysis, differential analysis, enrichment analysis, and drug sensitivity prediction on different subtypes. The results not only predict the prognosis and immune status of liver cancer patients, but also have certain guiding significance for the management of liver cancer. However, this study still has some limitations, particularly the mechanisms of the small molecule targeted drugs predicted in this study in liver cancer, which requires further investigation.

Author Contributions

Yuming Hu conceived and designed the study. Yuming Hu and Peipei Huang performed the experiments. Yuming Hu and Feizhao Jiang wrote the paper. Peipei Huang reviewed and edited the manuscript. All authors read and approved the manuscript.

Conflict of Interest

The authors declare no conflict of interest.

References
 
© 2024 Tohoku University Medical Press

This article is licensed under a Creative Commons [Attribution-NonCommercial-NoDerivatives 4.0 International] license.
This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License (CC-BY-NC-ND 4.0). Anyone may download, reuse, copy, reprint, or distribute the article without modifications or adaptations for non-profit purposes if they cite the original authors and source properly.
https://creativecommons.org/licenses/by-nc-nd/4.0/
feedback
Top