Integrated Multi-Omics Data Analysis Reveals Associations Between Glycosylation and Stemness in Hepatocellular Carcinoma

Background: Glycosylation plays an essential role in driving the progression and treatment resistance of hepatocellular carcinoma (HCC). However, its function in regulating the acquisition and maintenance of the cancer stemness-like phenotype in HCC remains largely unknown. There is also very little known about how CAD and other potential glycosylation regulators may influence stemness. This study explores the relationship between glycosylation and stemness in HCC.

Methods: Gene set variance analysis (GSVA) was used to assess the TCGA pan-cancer enrichment in glycosylation-related pathways. Univariate, LASSO, and multivariate COX regression were then used to identify prognostic genes in the TCGA-LIHC and construct a prognostic signature. HCC patients were classified into high- and low-risk subgroups based on the signature. The relationship between gene expression profiles and stemness was confirmed using bulk and single-cell RNA-sequencing data. The role of CAD and other genes in regulating the stemness of HCC was also validated by RT-qPCR, CCK-8, and colony formation assay. Copy number variation (CNV), immune infiltration, and clinical features were further analyzed in different subgroups and subsequent gene expression profiles. Sensitive drugs were also screened.

Results: In the pan-cancer analysis, HCC was shown to have specific glycosylation alterations. Five genes, CAD, SLC51B, LGALS3, B3GAT3, and MT3, identified from 572 glycosylation-related genes, were used to construct a gene signature and predict HCC patient survival in the TCGA cohort. The results demonstrated a significant positive correlation between patients in the high-risk group and both elevated gene expression and HCC dedifferentiation status. A significant reduction in the stemness-related markers, CD24, CD44, CD20, FOXM1, and EpCAM, was found after the knockdown of CAD and other genes in HepG2 and Huh7 cells. Frequent mutations increased CNVs, immune-suppressive responses, and poor prognosis were also associated with the high-risk profile. The ICGC-LIRI-JP cohort confirmed a similar relationship between glycosylation-related subtypes and stemness. Finally, 84 sensitive drugs were screened for abnormal glycosylation of HCC, and carfilzomib was most highly correlated with CAD.

Conclusions: Glycosylation-related molecular subtypes are associated with HCC stemness and disease prognosis. These results provide new directions for further research on the relationship between glycosylation and stemness phenotypes.

Introduction

Hepatocellular carcinoma (HCC) accounts for 75–85% of primary liver cancers and is the second leading cause of cancer death, with 5-year survival rates of only 4–17% (1). This disease is highly malignant and progresses rapidly, resulting in almost one million deaths each year (2). Surgical resection of HCC followed by chemotherapy is an ideal curative treatment strategy, but it is limited by advanced stage and metastasized tumor cells. While several new therapeutic agents, including checkpoint or tyrosine kinase inhibitors, have been approved by the FDA for patients who cannot undergo surgery or transplantation, their efficacy is unsatisfactory (3, 4), and HCC patients still have an average survival of only 6 months (5, 6). Liver cancer stem cells (LCSCs) are closely associated with the poor prognosis of HCC because they have more robust metastatic and tumorigenic properties than non-LCSCs. Thus, eliminating LCSCs or reducing tumor size is critical to improving HCC treatment efficacy (7). There is an urgent need to explore the molecular biological features of LCSCs.

Glycosylation is a complex process by which a carbohydrate is added to a protein or lipid carrier and is involved in many cellular mechanisms, including cell–cell adhesion, trafficking, motility, inflammation, signaling, host–pathogen interactions, and innate immune responses. All these processes play an essential role in the development and progression of HCC (8, 9). Aberrant glycosylation is frequently cited as a hallmark of malignancy, and from the perspective of epigenetics, glycosylation is thought to directly impact key processes supporting the stemness of HCC, including cell adhesion, motility, invasion, and evasion (9–11). With the presence of highly expressed glycosyltransferases, altered glycosylation is ubiquitous in HCC cells.

Carbamoyl-phosphate synthetase 2, aspartate transcarbamoylase, and dihydrooro-tase (CAD) are multifunctional proteins that play prominent roles in glycosylation. CAD mutations can significantly reduce glycosylation and angiogenesis (12). A product of CAD, uridine diphosphate (UDP), is a specific target for interventional tumors (13, 14). Numerous drugs and anti-tumor vaccines targeting glycosylation are currently in clinical trials. One drug, trastuzumab, is shown to increase the sensitivity of drug-resistant breast cancer cell lines by removing siglec ligands and boosting antibody-dependent natural killer (NK) cell cytotoxicity (15). More glycomic and glycoproteomic studies will help define novel targets and strategies for improving cancer treatment (16). While the B3GAT3, SLC51B, LGALS3, and MT3 genes also play an essential role in the glycosylation pathway, their relationship with stemness remains unclear.

Because of the Warburg effect, sustained high glucose levels in HCC can promote abnormal glycosylation reactions, activate particular signaling pathways, and produce irreversible toxic products, such as glyoxal, methylglyoxal, and 3-deoxyglucosone, that accelerate HCC proliferation and metastasis (17–20). Glycosylation of the stem cell markers, CD24, CD20, CD44, EpCAM, and FOXM1, plays an important role in regulating LCSCs (21–23). Tumor-related glycoprotein or glycan antigen alterations approved by the Food and Drug Administration (FDA), such as core fucosylated AFP (AFP-L3), are better targets for tumor diagnosis and prognosis than AFP alone (24). Further understanding of HCC-related glycosylation patterns will provide advances in treatment and prognosis and reduce mortality.

In this study, glycosylation-related gene sets were downloaded from the GSEA (http://www.gsea-msigdb.org) and a model of HCC prognosis was constructed. Bulk and single-cell RNA-sequencing revealed that the genes and model correlated closely with the stemness of HCC. A systematic analysis of the multi-omics results, including CNV mutations, transcription factors, immunity, and clinical characteristics, may inform further study of the relationships between glycosylation and the stemness of HCC.

Materials and Methods

Pan-Cancer Data Collection and Analysis

All glycosylation-related pathways and relevant gene sets were retrieved from the Molecular Signatures Database (MSigDB, http://www.gsea-msigdb.org/gsea/msigdb/) (Table S1). TCGA pan-cancer RNA-seq data (FPKM values) were downloaded from the genomic data common website (https://gdc.cancer.gov/about‐data/publications/pancanatlas). Gene set variance analysis (GSVA) was used to assess the enrichment of each pan-cancer sample in glycosylation-related pathways, and the distribution of scores was shown using the R “pheatmap” package. Principal component analysis (PCA) was performed after scaling by Z-score, and the pan-cancer RNA-seq data were projected into two dimensions.

HCC Data Collection

In total, 424 RNA-seq transcriptome cases with corresponding clinical information were extracted from the TCGA-LIHC database (https://portal.gdc.cancer.gov) using the Genomic Data Commons (GDC) tool. RNA-seq data were normalized by FPKM. The somatic mutation and CNV data were also downloaded from the TCGA. RNA-seq and clinicopathological data from a Japanese HCC cohort were obtained from the ICGC (LIRI-JP, https://dcc.icgc.org/projects/LIRI-JP) as a comparison.

Construction of a Potential Prognostic Signature

First, 143 differentially expressed genes (DEGs) were identified using LIMMA analysis (adjusted P <0.05 and |Log FC| >1). Univariate Cox regression analysis was then performed to identify potential prognostic DEGs (P <0.05). Least absolute shrinkage and selection operator (LASSO) regression analysis and stepwise Cox regression analysis were used to construct a prognostic model. The risk score of each sample was calculated using the following formula:

$$\text{Risk score}={\displaystyle {\sum }_i^{k}Xi\times Yi\left(X:coefficients, Y:gene expression level\right)}.$$

Based on risk scores, the samples were divided into high- and low-risk subgroups.

Calculation of Stemness-Associated Scores

The mRNA-si stemness score (25) and other stemness signatures (Ben-Porath ESC score (26), Wong ESC score (27), and Bhattacharya ESC score (28)) were used to assess stemness in both the TCGA-LIHC and ICGC-LIRI-JP cohorts.

Single-Cell RNA-Sequencing Data Analysis

Smart-seq2 data (GSE103866) from 55 HuH-1 cells, 63 HuH-7 cells, and 12 patient HCC cells were obtained from the Gene Expression Omnibus (GEO) website. After preprocessing, the scRNA-seq data were converted into a Seurat object, and further quality control was conducted. Data were excluded if 1) cells had >30% mitochondrial genes, 2) genes were identified in <3 cells, or 3) cells had <300 detectable genes. Log normalization, centralization, generation of hypervariable genes, PCA, and clustering analysis were used to perform dimensionality reduction. Uniform manifold approximation and projection (UMAP) was used to visualize the cell distribution using EpCAM, CD24, and CD133 as cell surface markers. “FeaturePlot” and “VlnPlot” were also used to visualize the glycosylation-related genes.

Experimental Validation of the Relationship Between Gene Expressions and Stemness Phenotype

The expression of the HepG2, PLC, Huh7, and Hep3B genes in HCC tumors and cells was determined by quantitative real-time PCR (RT-qPCR) according to the protocols of the manufacturer. Cell Counting Kit-8 (CCK-8) and colony formation assays were used to measure cell proliferation. Small interfering RNA (siRNA) was transfected using Lipofectamine 2000 (Invitrogen, Carlsbad, CA, USA) and mRNA levels were assessed after 48 h. The primer and specific siRNA sequences are included in Table S2.

Analysis of Somatic Mutations and Gene Copy Number Variations (CNVs)

Somatic mutation and CNV data were extracted from the TCGA. The somatic mutation data were visualized using the R “maftools” package. GISTIC2.0 was used to determine the significantly deleted or amplified genomic regions in low- and high-risk subgroups in TCGA-LIHC, referenced to the Consortium Human build 38 (GRCh38). Gene locations were obtained online (ftp://ftp.ensembl.org/pub/current_gtf).

Analysis of Tumor Immune Infiltration

Six previously reported immune subtypes of TCGA-LIHC were identified (29). Infiltrating immune cell fractions were analyzed in tumor samples using single-sample gene set enrichment analysis (ssGSEA) and cell-type identification. Relative subsets of RNA transcript (CIBERSORT) algorithms were estimated in both the TCGA-LIHC and ICGC-LIRI-JP cohorts. Intratumor heterogeneity, IFN-response, TGF-β response, proliferation, and wound healing scores were estimated in the low- and high-risk groups (29).

Analysis of Clinical Data

Kaplan–Meier survival curves were used to evaluate the overall survival (OS) of risk subgroups from the TCGA-LIHC and ICGC-LIRI-JP cohorts. The predictive sensitivity and specificity of the gene signatures were assessed by receiver operating characteristic (ROC). PCA was conducted using the R “prcomp” function of the “stats” package. Univariate and multivariate Cox regression analyses were performed to identify whether clinical characteristics and risk scores were independent risk factors. A nomogram was established to predict the 1-, 3-, and 5-year OS using the R “rms” package. Calibration plots, concordance index (C-index), and ROC were used to evaluate nomogram performance.

Analysis of Drug Sensitivity

The drug sensitivity inhibitory concentration (IC50) value and the mRNA profiles of NCI60 cell lines were extracted from the Cell-Miner™ database (https://discover.nci.nih.gov/cellminer/home.do). The sensitivity of the five selected genes was tested against 216 FDA-approved drugs. Co-expression network analysis of the genes and drugs was visualized using Cytoscape.

Results

HCC Has Specific Glycosylation Alterations in Human Pan-Cancer

A flowchart of this study is shown in Figure 1A. mRNA expression of glycosylation-related genes across 11 cancer types in TCGA was evaluated using PCA. HCC had the most distinct glycosylation patterns (Figure 1B). The GSVA value of glycosylation-related pathways was also calculated in different samples of 11 pan-cancers, and HCC had the most specific enrichment of all cancers (Figure 1C).

FIGURE 1

Figure 1 Glycosylation alterations in human pan-cancer. (A) Study flowchart. (B) PCA projection of paired tumor and normal tissue samples from 11 different cancer types in TCGA. Different colors represent different cancer types. Circles and triangles represent normal and tumor tissue, respectively. (C) Heatmap showing glycosylation-related pathway alterations across 11 cancer types. MsigDB, molecular signatures database; PCA, principal component analysis; HCC, hepatocellular carcinoma; GSVA, gene set variance analysis; CNV, copy number variation; BLCA, bladder urothelial carcinoma; BRCA, breast invasive carcinoma; CESC, cervical squamous cell carcinoma and endocervical adenocarcinoma; HNSC, head and neck squamous cell carcinoma; KIRP, kidney renal papillary cell carcinoma; LIHC, liver hepatocellular carcinoma; LUAD, lung adenocarcinoma; LUSC, lung squamous cell carcinoma; OV, ovarian serous cystadenocarcinoma; PRAD, prostate adenocarcinoma; UCEC, uterine corpus endometrial carcinoma.

Identification of the Candidate Genes and Construction of the Prognostic Signature of HCC

Using gene expression profiling and corresponding clinical information from 50 normal and 374 tumor samples in the TCGA-LIHC database, 143 glycosylation-related DEGs were selected (adjusted P <0.05 and |Log FC| >1). Univariate Cox regression, LASSO regression, and stepwise Cox regression analyses were used to further investigate the importance of these DEGs (Figures 2A, B; Table 1). Five glycosylation-related genes were selected to construct the prognostic signature. Risk score = (0.403 ∗ CAD exp.) + (0.371 ∗ B3GAT3 exp.) + (0.068 ∗ SLC51B exp.) + (0.124 ∗ LGALS3 exp.) + (0.0767 ∗ MT3 exp.).

FIGURE 2

Figure 2 Construction of the prognostic model of the TCGA-LIHC cohort. (A, B) Lasso regression analysis of survival-associated genes. (C) Multivariate Cox regression confirming five glycosylation-related genes. (D) Volcano plot of differentially expressed genes between HCC and non-tumor tissues. The five genes are marked. Red: significant upregulation; blue: significant downregulation; grey: no statistical significance. (E) Expression levels of five glycosylation-related genes between normal and HCC samples. *P <0.05; **P <0.01; *** P <0.001.

TABLE 1

Table 1 Multivariate Cox analysis results of glycosylation-related genes.

Based on the median risk score, HCC patients from the TCGA-LIHC and ICGC-LIRI-JP cohorts were classified into high- and low-risk groups. The five selected genes were identified as risk factors and had higher mRNA expressions in tumor tissues (Figures 2C–E). The protein expression patterns are presented in Supplementary Figure 1.

Gene Expression Profiles of High-Risk Patients Were Enriched With HCC Stemness Markers in TCGC-LIHC and ICGC-LIRI-JP Cohorts

As reported previously, glycosylation is highly correlated with stemness. This study aimed to compare gene expression profiles and risk scores with stemness markers. A mostly positive Spearman’s correlation with the stemness-associated transcriptome-based signatures (Ben-Porath, Wong, Bhattacharya, and mRNA-si) was observed in the TCGA-LIHC (Figure 3A) and ICGC-LIRI-JP cohorts (Supplementary Figure 2A). The stemness markers, CD24, CD44, CD20, FOXM1, and EpCAM, also correlated strongly with gene expression levels and risk scores in both the TCGA-LIHC (Figure 3B) and ICGC-LIRI-JP cohorts (Supplementary Figure 2B). The high-risk cohorts had significantly higher signature scores (Figure 3C; Supplementary Figure 2C) and the gene expression profiles in the TCGA-LIHC and ICGC-LIRI-JP high-risk subgroups were significantly enriched in the Bhattacharya ESC signature (Figure 3D; Supplementary Figure 2D). The TCGA-LIHC and ICGC-LIRI-JP high-risk groups also had a significantly larger proportion of higher stages and grades (Figure 3E; Supplementary Figure 2E), suggesting the presence of a dedifferentiated phenotype.

FIGURE 3

Figure 3 Relationship between gene expression profiles and HCC stemness using HCC bulk data from TCGA-LIHC. (A) Heatmap of Spearman’s correlation results of the gene expression profiles and four distinct stemness indices (Ben-Porath signature, Wong signature, Bhattacharya signature, and mRNA-si). A darker color represents a stronger correlation. (B) Correlation between the CSC markers, CD24, CD44, CD20, FOXM1, and EpCAM, and the gene expression profiles. (C) Different scores of the four distinct stemness indices between the low- and high-risk groups. (D) The transcriptome profiles of high-risk HCC patients were significantly enriched with stemness markers. (E) Among the high-risk patients, the frequency of higher stages and grades was more significantly elevated. Tumor stages and grades were color-coded as shown in the legend. NES, normalized enrichment score; CSC, cancer stem cells.

Expression Profiles of the Five Glycosylation-Related Genes Were Associated With Stemness Markers in LCSCs

To further define the relationship between gene expression profiles and stemness in terms of stem cells, a single-cell RNA-sequencing dataset (GSE103866) including 55 HuH-1 cells, 63 HuH-7 cells, and 12 patient-derived cancer stem cells (CSCs) was downloaded. The UMAP algorithm was adopted for the three cell types, three CSC markers, and the four gene expression distributions (Figures 4A, B). Vlnplot was used to visualize differences in marker gene expression in the distinct immunophenotypes (Figure 4C). Compared with the triple-negative CSCs, these genes showed high expression in other groups with stemness phenotypes (Figure 4D). These results confirm that high expression of these glycosylation-related genes may enhance the stemness of HCC.

FIGURE 4

Figure 4 The relationship between gene expression profiles and HCC stemness in single-cell RNA-sequencing data. (A) UMAP plots of 130 single cells using different classifications. (B) UMAP plots of the four glycosylation-related genes in CSC (MT3 was not detected in GSE 103866). (C) Violin plots of the four glycosylation-related gene expressions in CSC. (D) Expression of the four glycosylation-related genes in different immune phenotypes. CSC, cancer stem cells.

Experimental Validation of the Relationship Between the Five Glycosylation-Related Genes and Stemness

The relationship between the five genes and stemness was also verified in HCC samples. CAD expression was the most different between HCC and normal tissues in the TCGA-LIHC cohort. HCC patient tissues and cell lines also had significantly higher CAD mRNA levels than normal samples (Figures 5A, B). Gene knockdown was conducted in the HepG2 and Huh7 cell lines and the efficiencies were verified (Figure 5C). CAD knockdown resulted in a significant decline in tumor cell viability (P <0.05; Figure 5D) and a marked decrease in the proliferative capacity (Figure 5E). CD24, CD44, CD20, FOXM1, and EpCAM expression were also significantly lower in both si-CAD HepG2 and Huh7 cells than in control cells (Figure 5F). Knockdown of the remaining four genes also significantly inhibited CD24, CD44, CD20, FOXM1, and EpCAM expression (Supplementary Figure 3).

FIGURE 5

Figure 5 The influence of CAD on the HCC stemness phenotype. (A) Differences in CAD expression between normal and HCC tissues. (B) Differential expression of CAD in HCC and normal cells. (C) The efficiency of CAD knockdown in HepG2 and Huh7 cells. (D) CCK-8 experiment comparing the si-CAD and control groups in HepG2 cells. (E) Colony formation assay comparing the si-CAD and control group in HepG2 cells. (F) Significant decrease in stemness-related markers, CD24, CD44, CD20, FOXM1, and EpCAM, after CAD knockdown in HepG2 and Huh7 cells. ***P <0.001.

Somatic Mutation Alterations and CNVs in Different Gene Expression Profiles

Gene mutations in normal stem or progenitor cells may lead to the development and activation of LCSCs and are closely linked to the stemness of HCC. In mutation frequency analyses, a waterfall diagram showed the different status of somatic mutations in the TCGA-LIHC low- and high-risk groups. Sixty percent of genes, including TP53, TTN, MUC16, RYR2, LRP1B, OBSCN, CSMD3, XIRP2, FAT3, CACNA1E, HMCN1, and ARID1A, had a higher mutation frequency in the high-risk than in the low-risk group, while only 20%, including CTNNB1, APOB, ALB, and AXIN1, had a higher frequency in the low-risk group (Figure 6A). The chromosomal locations and CNV alterations of the glycosylation-related genes are shown in Figure 6B. High-risk patients also had a higher number of segments and some mutation scores (29) (Figure 6C). Many chromosomal regions showed significant copy amplification and gene deletion. CNVs were associated with a high-risk prognosis (Figure 6D). GISTIC showed that the scores of both amplification and deletion CNVs were considerably higher in the high-risk group than in the low-risk group (Figure 6E). While the low-risk group had a higher frequency of CNVs than the high-risk group, the high-risk group had more variability (Figure 6F).

FIGURE 6

Figure 6 The high-risk group had a higher frequency of somatic mutations and CNVs in TCGA data. (A) Differences in somatic mutation frequency in the high- and low-risk groups. (B) Five gene copy loss and copy amplification proportion distributions in the genome. (C) Scores of the five distinct mutation indices (number of segments, fraction altered, aneuploidy score, nonslinet mutation rate, and homologous recombination defects) in the low- and high-risk groups. (D–F) CNVs in different risk subgroups. Red and blue represent the two types of CNVs, amplification and deletion, respectively.

The Landscape of Immune Infiltration in Different Gene Expression Profiles

Immune suppression drives tumor evolution toward a stem cell-like phenotype. This study assessed the correlation between risk classification and immune molecular subtype and found that C4 (lymphocyte-depleted) was the most common subtype (Figure 7A). A higher proportion of C1 and C2 subtypes, and a lower level of the C3 subtype, were associated with a poorer prognosis (P <0.0001; Figure 7B). The correlation between immune cell infiltration and glycosylation-related gene expression also varied (Figure 7C). In the TCGA-LIHC cohort, correlation analyses further revealed a clear positive association between risk scores and infiltrating Th2 cells (R = 0.33; P = 3.6e−10) and M0 macrophages (R = 0.36; P = 1.2e−12). In contrast, both Th17 cells (R = −0.31; P = 3e−9) and M1 macrophages (R = −0.13; P = 0.018) were negatively associated with risk scores (Figure 7D). Similar results were also shown in the ICGC-LIRI-JP cohort (Supplementary Figure 4A). Differences in immune cells were also evaluated between the low- and high-risk groups in the TCGA-LIHC and ICGC-LIRI-JP cohorts. While aDCs and macrophages were significantly increased in the high-risk groups in both the TCGA-LIHC and ICGC-LIRI-JP cohorts, B cells, CD8⁺ T cells, mast cells, neutrophils, NK cells, and TIL were reduced in the high-risk score group in the TCGA-LIHC cohort (Figures 7E, F). Meanwhile, only NK cells were reduced in the high-risk group in the ICGC-LIRI-JP cohort (Supplementary Figures 4B, C). Some immune-related scores were also assessed in the low- and high-risk groups (Figure 7G). Except for a striking negative association between the IFN-γ response score and the high-risk subgroup, other scores were significantly elevated in the high-risk group. These results indicated that high-risk patients had a distinct stemness phenotype associated with the accumulation of immune-suppressive cells.

FIGURE 7

Figure 7 The immune landscape of HCC tumors in TCGA-LIHC. (A, B) Classification of low- and high-risk HCC tumors into C1–C6 classes. C1, wound healing; C2, IFN-γ dominant; C3, inflammatory; C4, lymphocyte depleted; C6, TGF-β dominant. (C) Heatmap of Spearman’s correlation between transcriptome profiles, risk scores, and immune infiltration. Only statistically significant correlations are shown (P <0.05). (D) Spearman’s correlation of Th2/Th17 infiltration, M0/M1 macrophage infiltration, and the glycosylation-based risk score. (E, F) Different relative proportions of immune cells in different groups. (G) Scores of the five distinct immune indices between the low- and high-risk groups. *P <0.05; **P <0.01; ***P < 0.001, ns, not significant.

Validation of the Prognostic Signature and Establishment of a Novel Nomogram

HCC patients in the high-risk group had a significantly poorer prognosis and lower OS than those in the low-risk group (Figure 8A; Supplementary Figure 5A). The area under the ROC curve (AUC) of the risk score in predicting 1-, 2-, and 3-year survival was 0.747, 0.741, and 0.730, respectively, in the TCGA dataset and 0.672, 0.642, and 0.647, respectively, in the ICGC dataset (Figure 8B; Supplementary Figure 5B). Compared with the risk subgroups, whole gene expression patterns were separated into two dispersion directions in both the ICGC and TCGA (Figure 8C; Supplementary Figure 5C). Race, gender, age, stage, fibrosis, and risk score were then included in univariate and multivariate Cox regression analyses. Notably, risk score was found to be an independent risk factor (HR >1, P <0.001) (Figure 8D). To further verify the prognostic value of the risk score, a nomogram that included gender, age, stage, and risk score was designed to illustrate patient survival more intuitively (Figure 8E). The calibration curves evaluated the predictive power of the nomogram at 1, 3, and 5 years, and the C-indexes were 0.704, 0.705, and 0.703, respectively (Figure 8F). The ICGC dataset was used to validate the nomogram and showed better discrimination and calibration ability, with C-indexes of 0.737, 0.739, and 0.739, respectively (Figure 8G). The AUC values of the nomogram at 1, 3, and 5 years reached 0.778, 0.748, and 0.741, respectively, which were better than those of the traditional HCC marker, AFP (Figure 8H).

FIGURE 8

Figure 8 Validation of the prognostic signature and establishment of a novel nomogram. (A) Kaplan–Meier curve analysis of the low- and high-risk groups. (B) ROC curve showing the prognostic risk model. (C) PCA plot of TCGA-LIHC patients in different risk groups. (D) Univariate and multivariate independent prognostic analysis of clinical features and risk scores. (E) Nomogram for survival prediction in TCGA-LIHC. (F, G) Nomogram calibration curves of 1-, 3-, and 5-year survival probabilities in TCGA-LIHC and ICGC-LIRI-JP, respectively. (H) ROC curves of the nomogram and AFP for the survival prediction at 1-, 3-, and 5-years. ROC, receiver operating characteristic; AUC, areas under the receiver operating characteristic curve; PCA, Principal component analysis; AFP, alpha-fetoprotein; OS overall survival.

Drug Sensitivity Analysis for the Five Glycosylation-Related Prognostic Genes

Glycosylation provides a range of potential targets for therapeutic intervention. However, potential drugs are still in clinical trials. Thus, the drug sensitivity of five selected glycosylation-related genes across diverse human cancer cell lines was further analyzed (correlation coefficient |R| >0.25, P <0.05). Consequently, 84 sensitive drugs were screened for abnormal glycosylation during HCC. Carfilzomib had the most obvious correlation with CAD (R = −0.341, P = 0.008) and mitoxantrone had the most obvious correlation with LGALS3 (R = −0.517, P = 2.36e−5). Dasatinib is associated with several genes. The top 16 strongest negative correlations between gene expression and IC50 are shown in Figure 9A. The link between the genes and sensitive drugs with a negative correlation of IC50 is shown in Figure 9B.

FIGURE 9

Figure 9 Screening for sensitive drugs. (A) The top 16 most relevant negative correlations. (B) The network between five genes and sensitive drugs with negative IC50s (R<−0.25, P <0.05). Stronger correlations are represented by a thick dark blue line.

Discussion

HCC involves complex architecture and heterogeneity and lacks effective individualized treatment targets (30). LCSCs with a specific phenotype are believed to promote HCC relapse, metastasis, and chemoresistance. Recent studies have shown that aberrant glycosylation of signaling pathways and LCSC markers directly impacts key processes that maintain cell survival, self-renewal, and extravasation properties (31–33). This study focused on the importance of glycosylation in promoting the stemness of HCC by assessing alterations in somatic mutations, immune cell infiltration, and clinical characteristics. Potential target drugs that could be used to improve HCC prognosis were selected.

To our knowledge, this is the first integrated multi-omics study in which the association between glycosylation and stemness was shown in HCC patients using bulk and single-cell RNA-sequencing data. A glycosylation-related prognostic signature was constructed consisting of five genes: CAD, B3GAT3, SLC51B, LGALS3, and MT3, using univariate, stepwise, and multivariate Cox regression analysis. Patients were then divided into low- and high-risk subgroups. Gene expression profiles of the high-risk group correlated positively with the upregulation of CSC markers, CD24, CD44, CD20, FOXM1, and EpCAM, and significant enrichment of other ESC signatures, indicating that these were tumor-promoting targets. Several studies have shown that stemness is associated with genetic mutations, epigenetic changes, and differences in the tumor microenvironment (34, 35). Similarly, our study found that the high-risk group had unique somatic mutations, CNVs, and immune patterns as well as enriched stem cell-like characteristics and high levels of gene expression. Clinical information was then combined with the prognostic signature to construct a calibrated nomogram. Finally, the sensitivity of glycosylation-related genes to particular drugs was evaluated, providing novel insight into tumor treatments and the prevention of drug resistance. Among the drug candidates, carfilzomib showed the most obvious correlation with CAD.

The five glycosylation-related genes were identified as significant in HCC. CAD is a multifunctional protein that takes part in de novo pyrimidine nucleotide synthesis, protein glycosylation, and phospholipid biosynthesis in mammals (12). This protein must initiate the di-(UDP)-dependent glycosylation process by producing UDP. UDP-N-acetylglucosamine (UDP-GlcNAc), for example, serves as an essential sugar donor substrate for O-GlcNAc, which is involved in a central intracellular protein modification in diverse metabolism, signaling, and disease processes (36). A recent study by Ching-Yu et al. (37) showed that UDP overexpression significantly promoted human hepatoma cell proliferation both in vitro and in vivo. UDP-13C-glucose flux is also responsible for metabolic reprogramming and the expression of key stemness genes in CSCs (38). The expression of LGALS3, a siglec-9 ligand in cancer cells, correlates with tumor progression and increased expression of β-catenin and CSC markers induced through Wnt signaling (39, 40). LGALS3 also helps maintain LSCS stemness, expansion, and aggressiveness and may thus serve as a target for HCC treatment (41). B3GAT3 is a glucuronosyltransferase involved in glycosylation and the proliferation and metastasis of HCC tissues and cells (42). SLC51B is involved in the intestinal reclamation of bile acids and steroids and eicosanoid metabolism, promoting liver cancer cell proliferation and suppressing apoptosis, which is associated with stemness (43). MT3 is a member of the metallothionein family that regulates protein glycosylation and is closely linked to HCC progression. Studies indicate that MT overexpression may induce tumor cell differentiation (44, 45), but the mechanism will require further investigation. This study is consistent with the previously reported findings while also furthering our understanding of the relationship between LCSCs and glycosylation.

The production and activation of LCSCs may be induced by a series of gene mutations in normal stem or progenitor cells that disrupt cell metabolism, immune escape, and drug resistance. CNV is a hallmark of amplifications or deletions in the cancer genome, which can inactivate tumor suppressor genes, induce high oncogene expression, and increase stemness. A recent study of CNV showed that increases in the E2f1 or E2f3b loci promoted spontaneous HCC, whereas decreases in these loci suppressed HCC (46). In a genome-wide analysis, TRIM35 was shown as a novel tumor suppressor (47), and COL4A1 on the 13q34 locus was a frequent amplification target, a finding consistent with our results (48). Our study found that the RB1 deletion was most significant in the high-risk group, which supports prior research by Sung-Min et al. (49).

Stemness is significantly associated with immunosuppression in the immune microenvironment. As a result, this study also compared the proportions of immune cells in the high- and low-risk score groups and their gene expression profiles. Interestingly, cancer stem cell-like characteristics were enriched in the high-risk groups. To characterize intratumoral immune states, HCC samples were clustered into six subtypes using definitions developed by Thorsson et al. (29). The C4 subtype, associated with an imbalanced ratio of Th2/Th17 and M0/M1 cells, was the most common. Th2-related cytokines such as IL-4, IL-5, and IL-10 are involved in stimulating B cell proliferation and mediating humoral immunity. Some studies indicate that Th2 cells induce HCC progression (50). Indeed, Th2 cytokines facilitate tumor escape by inhibiting Th1 cytokine production (51) and furthering the progression of liver damage. However, IL-17 production by Th17 cells is generally associated with higher cancer survival (52). After transcatheter arterial embolization (TAE), the number of Th17 cells in the tumor microenvironment increases significantly (53). The Th2/Th17 cytokine profile is unbalanced during HCC and the shift towards Th17 cells can promote anti-tumor effects. Macrophages are key antigen-presenting cells of the innate immune system, and different macrophage phenotypes have distinct functions in regulating the tumor microenvironment. In particular, M1 macrophages are involved in suppressing tumor growth and inducing liver tumor regression. Guerra et al. (54) found that hydrogel-embedded M1 macrophages induced apoptosis of HCC cells and promoted tumor regression. Localization of activated M1 macrophages may be a new novel HCC treatment strategy. Our study also found significant enrichment in intratumor heterogeneity, TGF-β mediated signaling, proliferation, and a cancer stemness phenotype in the high-risk group, while IFN signaling was increased in the low-risk group. Consistent with these results, TGF-β-mediated signaling has been strongly linked to a cancer stem-like phenotype in tumors (55), while IFN signaling has been negatively correlated with stemness and HCC progression (56, 57).

Current chemotherapy drugs remain unable to significantly promote the long-term survival of HCC patients. Thus, this study screened potential target drugs using the Cell-Miner database. Carfilzomib, a second-generation proteasome inhibitor, had the strongest correlation with CAD. This drug exhibits an anti-tumor effect by inhibiting MAPK signaling (58). CAD activation requires MAPK, so it is possible that carfilzomib indirectly inhibits CAD by suppressing MAPK signaling and promoting an anti-tumor response. Mitoxantrone, an antineoplastic drug developed in the 1980s, had the strongest correlation with genes. A recent study found that mitoxantrone inhibits HCC growth and proliferation by inducing autophagy (59). In the mitoxantrone-resistant cell line, LGALS3 was detected with an abundance ratio well above two (60). Dasatinib, a tyrosine kinase inhibitor that exerts anti-tumor effects by inhibiting glucuronosyltransferase, is associated with the expression of several genes (61). SLC51B is a subtype of the human solute carrier (SLC) and correlates negatively with drug sensitivity (62). MT3 is also involved in SLC synthesis and is negatively associated with Dasatinib. BAGAT3 is an important type of glucuronosyltransferase. However, whether these drugs could be used to treat abnormal glycosylation requires further study.

Glycosylation is a promising target for tumor treatment. However, there are some limitations that should be addressed in our study. First, few studies have assessed the relationship between stemness and glycosylation, especially those involving the development of HCC. This study relied heavily on TCGA, GEO, and ICGC data and lacked experimental evidence. Thus, we could not determine whether the selected genes have corresponding roles in glycosylation pathways. Further functional experiments are necessary to define the molecular mechanisms. Second, CNVs were not available in the ICGC database. A larger CNV sample size would help ensure the reliability of the study. Third, the mechanism of potentially sensitive drugs in abnormal glycosylation during HCC remains unknown and requires further study.

Conclusions

In conclusion, this study proposed a prognostic signature founded on glycosylation-related gene expression and classified HCC samples into low- and high-risk subgroups. According to bulk and single-cell RNA-sequencing, the high-risk group was most likely to have a stem-like phenotype. This was verified by assessing somatic mutations, CNVs, immune infiltration, and clinical characteristics. These findings will help inform future studies on the relationship between glycosylation and HCC stemness.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

Ethics Statement

The studies involving human participants were reviewed and approved by the Ethics Committee of Tianjin Second People’s Hospital. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

All authors contributed to the study conception and design. PL and QZ contributed equally to this study. PL and QZ designed and wrote the manuscript. PL and QZ screened the literature and collected data. PL and QZ were responsible for pictures and statistics. PL and QZ polished the article. JL critically revised the manuscript. All authors listed have made a substantial, direct, and intellectual contribution to the work and approved it for publication.

Funding

This work was funded by the Natural Science Foundation of Tianjin City (20JCYBJC01150); the Tianjin Health Science and Technology Project (Nos. TJWJ2021QN063, TJWJ2021ZD010, and TJWJ2021MS034), and the Tianjin Key Medical Discipline (Specialty) Construction Project.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

The authors acknowledge The Third Affiliated Hospital of Sun Yat-Sen University for their technical assistance.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2022.913432/full#supplementary-material

Supplementary Table 1 | Glycosylation-related pathways and genes retrieved from MSigDB. MSigDB, Molecular signatures Database.

Supplementary Table 2 | Primer and specific siRNA sequences were used in this study. siRNA, small interfering RNA.

Supplementary Figure 1 | The representative immunohistochemistry of the five glycosylation-related genes in normal and HCC tissues downloaded from HPA. HPA, The Human Protein Atlas.

Supplementary Figure 2 | Relationship between gene expression profiles and HCC stemness using HCC bulk data from ICGC-LIRI-JP. (A) Heatmap of Spearman’s correlation results of the gene expression profiles and four distinct stemness indices (Ben-Porath signature, Wong signature, Bhattacharya signature, and mRNA-si). A darker color represents a stronger correlation. (B) Correlation between the CSC markers, CD24, CD44, CD20, FOXM1, and EpCAM, and the gene expression profiles. (C) Different scores of the four distinct stemness indices between the low- and high-risk groups. (D) The transcriptome profiles of high-risk HCC patients were significantly enriched with stemness markers. (E) Among the high-risk patients, the frequency of higher stages and grades was more significantly elevated. Tumor stages and grades were color-coded as shown in the legend. NES, normalized enrichment score; CSC, cancer stem cells.

Supplementary Figure 3 | The influence of the remaining four genes on the stemness phenotype of HCC. (A) Differential expression in HCC and normal cells. (B) Significant decrease in stemness-related markers, CD24, CD44, CD20, FOXM1, and EpCAM, after SLC51B knockdown in HepG2 and Huh7 cells. (C) Significant decrease in stemness-related markers, CD24, CD44, CD20, FOXM1, and EpCAM, after LGALS3 knockdown in HepG2 and Huh7 cells. (D) Significant decrease in stemness-related markers, CD24, CD44, CD20, FOXM1, and EpCAM, after B3GAT3 knockdown in HepG2 and Huh7 cells. (E) Significant decrease in stemness-related markers, CD24, CD44, CD20, FOXM1, and EpCAM, after MT3 knockdown in HepG2 and Huh7 cells. *** P < 0.001.

Supplementary Figure 4 | The immune landscape of HCC tumors in ICGC-LIRI-JP. (A) Spearman’s correlation of Th2/Th17 infiltration, M0/M1 macrophage infiltration, and the glycosylation-based risk score. (B, C) Different relative proportions of immune cells in different groups. *P <0.05; **P <0.01; *** P <0.001.

Supplementary Figure 5 | Validation of the prognostic signature in ICGC-LIRI-JP. (A) Kaplan-Meier curve analysis of the low- and high-risk groups. (B) ROC curve showing the prognostic risk model. (C) PCA plot of TCGA-LIHC patients in different risk groups. ROC, receiver operating characteristic; PCA, principal component analysis.

References

1. Emal A, Ward EM, Johnson CJ, Cronin KA, Ma J, Ryerson B, et al. Annual Report to the Nation on the Status of Cancer, 1975-2014, Featuring Survival. J Natl Cancer Inst (2017) 109(9):djx030. doi: 10.1093/jnci/djx030

CrossRef Full Text | Google Scholar

2. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, et al. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin (2021) 71(3):209–49. doi: 10.3322/caac.21660

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Yu M, Chen Z, Zhou Q, Zhang B, Huang J, Jin L, et al. PARG InhibitionLimits HCC Progression and Potentiates the Efficacy of Immune CheckpointTherapy. J Hepatol (2022) 12:S0168–8278(22)00072. doi: 10.1016/j.jhep.2022.01.026

CrossRef Full Text | Google Scholar

4. Romito I, Porru M, Braghini MR, Pompili L, Panera N, Crudele A, et al. Focal Adhesion Kinase Inhibitor TAE226 Combined With Sorafenib Slows Down Hepatocellular Carcinoma by Multiple Epigenetic Effects. J Exp Clin Cancer Res (2021) 40(1):364. doi: 10.1186/s13046-021-02154-8

PubMed Abstract | CrossRef Full Text | Google Scholar

5. European Association for the Study of the Liver. EASL-EORTC Clinical Practice Guidelines: Management of Hepatocellular Carcinoma. J Hepatol (2012) 56(4):908–43. doi: 10.1016/j.jhep.2011.12.001

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Bruix J, Reig M, Sherman M. Evidence-Based Diagnosis, Staging, and Treatment of Patients With Hepatocellular Carcinoma. Gastroenterology (2016) 150(4):835–53. doi: 10.1053/j.gastro.2015.12.041

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Tricot T, Verfaillie CM, Kumar M. Current Status and Challenges of Human Induced Pluripotent Stem Cell-Derived Liver Models in Drug Discovery. Cells (2022) 11(3):442. doi: 10.3390/cells11030442

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Zhang S, Cao X, Gao Q, Liu Y. Protein Glycosylation in Viral Hepatitis-Related HCC: Characterization of Heterogeneity, Biological Roles, and Clinical Implications. Cancer Lett (2017) 406:64–70. doi: 10.1016/j.canlet.2017.07.026

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Jia L, Li J, Li P, Liu D, Li J, Shen J, et al. Site-Specific Glycoproteomic Analysis Revealing Increased Core-Fucosylation on FOLR1 Enhances Folate Uptake Capacity of HCC Cells to Promote EMT. Theranostics (2021) 11(14):6905–21. doi: 10.7150/thno.56882

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Zhang Q, Jiang K, Li Y, Gao D, Sun L, Zhang S, et al. Histidine-Rich Glycoprotein Function in Hepatocellular Carcinoma Depends on its N-Glycosylation Status, and it Regulates Cell Proliferation by Inhibiting Erk1/2 Phosphorylation. Oncotarget (2015) 6(30):30222–31. doi: 10.18632/oncotarget.4997

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Zhang C, Niu Y, Wang Z, Xu X, Li Y, Ma L, et al. Corosolic Acid Inhibits Cancer Progression by Decreasing the Level of CDK19-Mediated O-GlcNAcylation in Liver Cancer Cells. Cell Death Dis (2021) 12(10):889. doi: 10.1038/s41419-021-04164-y

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Coxam B, Neyt C, Grassini DR, Le Guen L, Smith KA, Schulte-Merker S, et al. Carbamoyl-Phosphate Synthetase 2, Aspartate Transcarbamylase, and Dihydroorotase (Cad) Regulates Notch Signaling and Vascular Development in Zebrafish. Dev Dyn (2015) 244(1):1–9. doi: 10.1002/dvdy.24209

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Teoh ST, Ogrodzinski MP, Lunt SY. UDP-Glucose 6-Dehydrogenase Knockout Impairs Migration and Decreases In Vivo Metastatic Ability of Breast Cancer Cells. Cancer Lett (2020) 492:21–30. doi: 10.1016/j.canlet.2020.07.031

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Wolfe AL, Zhou Q, Toska E, Galeas J, Ku AA, Koche RP, et al. UDP-Glucose Pyrophosphorylase 2, a Regulator of Glycogen Synthesis and Glycosylation, is Critical for Pancreatic Cancer Growth. Proc Natl Acad Sci U.S.A. (2021) 118(31):e2103592118. doi: 10.1073/pnas.2103592118

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Xiao H, Woods EC, Vukojicic P, Bertozzi CR. Precision Glycocalyx Editing as a Strategy for Cancer Immunotherapy. Proc Natl Acad Sci U.S.A. (2016) 113(37):10304–9. doi: 10.1073/pnas.1608069113

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Pinho SS, Reis CA. Glycosylation in Cancer: Mechanisms and Clinical Implications. Nat Rev Cancer (2015) 15(9):540–55. doi: 10.1038/nrc3982

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Ashraf JM, Shahab U, Tabrez S, Lee EJ, Choi I, Aslam yusuf M, et al. DNA Glycation From 3-Deoxyglucosone Leads to the Formation of AGEs: Potential Role in Cancer Auto-Antibodies. Cell Biochem Biophys (2016) 74(1):67–77. doi: 10.1007/s12013-015-0713-6

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Khan MS, Rabbani N, Tabrez S, Islam BU, Malik A, Ahmed A, et al. Glycation Induced Generation of Amyloid Fibril Structures by Glucose Metabolites. Protein Pept Lett (2016) 23(10):892–7. doi: 10.2174/0929866523666160831153858

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Ahmad S, Akhter F, Shahab U, Rafi Z, Khan MS, Nabi R, et al. Do All Roads Lead to the Rome? The Glycation Perspective! Semin Cancer Biol (2018) 49:9–19. doi: 10.1016/j.semcancer.2017.10.012

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Ashraf JM, Ahmad S, Rabbani G, Hasan Q, Jan AT, Lee EJ, et al. 3-Deoxyglucosone: A Potential Glycating Agent Accountable for Structural Alteration in H3 Histone Protein Through Generation of Different AGEs. PloS One (2015) 10(2):e0116804. doi: 10.1371/journal.pone.0116804

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Shen YA, Wang CY, Chuang HY, Hwang JJ, Chi WH, Shu CH, et al. CD44 and CD24 Coordinate the Reprogramming of Nasopharyngeal Carcinoma Cells Towards a Cancer Stem Cell Phenotype Through STAT3 Activation. Oncotarget (2016) 7(36):58351–66. doi: 10.18632/oncotarget.11113

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Yamashita T, Honda M, Nakamoto Y, Baba M, Nio K, Hara Y, et al. Discrete Nature of EpCAM+ and CD90+ Cancer Stem Cells in Human Hepatocellular Carcinoma. Hepatology (2013) 57(4):1484–97. doi: 10.1002/hep.26168

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Yuan B, Liu Y, Yu X, Yin L, Peng Y, Gao Y, et al. FOXM1 Contributes to Taxane Resistance by Regulating UHRF1-Controlled Cancer Cell Stemness. Cell Death Dis (2018) 9(5):562. doi: 10.1038/s41419-018-0631-9

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Ahn KS, O'Brien DR, Kim YH, Kim TS, Yamada H, Park JW, et al. Associations of Serum Tumor Biomarkers With Integrated Genomic and Clinical Characteristics of Hepatocellular Carcinoma. Liver Cancer (2021) 10(6):593–605. doi: 10.1159/000516957

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Malta TM, Sokolov A, Gentles AJ, Burzykowski T, Poisson L, Weinstein JN, et al. Machine Learning Identifies Stemness Features Associated With Oncogenic Dedifferentiation. Cell (2018) 173(2):338–354.e15. doi: 10.1016/j.cell.2018.03.034

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Ben-Porath I, Thomson MW, Carey VJ, Ge R, Bell GW, Regev A, et al. An Embryonic Stem Cell-Like Gene Expression Signature in Poorly Differentiated Aggressive Human Tumors. Nat Genet (2008) 40(5):499–507. doi: 10.1038/ng.127

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Wong DJ, Liu H, Ridky TW, Cassarino D, Segal E, Chang HY. Module Map of Stem Cell Genes Guides Creation of Epithelial Cancer Stem Cells. Cell Stem Cell (2008) 2(4):333–44. doi: 10.1016/j.stem.2008.02.009

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Bhattacharya B, Miura T, Brandenberger R, Mejido J, Luo Y, Yang AX, et al. Gene Expression in Human Embryonic Stem Cell Lines: Unique Molecular Signature. Blood (2004) 103(8):2956–64. doi: 10.1182/blood-2003-09-3314

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Thorsson V, Gibbs DL, Brown SD, Wolf D, Bortone DS, Ou Yang TH, et al. The Immune Landscape of Cancer. Immunity (2018) 48(4):812–830.e14. doi: 10.1016/j.immuni.2018.03.023

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Losic B, Craig AJ, Villacorta-Martin C, Martins-Filho SN, Akers N, Chen X, et al. Intratumoral Heterogeneity and Clonal Evolution in Liver Cancer. Nat Commun (2020) 11(1):291. doi: 10.1038/s41467-019-14050-z

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Li C, Du Y, Yang Z, He L, Wang Y, Hao L, et al. GALNT1-Mediated Glycosylation and Activation of Sonic Hedgehog Signaling Maintains the Self-Renewal and Tumor-Initiating Capacity of Bladder Cancer Stem Cells. Cancer Res (2016) 76(5):1273–83. doi: 10.1158/0008-5472.CAN-15-2309

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Venetz D, Hess C, Lin CW, Aebi M, Neri D. Glycosylation Profiles Determine Extravasation and Disease-Targeting Properties of Armed Antibodies. Proc Natl Acad Sci U.S.A. (2015) 112(7):2000–5. doi: 10.1073/pnas.1416694112

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Zhou B, Zhao YC, Liu H, Luo S, Amos CI, Lee JE, et al. Novel Genetic Variants of ALG6 and GALNTL4 of the Glycosylation Pathway Predict Cutaneous Melanoma-Specific Survival. Cancers (Basel) (2020) 12(2):288. doi: 10.3390/cancers12020288

CrossRef Full Text | Google Scholar

34. Weng CC, Ding PY, Liu YH, Hawse JR, Subramaniam M, Wu CC, et al. Mutant Kras-Induced Upregulation of CD24 Enhances Prostate Cancer Stemness and Bone Metastasis. Oncogene (2019) 38(12):2005–19. doi: 10.1038/s41388-018-0575-7

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Miranda A, Hamilton PT, Zhang AW, Pattnaik S, Becht E, Mezheyeuski A, et al. Cancer Stemness, Intratumoral Heterogeneity, and Immune Response Across Cancers. Proc Natl Acad Sci U.S.A. (2019) 116(18):9020–9. doi: 10.1073/pnas.1818210116

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Winans KA, Bertozzi CR. An Inhibitor of the Human UDP-GlcNAc 4-Epimerase Identified From a Uridine-Based Library: A Strategy to Inhibit O-Linked Glycosylation. Chem Biol (2002) 9(1):113–29. doi: 10.1016/s1074-5521(02)00093-5

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Lai CY, Liu H, Tin KX, Huang Y, Yeh KH, Peng HW, et al. Identification of UAP1L1 as a Critical Factor for Protein O-GlcNAcylation and Cell Proliferation in Human Hepatoma Cells. Oncogene (2019) 38(3):317–31. doi: 10.1038/s41388-018-0442-6

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Bi Z, Zhang Q, Fu Y, Wadgaonkar P, Zhang W, Almutairy B, et al. Nrf2 and HIF1α Converge to Arsenic-Induced Metabolic Reprogramming and the Formation of the Cancer Stem-Like Cells. Theranostics (2020) 10(9):4134–49. doi: 10.7150/thno.42903

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Ilmer M, Mazurek N, Gilcrease MZ, Byrd JC, Woodward WA, Buchholz TA, et al. Low Expression of Galectin-3 is Associated With Poor Survival in Node-Positive Breast Cancers and Mesenchymal Phenotype in Breast Cancer Stem Cells. Breast Cancer Res (2016) 18(1):9. doi: 10.1186/s13058-016-0757-6

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Chung LY, Tang SJ, Wu YC, Sun GH, Liu HY, Sun KH. Galectin-3 Augments Tumor Initiating Property and Tumorigenicity of Lung Cancer Through Interaction With β-Catenin. Oncotarget (2015) 6(7):4936–52. doi: 10.18632/oncotarget.3210

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Tummala KS, Brandt M, Teijeiro A, Graña O, Schwabe RF, Perna C, et al. Hepatocellular Carcinomas Originate Predominantly From Hepatocytes and Benign Lesions From Hepatic Progenitor Cells. Cell Rep (2017) 19(3):584–600. doi: 10.1016/j.celrep.2017.03.059

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Zhang YL, Ding C, Sun L. High Expression B3GAT3 Is Related With Poor Prognosis of Liver Cancer. Open Med (Wars) (2019) 14:251–8. doi: 10.1515/med-2019-0020

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Cheng J, Li Y, Wang X, Dong Z, Chen Y, Zhang R, et al. Response Stratification in the First-Line Combined Immunotherapy of Hepatocellular Carcinoma at Genomic, Transcriptional and Immune Repertoire Levels. J Hepatocell Carcinoma (2021) 8:1281–95. doi: 10.2147/JHC.S326356

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Jin R, Bay BH, Chow VT, Tan PH, Lin VC. Metallothionein 1e mRNA is Highly Expressed in Oestrogen Receptor-Negative Human Invasive Ductal Breast Cancer. Br J Cancer (2000) 83(3):319–23. doi: 10.1054/bjoc.2000.1276

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Sens MA, Somji S, Garrett SH, Beall CL, Sens DA. Metallothionein Isoform 3 Overexpression is Associated With Breast Cancers Having a Poor Prognosis. Am J Pathol (2001) 159(1):21–6. doi: 10.1016/S0002-9440(10)61668-9

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Kent LN, Bae S, Tsai SY, Tang X, Srivastava A, Koivisto C, et al. Dosage-Dependent Copy Number Gains in E2f1 and E2f3 Drive Hepatocellular Carcinoma. J Clin Invest (2017) 127(3):830–42. doi: 10.1172/JCI87583

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Jia D, Wei L, Guo W, Zha R, Bao M, Chen Z, et al. Genome-Wide Copy Number Analyses Identified Novel Cancer Genes in Hepatocellular Carcinoma. Hepatology (2011) 54(4):1227–36. doi: 10.1002/hep.24495

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Lau SH, Guan XY. Cytogenetic and Molecular Genetic Alterations in Hepatocellular Carcinoma. Acta Pharmacol Sin (2005) 26(6):659–65. doi: 10.1111/j.1745-7254.2005.00126.x

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Ahn SM, Jang SJ, Shim JH, Kim D, Hong SM, Sung CO, et al. Genomic Portrait of Resectable Hepatocellular Carcinomas: Implications of RB1 and FGF19 Aberrations for Patient Stratification. Hepatology (2014) 60(6):1972–82. doi: 10.1002/hep.27198

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Gramenzi A, Andreone P, Loggi E, Foschi FG, Cursaro C, Margotti M, et al. Cytokine Profile of Peripheral Blood Mononuclear Cells From Patients With Different Outcomes of Hepatitis C Virus Infection. J Viral Hepat (2005) 12(5):525–30. doi: 10.1111/j.1365-2893.2005.00634.x

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Ji L, Gu J, Chen L, Miao D. Changes of Th1/Th2 Cytokines in Patients With Primary Hepatocellular Carcinoma After Ultrasound-Guided Ablation. Int J Clin Exp Pathol (2017) 10(8):8715–20. doi: 10.4161/2162402X.2014.984547

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Punt S, Langenhoff JM, Putter H, Fleuren GJ, Gorter A, Jordanova ES. The Correlations Between IL-17 vs Th17 Cells and Cancer Patient Survival: A Systematic Review. Oncoimmunology (2015) 4(2):e984547. doi: 10.4161/2162402X.2014.984547

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Avritscher R, Jo N, Polak U, Cortes AC, Nishiofuku H, Odisio BC, et al. Hepatic Arterial Bland Embolization Increases Th17 Cell Infiltration in a Syngeneic Rat Model of Hepatocellular Carcinoma. Cardiovasc Intervent Radiol (2020) 43(2):311–21. doi: 10.1007/s00270-019-02343-1

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Guerra AD, Yeung OWH, Qi X, Kao WJ, Man K. The Anti-Tumor Effects of M1 Macrophage-Loaded Poly (Ethylene Glycol) and Gelatin-Based Hydrogels on Hepatocellular Carcinoma. Theranostics (2017) 7(15):3732–44. doi: 10.7150/thno.20251

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Johnston CJ, Smyth DJ, Dresser DW, Maizels RM. TGF-β in Tolerance, Development and Regulation of Immunity. Cell Immunol (2016) 299:14–22. doi: 10.1016/j.cellimm.2015.10.006

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Shirasaki T, Honda M, Yamashita T, Nio K, Shimakami T, Shimizu R, et al. The Osteopontin-CD44 Axis in Hepatic Cancer Stem Cells Regulates IFN Signaling and HCV Replication. Sci Rep (2018) 8(1):13143. doi: 10.1038/s41598-018-31421-6

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Xing M, Li J. A New Inflammation-Related Risk Model for Predicting Hepatocellular Carcinoma Prognosis. BioMed Res Int (2022) 5:5396128. doi: 10.1155/2022/5396128

CrossRef Full Text | Google Scholar

58. Jiang C, Xu R, Li XX, Zhou YF, Xu XY, Yang Y, et al. Sorafenib and Carfilzomib Synergistically Inhibit the Proliferation, Survival, and Metastasis of Hepatocellular Carcinoma. Mol Cancer Ther (2018) 17(12):2610–21. doi: 10.1158/1535-7163.MCT-17-0541

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Xie B, He X, Guo G, Zhang X, Li J, Liu J, et al. High-Throughput Screening Identified Mitoxantrone to Induce Death of Hepatocellular Carcinoma Cells With Autophagy Involvement. Biochem Biophys Res Commun (2020) 521(1):232–7. doi: 10.1016/j.bbrc.2019.10.114

PubMed Abstract | CrossRef Full Text | Google Scholar

60. Wang J, Gutierrez P, Edwards N, Fenselau C. Integration of 18O Labeling and Solution Isoelectric Focusing in a Shotgun Analysis of Mitochondrial Proteins. J Proteome Res (2007) 6(12):4601–7. doi: 10.1021/pr070401e

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Virakul S, Dalm VA, Paridaens D, van den Bosch WA, Hirankarn N, van Hagen PM, et al. The Tyrosine Kinase Inhibitor Dasatinib Effectively Blocks PDGF-Induced Orbital Fibroblast Activation. Graefes Arch Clin Exp Ophthalmol (2014) 252(7):1101–9. doi: 10.1007/s00417-014-2674-7

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Polillo M, Galimberti S, Baratè C, Petrini M, Danesi R, Di Paolo A. Pharmacogenetics of BCR/ABL Inhibitors in Chronic Myeloid Leukemia. Int J Mol Sci (2015) 16(9):22811–29. doi: 10.3390/ijms160922811

PubMed Abstract | CrossRef Full Text | Google Scholar