Deciphering the Prognostic Implications of T Cells Marker Genes Signature in the tumor Microenvironment of Hepatocellular Carcinoma

DOI: https://doi.org/10.21203/rs.3.rs-1535216/v1

Abstract

Background: Hepatocellular carcinoma is one of the most aggressive gastrointestinal tumor, with a high recurrence and mortality rate. Novel immunotherapy targeting the tumor microenvironment are innovative and promising therapeutic approaches for some tumor, including hepatocellular carcinoma. However, patr of patients with HCC remain a diminished response to immunotherapy on account of the insufficient understanding on the tumor microenvironment.

Methods: The correlation between infiltrated immune cells and hepatocellular carcinoma patients prognosis were confirmed via MCP-counter algorithms. The single-cell RNA sequencing dataset (GSE146115) was performed to identify T cells marker genes in hepatocellular carcinoma. GSE10141, GSE14520 dataset and TCGA-liver hepatocellular carcinoma were used to construct and validate the T cell marker genes signature (TCMS).

Results: The survival analysis revealed the correlation between infiltrating T cells and the overall survival of patients with hepatocellular carcinoma. Four T cell marker genes included in the TCMS model were: HELZ, GZMA, SLC2A2, JAK3. The TCMS risk score could stratify patients with hepatocellular carcinoma into high- and low-risk score subgroups. TCMS risk score remained a influential feature of overall survival and one of the independent prognostic risk factors in multivariate analysis in TCGC validation cohort. Besides, correlation analysis indicated increased expression of HELZ, GZMA and JAK3 were significantly correlated with high degree of CD8+ T cell, CD4+ T cell infiltration.

Conclusion: We successfully constructed the T cell marker genes signature with powerful predictive function and provided a new understanding of T cell infiltration in tumor microenvironment which might offer practices instructions for HCC immunotherapy.

1 Introduction

Liver cancer, predominantly hepatocellular carcinoma (HCC), is a life-threatening malignancy with the third leading cause of cancer-related mortality nowadays.1 Although clinical management modalities for HCC have made great progress in recent years, the five-year overall survival rate of liver tumor patients in most countries is only 10–19%, which remains a growing concern.2,3 Accumulating scientific evidence indicate tumor-immunity interactions along with tumor microenvironment (TME) that develop during disease evolution play a pivotal role in patients’ prognosis and therapeutic response.4 Given that HCCs are strongly and wildly resistant to conventional chemotherapy, immunotherapy with immune checkpoint inhibitors (ICIs) had revolutionized the treatment modalities in the clinic for HCC.5 The combination of anti-programmed death-ligand 1 (PD-L1) antibody (atezolizumab) and vascular endothelial growth factor (VEGF) antibody (bevacizumab) were recommended as the first-line therapy for advanced HCC.6,7

T cells repertoire of HCC represent a complex interplay between immune-activating and immune-suppressive factors in TME, that reflected in the infiltration degree of T cell subsets, including CD4 + T cell (T-helper cells) and CD8 + T cells.8,9 There were no significant difference noted in immune contexture of HCCs with different etiologies, however the selective reduction of CD8 + T cells and the decrease in co-stimulatory activity and cytolytic were observed in more advanced HCC.10 CD4 + T cells are responsible for modulating the activity of cytotoxic T-cell and CD8 + T cells, also called cytotoxic T cells, are indispensable for immunosurveillance. The exhausted of CD8 + T cells impaired tumor immune surveillance and form a tumor-promoting milieu.11 Indeed, HCC patients with higher CD8 + T cells infiltration overexpressing coinhibitory receptors like PD-1, cytotoxic T-lymphocyte-associated protein 4 (CTLA-4), lymphocyte-activation gene 3 (LAG-3) and T cell immunoglobulin domain and mucin domain-3 (TIM-3) show an exhausted phenotype and poor clinical outcome.12 Hence, elucidating the mechanism of components and signatures, especially T cells, in HCC-TME were significant and imperative.

Advances in single-cell RNA sequencing (scRNA-seq) technology provided high-dimensional single-cell data to present the molecular characteristics of varied cell components and afforded a high-throughput method for researchers to interpret intratumoral heterogeneity. Therefore, scRNA-seq analysis of immune cells in HCC-TME helped to discover specific signatures served as potential targets of immunotherapy.13,14 The stringent limitation of scRNA-seq analysis was the difficulty of combining the sequencing findings with clinical information, such as tumor staging and survival expectancy. Appropriate combination with the merit of scRNA-seq and bulk sequencing would optimize the utilization of scRNA-seq data, giving the availability and complete clinical information in sequencing cohorts. Here, we annotated the T cell clusters in HCC using scRNA-seq data and combined sequencing cohorts from Gene Expression Omnibus (GEO) and The Cancer Genome Atlas (TCGA) to construct and validate the T cell marker genes signature. The correlation with T cells infiltration were further confirmed in Tumor IMmune Estimation Resource database.

2 Methods

2.1 Data source

The single-cell RNA sequencing dataset (GSE146115) counts from four hepatocellular carcinoma was downloaded from GEO. The sequencing datasets were obtained from GEO (GSE10141, GSE14520) and (TCGA-LIHC). Concurrently,both the individual patient files and transcriptome raw data of each dataset were derived from the respective websites, Only HCC tissues with completed follow-up information and survival status were included in the further analysis, while indistinct clinical information and non tumor tissues were excluded. All of the RNA expression raw data including GSE10141, GSE14520 and TCGA-LIHC dataset were transformed to the format of Log2(count + 1) and normalized for subsequent research.

2.2 Bioinformatics Analysis

2.2.1 Estimation of Tumor Infiltrated Immune Cells

Based on Tumor IMmune Estimation Resource (TIMER) 2.0, we estimated the proportions of tumoral infiltrated immune cells in GSE10141, GSE14520 dataset using six algorithms, including XCELL, CIBERSORT, TIMER, QUANTISEQ, MCP-counter and EPIC algorithms.15 The count of infiltrated immune cells, namely, T cell, B cells, macrophages/monocyte, neutrophil, natural killer cell and dendritic cells, were estimated by MCP-counter.16 The median was set to divide samples into high- and low- infiltration groups. The “survival” R package was applied to exhibited the prognosis of high- and low- infiltration groups.

2.2.2 Analysis of single-Cell RNA Sequencing Data

The “Seurat” and ”SingleR” R package were used to perform the scRNA-seq data analysis. We controlled the quality of scRNA-seq data by excluding genes detected in < 3 cells, cells with < 200 total mapped genes, and 5% or more of unique molecular identifiers mapped in mitochondrial genes.. Afterwards, the top 2000 highly variable genes were screened for following analysis by the “FindVariableFeatures” function. Based on these 2000 variable genes, principal component analysis (PCA) were performed to visualize transcriptional variability of the scRNA-seq. T-distributed Stochastic Neighbour Embedding (t-SNE) algorithm was performed for further dimensional reduction with the primal principal components. Uniform manifold approximation and projection (UMAP) was also applied for non-linear dimensional reduction. According to the composition patterns of the feature markers, we annotated different cell clusters through “singleR” R package and then manually verified and proofread in the CellMarker database. The cell cluster associated with prognosis in survival analysis were singled out. The “clusterProfiler”, “GO plot”, “org.Hs.eg.db”, “DOSE” and “enrichplot” R package were used to conduct Gene Ontology (GO) and Kyoko Encyclopaedia of Genes and Genomes (KEGG) pathway enrichments analysis.

2.2.3 Construction and validation of prognosis risk model

The GSE10141 and GSE14520 dataset were combined to form a training set including 301 HCC tissues samples. The validation set consisted of 343 HCC cases from TCGA-LIHC dataset. The marker genes expression data singled out from training set were merged with survival time and status for ecah sample. We performed the least absolute shrinkage and selection operator (LASSO) regression by “glmnet” R package to obtain the most significant genes account for overall survival. The optimal lambda (λ) parameter was estimated with 10-fold cross-validation approach to avoid over-fitting. The significant genes in LASSO regression were included to construct multivariate Cox regression model. The optimum cut-off value was defined by X-tile software (version 3.6.1). We constructed the Kaplan–Meier survival curves to evaluate prognosis differences between high-risk and low-risk groups by ”Survival” and “Survminer” R packages both in training and validation set. The capacity of the risk score to predict the prognosis was assessed by adjusting for clinical variables in multivariate Cox analysis in validation set. Furthermore, receiver operating characteristic (ROC) curve, under the curve (AUC) and decision curve analysis (DCA) and predictiveness curves were implemented to the discrimination power of risk score. The characteristics of distribution of the risk score in different classical clinical variables groups were revealled by “ggstatsplot” R package.

2.2.4 Correlation between T cell infiltration and gene expression

We applied TIMER 2.0 to analysis the correlation between T cell infiltration and the expression level of gene included in the risk model. The correlation between marker gene expressions and infiltrated T cell counts (MCP-counter) in validation set were exhibited by “psych” R backage. T cell mainly contained two different subtype CD8 + T cell and CD4 + T cell according to the cluster of differentiation. The infiltration degree of CD8 + T cell and CD4 + T cell were calculated via TIMER algorithms. Spearmans correlation coefficient and P value revealed the strength and significance of interaction, respectively.

2.2.5 Statistical analysis

All of the figure in this study were constructed by using R (version 4.0.3) software. Statistical significance was set at P value of less than 0.05.

3 Result

3.1 The survival analysis revealed the correlation between infiltrating T cells and the OS of patients with HCC

Over the past few decades, many algorithms have been developed to quantify or estimate the infiltration percentage of immune cells in tumor tissues, like transcriptome or immunohistochemical staining-based estimation. We performed six algorithms, XCELL, CIBERSORT, TIMER, quanTIseq, MCP-counter and EPIC algorithms, to quantify the infiltration level of immune cells.15,17,18 Then, we integrated the survival data of HCC patients in two GEO cohorts (GSE10141, GSE14520) with the infiltration level of immune cell. Six basic immune cell types (T cells, B cells, natural killer cells, dendritic cells, macrophages, and neutrophils) quantified by the MCP-counter algorithm were selected for screening overall survival-related components.16 The results showed that B cell infiltration was associated with worse OS only in GSE14520 cohorts, while T cell infiltration was related to prolonged OS for HCC in both GEO cohorts with marginal statistical significance. (Fig. 1)

3.2 Identification of T cell marker genes based on single-cell sequencing data

After performing the quality control specification sprocedures described in the Methods section, the top 2000 variable genes of the scRNA-seq data (GSE146115) were used in cells clustering. We identified 7 clusters after carrying out the tSNE algorithm for dimensionality reduction with primal PCs. Each clusters possessed a subset of differentially expressed marker genes distinguished from other 6 cluster (with the thresholds |log2 (fold change)| > 1.0 and adjusted p value < 0.05 ) (Fig. 2B) According to this subset of differentially expressed marker genes, we annotated the cell-type of 7 clusters using “singleR” R package and manually verified in the CellMarker database. (Fig. 2A) The previously reported cell-type-specific marker genes were illuminated in specific cell cluster, such as CD79A, DERL3, CD38, and MZB1 in B cells cluster, CD3D, CD3E, CD3G, CD8A, CD8B, and PTPRC in T cells cluster, CD163, CD14, and CSF1R in Macrophage. (Fig. 2C) Cells in cluster 3 expressing T cell-specific markers were classified as T cells, cluster 5 as macrophages, and cluster 6 as B cells. All of the significant marker genes of 7 clusters are presented in Supplementary Table. T cells cluster contained 511 T cell marker genes. Subsequently, the GO functional analysis suggested these marker genes primarily participated in the process of T cell activation, antigen receptor-mediated signaling pathway, leukocyte proliferation, mononuclear cell differentiation, leukocyte cell-cell adhesion, regulation of T cell activation, regulation of immune effector process, immune response-activating cell surface receptor signaling pathway, T cell receptor signaling pathway, and immune response-activating signal transduction. The KEGG pathway analysis demonstrated that they were mostly involved in the T cell receptor signaling pathway, PD-L1 expression and PD-1 checkpoint pathway in cancer, leishmaniasis, Th1 and Th2 cell differentiation, Th17 cell differentiation, allograft rejection, human T cell leukemia virus 1 infection, primary immunodeficiency, etc. (Fig. 2D)

3.3 Construction of prognostic T cell marker genes signature (TCMS)

The training set obtained by merging 80 HCC tissues samples from GSE10141 and 221 HCC tissues samples from GSE14520 was used to identify a T cell gene prognostic signature. The LASSO regression analysis of training set revealed that 14 T cell marker genes ( RPL10, PSME1, TBC1D4, HELZ, PEBP1, CD3E, GZMA, F12, SLC2A2, JAK3, MACF1, H3F3B, CES1, TTR ) were associated with OS. (Fig. 3A, 3B) These 14 marker genes regarded as candidate genes of TCMS were incorporated into the multivariate Cox regression to determine the hub genes and their relative regression coefficients. Ultimately, the TCMS model was established including 4 T cell marker genes (HELZ, GZMA, SLC2A2, JAK3). (Fig. 3C) We calculated the risk score of each HCC case by combining relative expression levels and relative regression coefficients of the T cell marker genes in the classifier as follows: TCMS risk score = HELZ expression * 0.794 + GZMA expression * (-0.472) + SLC2A2 expression * (-0.405) + JAK3 expression * (-0.472). Moreover, among the four T cell marker genes, HELZ had positive coefficients meanwhile GZMA, SLC2A2, and JAK3 had negative coefficients.

TCMS risk score was calculated for each HCC patients from the training cohort. (Fig. 3D) According to cut-off value of risk score (-3.25) generated by X-tile software, HCC patients were divide into high-risk group (n = 229) and low-risk group (n = 72). (Fig. 3E) Survival analysis indicated that patients in high-risk group had a significantly poor OS compared with low-risk group. (P < 0.001). (Fig. 3F) The ROC analysis demonstrated that the AUC of this classifier was 0.769 indicating TCMS possesses relatively high sensitivity and specificity to predicte prognosis in training set.(Fig. 4A) The DCA curve indicated that TCMS showed good potential for clinical application as shown in Fig. 4B.

3.4 Verification of the prognostic TCMS in TCGA cohort

We used TCGA-LIHC cohort (343 HCC patients) as the validation set to verify the predictive ability of the TCMS. The TCMS risk score and survival state of each patients in validation set were exhibitd in (Fig. 3G). Based on the optimal cut-off value generated by X-tile software, the patients were classified into high- and low-risk score subgroups. (Fig. 3H) Kaplan-Meier analysis demonstrated that OS in low-risk socre group was significantly better than those of high-risk score group in validation set. (Fig. 3I) AUC of the TCMS in validation set was 0.620. (Fig. 4A) DCA curves showed that TCMS might be a potential classifier for stable clinical applications in predicting HCC prognosis. (Fig. 4C). Furthermore, the TCMS and classical clinical factors were included in the univariate and multivariate Cox analysis. As shown in Table 1, risk score of TCMS (Hazard Ratio(HR) = 1.919, 95%CI: 1.280–2.878, P = 0.002) and TNM stage (Hazard Ratio(HR) = 1.808, 95%CI: 1.463–2.234, P < 0.001) were influential features of OS. Further multivariate Cox analysis suggested that the risk score of TCMS was an independent prognostic factor for HCC patients. (Table 2).

To investigate the correlation of TCMS risk score with classical clinical factor including age, gender, tumor grade and Tumor Node Metastasis (TNM) stage, we calculated TCMS risk scores distribution in patients from TCGA cohort stratified by each clinical risk variables. We found the TCMS risk score steadily increased from TNM stage I to TNM stage III, with a statistically significant increase in TNM stage III (P < 0.01) (Fig. 4G). The TCMS risk score of TNM stage IV showed a downward trend, probably due to the small sample size. No statistical significant difference of TCMS risk score were found in age, gender and tumor grade group. (Fig. 4D, 4E, 4F)

3.5 The correlation between T cell infiltration and expression in TCMS

To lay the foundation for further investigation of these four genes in tumor immunology, we analysed the the association between T cell infiltration and expression levels in the TCMS. As shown in Fig. 5, patients in high-risk score group showed a reduction infiltration of T cell and low-risk score group patients possessed high degree of T cell infiltration. In TCGA validation set, HELZ, GZMA and JAK3 were positive correlation with T cell infiltration, oppositely, SLC2A2 were negative correlation with T cell infiltration.

Further immune analysis were carried out in TIMER 2.0 database. The result of correlation analysis suggested that increased expression of HELZ, GZMA and JAK3 were significantly correlated with high degree of CD8 + T cell, CD4 + T cell infiltration. In contrast, with the increased expression of SLC2A2, the degree of CD8 + T cell, CD4 + T cell infiltration decreased. Among them, the association between GZMA and CD8 + T cell were strong (Rho = 0.571, P = 3.49e-31) and SLC2A2 were weakly correlated with CD8 + T cell and CD4 + T cell infiltration. The other correlation were shown in Fig. 6.

4 Discussion

The treatment modalities for several solid tumors have entered a new era of immunological therapy.19 But for HCC, a highly devastating malignancy with high mortality rates and poor prognosis, the efficacy of immunological treatment was not satisfactory yet.20 A stringent barrier in optimizing HCC treatment efficacy was its immunosuppressive microenvironment, which causes low response rates to targeted therapy and limited improvements in OS, including immune checkpoint inhibitors (ICIs) based immunotherapy.21,22 The comprehension of TME in HCC was the key to develop optimized treatment modalities. The development of scRNA-seq and transcriptomes sequencing techniques provided researchers with a potential opportunity to explore and reveal TME. The infiltration of specific immune cell types could be reflected via the expression levels of markers by using algorithms such MCP-counter and TIMER. The correlation between expression levels of specific immune cell type markers and HCC patient prognosis could be established with complete survival data in multiple HCC cohorts.

In this study, we first observed that the infiltration of T cell types were associated with patient

survival. Then, we tried to explore the marker genes related to T cell infiltration by analyzing the intratumoral heterogeneity in HCC using sc-RNA-seq profiles. The functions of T cell marker genes were enriched in T cell activation, regulation of immune effector process, Th1, Th2 and Th17 cell differentiation, etc. Subsequently, we applied LASSO cox regression to constructe TCMS for higher accuracy based on gene expression matrix and follow-up data in training set and validated the TCMS with TCGA-LIHC dataset. In validation set, TCMS were still serving as an independently prognostic risk factor of OS. Besides, immune analysis in TIMER 2.0 database further confirmed the strong correlation between merker genes and T cell infiltration, which made the result more convincing.

Over the years a greater emphasis of TME research were placed on immune cell infiltration, in particular, the correlation between T cell infiltration and tumor prognosis.23–25 A comprehensive transcriptomic analysis comparing data from the Genotype-Tissue Expression project and TCGA-LIHC defined the increased T cell markers expression and T cell infiltration associated with a good prognosis.26,27 Consistent with this, the real-word cohort study including 65 HCC tissue specimens from stage I to IV data indicated high infiltration of CD8 + T cells in HCC were generally associated with better prognoses, lower recurrence and a more prolonged disease free survival (RFS). 11,28 Our results also supported the view that T cell infiltration in HCC were is an excellent classifier for predicting prognoses. Specifically, CD8 + T cells could develop powful antitumoral immune responses via scereting perforin, granzyme eradicating cancer cells.29,30 But the inhibitory immune checkpoints generating the immunosuppressive tumor microenvironment substantially exhausted CD8 + cytotoxic T lymphocytes (CTLs) and attenuated antitumoral immune responses, e.g., programmed death-ligand 1 (PD-L1), cytotoxic T-lymphocyte antigen 4 (CTLA-4) and lymphocyte activating gene 3 protein (LAG-3).31–33 The immune checkpoint inhibitors (ICIs), one of the innovative immunotherapy therapeutic approach, selectively block the inhibitory immune checkpoints to repair antitumoral immune responses mediated by T cells.

Previous studies concentrated on studying HCC-TME at transcriptome level and developed the biomarkers of HCC prognosis.34,35 However, few of them explored the heterogeneity of HCC-TME and prognostic marker genes by analysing scRNA-seq. Here, we tried to explore prognostic T cells marker genes and further constructed a prognostic model in HCC (TCMS) using scRNA-seq data and 3 HCC transcriptome cohorts with bulk sequencing. The TCMS model might reflect T cell function in HCC-TME somehow on account of the specific T cell marker genes expression profile identified by analysing scRNA-seq data. Regrettably, further analysis of subtype in T cell cluster were limited by the small scRNA-seq sample and T cell counts and could not be performed. Nevertheless. TCMS model was still a powerful predictor of OS in both training and validation dataset. And the T cell marker genes we screened out provided a basis for further research of the role of T cell in TME and immunotherapy.

Four T cell marker genes (HELZ, GZMA, SLC2A2, JAK3) included in TCMS model indicated pivotal roles in T cells’ biological behaviours in HCC-TME. In correlation analysis, GZMA was

strongly associated with CD8 + T cell infiltration (Rho = 0.571, P = 3.49e-31). GZMA (Granzyme A) encoded a specific serine protease which may function as a necessary component for cytolytic T lymphocytes (CTL) to recognize and lyse specific target cells. The CTL-derived GZMA catalyzed the cleavage of gasdermin B, unleashed pore-forming activity and promoted pyroptosis.30 Therefore, GZMA may serve as an noteworthy component in antitumor immunity. Previous studies have indicated the correlation between GZMA and outcomes of HCC.36 In this study, we found high expression of GZMA predicted high degree of T cell infiltration and good prognosis, indicating GZMA functions anti-neoplastic roles in HCC. SLC2A2 was the member of solute carrier family 2 and encoded the integral plasma menbrance glycoprotein which was responsible for glucose uptake and mediated the bidirectional transfer of glucose across the hepatocytes plasma membrane. Cancer cells necessitated high glucose transport rates for satisfying the high metabolic demands, and because of this property cancer cells up-regulate the glucose transporters expressions.37 But the study of metabolic genes in HCC cell lines reported the expression of SLC2A2 were down-regulated, 38 and another research form Korea showed high SLC2A2 HCC patients had a good prognosis,37 which were consistent with the regression coefficients of SLC2A2 in TCMS model. In our scRNA-seq analysis, SLC2A2 was high expression in T cell cluster compared with other cell clusters and identified as the marker gene of T cell. However, in immuno-infiltration analyse based on the RNA sequencing data and TIMER database SLC2A2 were negative correlation with T cell infiltration (Fig. 5D, Fig. 6), which were not in conformity with the aforementioned results, even though the spearmans correlation coefficient were weak. Further researches were necessary to explain SLC2A2 in HCC-TME. JAK3 was primarily expressed in immune cells and mediated essential intracellular signal transduction via tyrosine phosphorylation activated by interleukin receptors. Non-receptor tyrosine kinase (JAK) involved in various biological processes like cell growth, development, or differentiation, particularly, played a crucial role during T-cells development. Previous studies have reported that mutations of JAK3 were identified in HCC patients and cell lines. 39 HELZ was the member of the superfamily I class of RNA helicases and implicated in various aspects of RNA metabolism. Nevertheless, functions of JAK3 in T cell were largely unknown and few research mentioned the relationship between HELZ and HCC patients.40 In this study, HELZ was identified a T cell marker gene and associated with HCC patients outcome, suggesting a novel role of HELZ.

Although our conclusions have been vcalidated by independent data set, we have not implemented futher prospective clinical study or molecular biology experiments to validate the results, which were the largest limitation in this study. So we would conduct bioinformatics research, subsequently, to confirm our conclusions and underlying molecular mechanism.

In conclusion, we successfully constructed the T cell marker genes signature with powerful predictive function and provided a new understanding of T cell infiltration in tumor microenvironment which might offer practices instructions for HCC immunotherapy.

Declarations

5  Author Contributions

Yaojian Shao and Chenyu Jiang: study design and manuscript drafting; Jiangmin Liang, Yiwan Wang, Chaohui Wang: revision of the manuscript; Shenyu Wei and Junjie Ni: statistical analysis. All authors reviewed the manuscript.

6  Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

7  Data Availability Statement

The datasets (GSE146115, GSE10141, GSE14520) for this study can be found in the Gene Expression Omnibus (GEO) database. GSE146115 can be accessed in the website: https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE146115; GSE10141 can be accessed in website: https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE10141;   and GSE14520 can be accessed in the website https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE14520. 

The TCGA-LIHC dataset can be accessed in The Cancer Genome Atlas (TCGA) database (https://portal.gdc.cancer.gov/).

8  Acknowledgements

Not applicable.

Ethics approval and Consent to participate

The data of this study is from publicly available database, and all methods were carried out in accordance with relevant guidelines and regulations. This study was carried out in compliance with the MIAME compliant.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

References

  1. Sung H, Ferlay J, Siegel RL. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. May 2021;71(3):209–249.
  2. Villanueva A. Hepatocellular Carcinoma. The New England journal of medicine. Apr 11 2019;380(15):1450–1462.
  3. Allemani C, Matsuda T, Di Carlo V, et al. Global surveillance of trends in cancer survival 2000-14 (CONCORD-3): analysis of individual records for 37 513 025 patients diagnosed with one of 18 cancers from 322 population-based registries in 71 countries. Lancet (London, England). Mar 17 2018;391(10125):1023–1075.
  4. Egen JG, Ouyang W, Wu LC. Human Anti-tumor Immunity: Insights from Immunotherapy Clinical Trials. Immunity. Jan 14 2020;52(1):36–54.
  5. Sangro B, Sarobe P. Advances in immunotherapy for hepatocellular carcinoma. Aug 2021;18(8):525–543.
  6. Galle PR, Finn RS, Qin S, et al. Patient-reported outcomes with atezolizumab plus bevacizumab versus sorafenib in patients with unresectable hepatocellular carcinoma (IMbrave150): an open-label, randomised, phase 3 trial. The Lancet Oncology. Jul 2021;22(7):991–1001.
  7. Rizzo A, Ricci AD, Brandi G. Immune-based combinations for advanced hepatocellular carcinoma: shaping the direction of first-line therapy. Future oncology (London, England). Mar 2021;17(7):755–757.
  8. Bian J, Lin J, Long J, et al. T lymphocytes in hepatocellular carcinoma immune microenvironment: insights into human immunology and immunotherapy. American journal of cancer research. 2020;10(12):4585–4606.
  9. Sia D, Jiao Y, Martinez-Quetglas I, et al. Identification of an Immune-specific Class of Hepatocellular Carcinoma, Based on Molecular Features. Gastroenterology. Sep 2017;153(3):812–826.
  10. Garnelo M, Tan A, Her Z, et al. Interaction between tumour-infiltrating B cells and T cells controls the progression of hepatocellular carcinoma. Gut. Feb 2017;66(2):342–351.
  11. Gabrielson A, Wu Y, Wang H, et al. Intratumoral CD3 and CD8 T-cell Densities Associated with Relapse-Free Survival in HCC. Cancer immunology research. May 2016;4(5):419–30.
  12. Ma J, Zheng B, Goswami S, et al. PD1(Hi) CD8(+) T cells correlate with exhausted signature and poor clinical outcome in hepatocellular carcinoma. Journal for immunotherapy of cancer. Nov 29 2019;7(1):331.
  13. Wang T, Dang N. Integrating bulk and single-cell RNA sequencing reveals cellular heterogeneity and immune infiltration in hepatocellular carcinoma. Feb 6 2022.
  14. Heinrich B, Gertz EM, Schäffer AA, Craig A, Ruf B. The tumour microenvironment shapes innate lymphoid cells in patients with hepatocellular carcinoma. Aug 2 2021.
  15. Li T, Fan J, Wang B, et al. TIMER: A Web Server for Comprehensive Analysis of Tumor-Infiltrating Immune Cells. Cancer research. Nov 1 2017;77(21):e108-e110.
  16. Becht E, Giraldo NA, Lacroix L, et al. Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome biology. Oct 20 2016;17(1):218.
  17. Aran D, Hu Z, Butte AJ. xCell: digitally portraying the tissue cellular heterogeneity landscape. Nov 15 2017;18(1):220.
  18. Chen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling Tumor Infiltrating Immune Cells with CIBERSORT. Methods in molecular biology (Clifton, NJ). 2018;1711:243–259.
  19. Kole C, Charalampakis N, Tsakatikas S, et al. Immunotherapy for Hepatocellular Carcinoma: A 2021 Update. Cancers. Oct 4 2020;12(10)doi:10.3390/cancers12102859
  20. Xing R, Gao J, Cui Q, Wang Q. Strategies to Improve the Antitumor Effect of Immunotherapy for Hepatocellular Carcinoma. Frontiers in immunology. 2021;12:783236.
  21. Leone P, Solimando AG. The Evolving Role of Immune Checkpoint Inhibitors in Hepatocellular Carcinoma Treatment. May 20 2021;9(5).
  22. Zhou G, Sprengers D, Boor PPC, et al. Antibodies Against Immune Checkpoint Molecules Restore Functions of Tumor-Infiltrating T Cells in Hepatocellular Carcinomas. Gastroenterology. Oct 2017;153(4):1107–1119.e10.
  23. Kumar-Sinha C, Sahai V. T-Cell Subsets as Potential Biomarkers for Hepatobiliary Cancers and Selection of Immunotherapy Regimens as a Treatment Strategy. Journal of the National Comprehensive Cancer Network: JNCCN. Feb 2022;20(2):203–214.
  24. Ott PA, Bang YJ, Piha-Paul SA, et al. T-Cell-Inflamed Gene-Expression Profile, Programmed Death Ligand 1 Expression, and Tumor Mutational Burden Predict Efficacy in Patients Treated With Pembrolizumab Across 20 Cancers: KEYNOTE-028. Journal of clinical oncology: official journal of the American Society of Clinical Oncology. Feb 1 2019;37(4):318–327.
  25. Zheng C, Zheng L, Yoo JK, et al. Landscape of Infiltrating T Cells in Liver Cancer Revealed by Single-Cell Sequencing. Cell. Jun 15 2017;169(7):1342–1356.e16.
  26. Foerster F, Hess M, Gerhold-Ay A, et al. The immune contexture of hepatocellular carcinoma predicts clinical outcome. Scientific reports. Mar 29 2018;8(1):5351.
  27. Ostroumov D, Fekete-Drimusz N, Saborowski M, Kühnel F, Woller N. CD4 and CD8 T lymphocyte interplay in controlling tumor growth. Feb 2018;75(4):689–713.
  28. Pfister D, Núñez NG. NASH limits anti-tumour surveillance in immunotherapy-treated HCC. Apr 2021;592(7854):450–456.
  29. Oura K, Morishita A. Tumor Immune Microenvironment and Immunosuppressive Therapy in Hepatocellular Carcinoma: A Review. May 28 2021;22(11).
  30. Zhou Z, He H. Granzyme A from cytotoxic lymphocytes cleaves GSDMB to trigger pyroptosis in target cells. May 29 2020;368(6494).
  31. Kim HD, Song GW, Park S, et al. Association Between Expression Level of PD1 by Tumor-Infiltrating CD8(+) T Cells and Features of Hepatocellular Carcinoma. Gastroenterology. Dec 2018;155(6):1936–1950.e17.
  32. Han Y, Chen Z, Yang Y, et al. Human CD14 + CTLA-4 + regulatory dendritic cells suppress T-cell response by cytotoxic T-lymphocyte antigen-4-dependent IL-10 and indoleamine-2,3-dioxygenase production in hepatocellular carcinoma. Hepatology (Baltimore, Md). Feb 2014;59(2):567–79.
  33. Ruffo E, Wu RC, Bruno TC, Workman CJ, Vignali DAA. Lymphocyte-activation gene 3 (LAG3): The next immune checkpoint receptor. Seminars in immunology. Apr 2019;42:101305.
  34. Zhu Y, Shan D, Guo L, Chen S. Immune-Related lncRNA Pairs Clinical Prognosis Model Construction for Hepatocellular Carcinoma. 2022;15:1919–1931.
  35. Huang S, Zhang J, Lai X, Zhuang L, Wu J. Identification of Novel Tumor Microenvironment-Related Long Noncoding RNAs to Determine the Prognosis and Response to Immunotherapy of Hepatocellular Carcinoma Patients. Frontiers in molecular biosciences. 2021;8:781307.
  36. Mauriello A, Zeuli R, Cavalluzzo B, et al. High Somatic Mutation and Neoantigen Burden Do Not Correlate with Decreased Progression-Free Survival in HCC Patients not Undergoing Immunotherapy. Cancers. Nov 20 2019;11(12).
  37. Kim YH, Jeong DC, Pak K, et al. SLC2A2 (GLUT2) as a novel prognostic factor for hepatocellular carcinoma. Oncotarget. Sep 15 2017;8(40):68381–68392.
  38. Nwosu ZC, Battello N, Rothley M, et al. Liver cancer cell lines distinctly mimic the metabolic gene expression pattern of the corresponding human tumours. Journal of experimental & clinical cancer research: CR. Sep 3 2018;37(1):211.
  39. Lu J, Yin J, Dong R, et al. Targeted sequencing of cancer-associated genes in hepatocellular carcinoma using next generation sequencing. Molecular medicine reports. Sep 2015;12(3):4678–4682.
  40. Hanet A, Räsch F. HELZ directly interacts with CCR4-NOT and causes decay of bound mRNAs. Oct 2019;2(5).

Tables

Table 1

Multivariate analysis of hepatocellular carcinoma in training set.

Variables

Coef

HR

HR.95L

HR.95H

p Value

HELZ

0.794 

2.213 

1.318 

3.715 

0.003 

GZMA

-0.472 

0.624 

0.381 

1.021 

0.061 

SLC2A2

-0.405 

0.667 

0.426 

1.044 

0.077 

JAK3

-0.472 

0.624 

0.441 

0.881 

0.007 

Coef : hepatocellular carcinoma, HR : Hazard ratio.

 

Table 2

Univariate and multivariate analysis of hepatocellular carcinoma in validation set

Variables

Univariate analysis

　

Multivariate analysis

HR (95% CI)

p Value

　

HR (95% CI)

p Value

Age (>60 /≤60)

1.172(0.800-1.717)

0.414 

　

1.137(0.768-1.682)

0.522 

Gender (female /male)

0.758(0.513-1.118)

0.162 

0.908(0.608-1.356)

0.637 

Grade (G1 /G2 /G3 /G4)

1.121(0.865-1.454)

0.388 

1.170(0.888-1.542)

0.263 

TNM Stage (I /II /III /IV)

1.808(1.463-2.234)

<.001

1.771(1.423-2.203)

<.001*

Risk score (high /low)

1.919(1.280-2.878)

0.002 

　

1.664(1.104-2.507)

0.015* 

Coef : hepatocellular carcinoma, HR : Hazard ratio.

p values < 0.05 were marked with * and highlighted in bold.