- 1Department of Gastroenterology, Daping Hospital, Army Medical University, Chongqing, China
- 2Department of General Surgery, The Southern Theater Air Force Hospital, Guangzhou, China
- 3Department of Neurology, Xijing Hospital, Fourth Military Medical University, Xi’an, China
- 4Department of General Surgery, Tangdu Hospital, Fourth Military Medical University, Xi’an, China
- 5Department of Surgery, Tangdu Hospital, Fourth Military Medical University, Xi’an, China
- 6Department of Cell Biology and Genetics, School of Basic Medical Sciences, Xi’an Jiaotong University Health Science Center, Xi’an, China
- 7Centre for Computational Systems Medicine, School of Biomedical Informatics, The University of Texas Health Science Centre at Houston, Houston, TX, United States
Necroptosis is a programmed form of necrotic cell death in regulating cancer ontogenesis, progression, and tumor microenvironment (TME) and could drive tumor-infiltrating cells to release pro-inflammatory cytokines, incurring strong immune responses. Nowadays, there are few identified biomarkers applied in clinical immunotherapy, and it is increasingly recognized that high levels of tumor necroptosis could enhance the response to immunotherapy. However, comprehensive characterization of necroptosis associated with TME and immunotherapy in Hepatocellular carcinoma (HCC) remains unexplored. Here, we computationally characterized necroptosis landscape in HCC samples from TCGA and ICGA cohorts and stratified them into two necroptosis clusters (A or B) with significantly different characteristics in clinical prognosis, immune cell function, and TME-landscapes. Additionally, to further evaluate the necroptosis levels of each sample, we established a novel necroptosis-related gene score (NRGscore). We further investigated the TME, tumor mutational burden (TMB), clinical response to immunotherapy, and chemotherapeutic drug sensitivity of HCC subgroups stratified by the necroptosis landscapes. The NRGscore is robust and highly predictive of HCC clinical outcomes. Further analysis indicated that the high NRGscore group resembles the immune-inflamed phenotype while the low score group is analogous to the immune-exclusion or metabolism phenotype. Additionally, the high NRGscore group is more sensitive to immune checkpoint blockade-based immunotherapy, which was further validated using an external HCC cohort, metastatic melanoma cohort, and advanced urothelial cancer cohort. Besides, the NRGscore was demonstrated as a potential biomarker for chemotherapy, wherein the high NRGscore patients with more tumor stem cell composition could be more sensitive to Cisplatin, Doxorubicin, Paclitaxel-based chemotherapy, and Sorafenib therapy. Collectively, a comprehensive characterization of the necroptosis in HCC suggested its implications for predicting immune infiltration and response to immunotherapy of HCC, providing promising strategies for treatment.
Introduction
Hepatocellular carcinoma (HCC) is acknowledged to be one of the most common malignant tumors globally, accounting for one-third of cancer mortalities (1). Risk factors of HCC progression contain metabolic disorders, viral infections by hepatitis B virus and hepatitis C virus, absorption of the aflatoxin-contaminated meal, and cirrhosis out of alcoholic hepatitis (2–4). Though early examination and intervention have achieved success, the ratio of diagnosis at the early stage remains low (5). The majority of HCC patients are diagnosed at an advanced stage (6). Although plentiful efforts were made in both diagnosis and treatment of HCC, its overall survival (OS) rate remains unfavorable. Seventy percent of HCC patients have recurrent neoplasm 5 years after resection (7). Consequently, the 5-year survival rate of HCC patients turns out unsatisfactory (8).
Chemotherapy failure has been a major obstacle during cancer treatment, among which apoptosis resistance (innate or acquired) is widely accepted mechanistically. How to bypass the apoptosis pathway and induce effective cell death pathways are becoming crucial in the treatments of HCC. Like apoptosis and necrosis, necroptosis is a novel programmed form of necrotic cell death that is mainly manipulated by receptor-interacting Protein Kinase 1 (RIPK1), RIPK3, and Mixed Lineage Kinase domain-like (MLKL) (9). Accumulating data demonstrate that necroptosis plays a crucial role in the regulation of cancer biology, including cancer progression (10, 11), cancer metastasis (12), cancer Immunosurveillance (13), and cancer subtypes (14, 15).
With the advent of immunotherapy, such as immune checkpoint blockades (ICBs), only a few HCC patients have been reported to gain clinical benefits from ICBs (16, 17). To date, Some predictors of response to ICBs were identified, including PD-L1 expression (18), the degree of cytotoxic T cell infiltration (19), tumor mutational burden (TMB) (20), mismatch repair deficiency (21), and activated Wnt/β-catenin signaling (22), et al. However, to our best knowledge, these above biomarkers have been confirmed in HCC (23). There is a rationale supporting the development of some novel biomarkers for ICBs in HCC patients. Recently, a few studies revealed crosstalk between necroptosis and antitumor immunity (24). Necroptosis not only has direct interaction with immune cells like dendritic cells (DCs) and natural killer T cells (NKT) but also could initiate an adaptive immune system by promoting DC cells and macrophages to release pro-inflammatory cytokines into the tumor microenvironment (TME), incurring strong immune responses (25). In in vivo and in vitro experiments, necroptotic tumor cells were shown to induce antitumor immunogenicity through the cross-priming and proliferation of CD8+ T cells (26). Furthermore, with the advance in nanomedicine, necroptotic cancer cells could boost antitumor immunity by administrating cell-mimicry nanovaccine with a tailored immunostimulatory modality (26). In addition, the effect of tumor regression by nanovaccine could be enhanced by combination with ICBs (26). However, the discovery of a specific necroptosis marker, a thorough investigation of the molecular mechanism, and a clarification of its crosstalk with other cell death machinery and its interaction with the immune system should be urgently further investigated. Taken together, necroptosis has a close relationship with TME and antitumor immunity, suggesting necroptosis could become a novel biomarker for chemotherapy and immunotherapy of HCC. However, it’s rarely depicted for necroptosis characteristics from a multi-omics perspective and its correlation with ICBs and chemotherapy in HCC.
In this study, we uncovered that necroptosis regulators were clinically predictive and independent factors for HCC patients by survival analysis and unsupervised clustering analysis. Necroptosis regulators could distinguish patients into two necroptosis clusters (Nclusters A and B) and different Nclusters were enriched into several TME-related pathways and various immune infiltration pathways. Therefore, we proposed a hypothesis: necroptosis could become a new and indispensable factor for HCC prognosis and TME characterization. A novel necroptosis-related gene score (NRGscore) was established by the Lasso algorithm and multiple Cox regression analysis. In TME cell infiltration of HCC, NRGscore had a negative relation with activated NK cells and macrophage M1 cells while positively associated with Treg cells. NRGscore had meaningful guidance for patients, where the high NRGscore group could be more sensitive to ICB treatment, which was validated in other independent cancer cohorts, as well. As for chemotherapy and targeted therapy, due to more stem cells and enhanced cell proliferation, the high NRGscore group could be more sensitive to Cisplatin, Doxorubicin, Paclitaxel-based chemotherapy, and Sorafenib targeted therapy than its counterparts. Taken together, the prognostic NRGscore system could help to dissect the TME characterization of HCC and to interpret the clinical responses to chemotherapies and immunotherapies, providing new target molecules for the treatment of cancers.
Methods
Data Sources
The workflow of this study is shown in Figure S1. Necroptosis-related genes were collected based on the GESA necroptosis gene set and published literature (Table S1). Clinical information and mRNA expression matrixes of HCC patients were acquired from The Cancer Genome Atlas (TCGA, https://portal.gdc.cancer.gov/ ) database and International Cancer Genome Consortium (ICGC, https://dcc.icgc.org/) database. Then we transformed FPKM values into transcripts per kilobase million (TPM) values. The basic clinical information of HCC datasets in the study is summarized in Table S2. The somatic mutation data of the TCGA and ICGC cohort were downloaded from the UCSC Xena (https://gdc.xenahubs.net/). Additionally, the TCGA cohort was applied to Copy Number Variation (CNV) analysis. Another HCC validated cohort (GSE54236) (27, 28) and corresponding clinical features were collected to assess predictive power in our study.
Identification of Significant Mutational Genes
The “maftools” R package was applied to process the mutation annotation format (maf) data, and the “MutSigCV” algorithm was implemented to screen the significant mutational genes (SMGs) (29). The significance of nonsilent somatic mutations in a gene was measured based on the background mutation rates by silent mutation. The false discovery rates (FDR (30)) were then calculated, and genes with statistical significance (FDR ≤ 0.1) were set as SMGs. Then, waterfall plots were employed to visualize the mutation information of these significant SMGs in the TCGA cohort. Besides, Fisher’s test was applied to detect the mutually exclusive or co-occurring ratio of necroptosis-related genes. By adopting the “ExtractSignatures” function that applies the Bayesian nonnegative matrix factorization-based framework, we determined the mutational signatures using the genomic data. The optimal number of mutational signatures for the TCGA cohort could be detected by the “SignatureEnrichment” function and then it automatically assigned a given signature to each sample.
Unsupervised Cluster Analysis
Based on the expression of necroptosis-related genes, unsupervised clustering analysis was performed to stratify patients into different clusters. The “ConsensuClusterPlus” R package was adopted to determine the number of clusters and guarantee the stability of classification (31). The principal component analysis (PCA) was utilized to investigate gene-expression arrays among distinct clusters.
Tumor Infiltration Cell and Immune-Related Function Analysis
To estimate the abundance and activity of tumor infiltration cells (TICs) in HCC, the single-sample gene set enrichment analysis (ssGSEA) algorithm was implemented by using the “GSVA” R package (32). In addition, CIBERSORT (33), an analysis algorithm based on the immune gene set, was also used to evaluate the TICs levels of HCC. The algorithm was run for 1000 permutations and HCC samples with an output P < 0.05 were selected as previously reported (34, 35). Twenty-three types of TICs were comprised of adaptive immune cells (B cells, T cells, CD8 T cells, T follicular helper (Tfh), Th1, Th2, Th17, and Treg cells), and innate immune cells (NK cells, CD56 dim NK cells, CD56 bright NK cells, DCs, plasmacytoid DCs, immature DCs, neutrophils, mast cells, and macrophages). Besides, immune-related pathways (such as cytolytic activity, T-cell costimulation, inflammation-promoting, and para-inflammation) were also calculated via ssGSEA. The biosimilarity of infiltrating immune cells and immune-related functions were estimated by the Gaussian fitting model.
Functional Annotation Analysis
To investigate the enrichments of biological processes, cell components, and molecular function pathways, Gene Ontology (GO) analysis, and Kyoto Encyclopedia of Genes and Genomes (KEGG) signaling pathway analysis were conducted using the “clusterProfiler” R package. Additionally, Gene Set Enrichment Analysis (GSEA) (36) was adopted to identify the differences in DEGs between the distinct clusters in the enrichment of the KEGG pathway. Permutations were performed 1000 times for each analysis. P-value < 0.05 and adjusted P-value (Q value) < 0.05 were considered statistically significant.
Differentially Expressed Genes Analysis
To identify necroptosis-related genes, patients were classified into distinct groups according to sample types, necroptosis clusters, and NRGscore, respectively. The “limma” R package was utilized to determine DEGs between different groups (37). The significance filtering criteria for determining DEGs were set as adjusted P-value < 0.001 and fold change > 1.5.
Dimension Reduction and Construction of NRGscore
Prognostic genes were identified from DEGs by performing the univariate Cox regression. The least absolute shrinkage and selection operator (Lasso) regression was conducted to necroptosis gene signature via utilizing the “glmnet” R package. Responding coefficients (β) of the signature were verified. Additionally, the signature was calculated by the following equation: NRGscore = ∑(exp(gene)∗β), where exp indicated RNA expression of HCC samples.
Clinical Characteristic Evaluation of NRGscore
The survival curves for the different subgroups were generated by Kaplan-Meier (K-M) methods. The areas under the curve (AUC) of the receiver operating characteristic (ROC) curve were applied to assess the predictive value of gene signature. The nomogram was built based on the NRGscore and clinicopathologic characteristics including age, gender, T stage, and N stage to predict the survival probability of 1-, 3-, and 5-year OS of HCC patients. The calibration curve of the nomogram was plotted to estimate the prediction possibilities according to the observed survival rates. The nomogram and calibration plots were generated based on the “rms” R package.
Estimation of Tumor Immune Microenvironment and Tumor Mutational Burden
To further dissect the immune landscape of HCC, the “ESTIMATE” package in R was used to evaluate the immune, stromal, and ESTIMATE scores, which reflect the ratio of the immune/stromal components of the tumor immune microenvironment (TME). To determine the tumor mutational burden (TMB) of each patient, we also counted the nonsynonymous and synonymous mutation counts in the TCGA cohort.
Assessment of Clinical Response to the Immunotherapy
TIDE (38) algorithm and immunophenoscore (IPS) function (39) were prevalently recognized to be effective methods to predict cancer patients’ response to immunotherapy. There are two main mechanisms to immune escape and resistance to ICB-based immunotherapy (40, 41): (1) high levels of dysfunctional CTL; (2) immunosuppressive factors to exclude T cells from the tumor region. TIDE algorithm integrates these two immune-escape mechanisms and was used to predict patients’ response to ICB-based therapy based on the transcriptome. According to characterizing the determinant factors of cancer immunogenicity and antigenomes, we stratified HCC patients into different IPS groups. HCC sample is more immunogenic when the z-score of IPS is higher (42). Two immunotherapy cohorts, metastatic melanoma cohort treated with Nivolumab (19), and advanced urothelial cancer cohort with the intervention of atezolizumab (IMvigor210, http://research-pub.gene.com/IMvigor210CoreBiologies/) (43) were downloaded to evaluate the predictive value of the necroptosis score system for immunotherapy. In addition, the TIDE website (http://tide.dfci.harvard.edu/login/) was used to further evaluate the predictive power of other cancer cohorts by inputting gene coefficients. The raw gene expression data of all cohorts were normalized according to previous literature (44, 45).
Dissection of Cancer Stem Cell And Prediction of Chemotherapeutic Drugs Sensitivity
To evaluate the cancer stem cell ratio and differentiation degree of HCC, the “limma” and “corrplot” R packages were applied. We adopted Spearman’s method to explore the correlation between NRGscore and cancer stem cells. In addition, based on the information retrieved from the Genomics of Drug Sensitivity in Cancer database (https://www.cancerrxgene.org/ ), we estimate the sensitivity of different chemotherapeutic drugs between high and low NRGscore subgroups. The prediction process used was the “pRRophetic” package (46) where the half-maximal inhibitory concentration was calculated by ridge regression model based on gene expression profiles.
Statistical Analysis
All statistical analyses were conducted by R software (http://www.R-project.org, version 4.1.2), PERL programming language (version 5.32.1.1, https://www.perl.org/). K-M curve analysis with a log-rank test was utilized to compare OS between diverse subgroups. Mann-Whitney and Kruskal-Wallis test with adjusted P-values were employed to compare either ssGSEA scores of immune cells or functions of the distinct clusters as indicated in the article. P-value < 0.05 was statistically significant.
Results
Hepatocellular Carcinoma Was Well Characterized by Necroptosis Genes
Firstly, we investigated 67 necroptosis genes (Table S1) regarding their expression and genetic variations between tumor and paratumor samples in TCGA cohorts. Fifty-one necroptosis genes presented significantly differential expression, of which 42 genes were upregulated and nine genes were downregulated in tumors (Figure 1A). Many differentially expressed necroptosis genes have been reported to be involved in modifying the microenvironment of HCC (15), including upregulated SIRT1/2 (47) and RIPK1 (48). It is well known that genetic variations, such as copy number variations (CNV) and mutations, result in perturbations of gene expression during tumorigenesis of HCC. The investigation of CNV alternation frequency showed a prevalent change in 67 necroptosis genes, and we observed that 44 genes exhibited the amplification in copy numbers and 23 genes showed depletion of CNVs (Figures 1B, C). Additionally, a total of 108 of 364 (29.67%) patients exhibited at least one type of mutation (Figure 1D). The most highly mutated genes, such as MYC, TERT, FASLG, RIPK1, TNFSF10, TARDBP, and CDKN2A, have been well characterized in liver cancer (15, 48), diffuse large B-cell lymphoma (49), neuroblastoma (50), and kidney renal clear cell carcinoma (51, 52). Especially, as the three most important necroptosis-driving molecules (53), RIPK1 and MLKL were upregulated in tumors and CNV of RIPK1 revealed amplification more than depletion, while RIPK3 presented no significant alternation in genomic expression. Extrapolating from all the above results, necroptosis genes present high heterogeneity among RNA expression, CNV, and mutations in HCC samples, which shows promise in tumorigenesis and development in HCC.
Figure 1 Genetic landscape and expression variation of necroptosis regulators in HCC. (A) The differential expression level of necroptosis regulators between tumor and normal tissues. (B) CNV alternation of necroptosis regulators in tumor tissues. Column represented the frequency of the variations. The green dot represented the deletion of CNV. The red dot represented the amplification of CNV. (C) Location of CNV of necroptosis regulators in chromosomes. Red dots represent genes gain than loss, blue dot presents genes loss than gain and the black dot means loss equal gain. (D) One hundred eight patients (29.67%) exhibited various genetic alterations, including missense, nonsense, splice, frameshift, and multiple mutations. Each column indicated individual HCC patients, and the upper bar diagram exhibited the TMB of HCC patients. The right number represented the mutation frequency, and the bar diagram on the right exhibited the proportion of each genetic alteration. P < 0.05 *; P < 0.01**; P < 0.001***. HCC, hepatocellular carcinoma; CNV, copy number variation; TMB, tumor mutation burden.
Necroptosis Clusters Exhibit Distinct Prognosis, Biological Functions, and Immune Characteristics
To integrate the prognostic value of these necroptosis genes, we sought to construct a prognosis and correlation network according to the clinical outcomes of HCC patients (Figure 2A). Three necroptosis genes were favorable factors while 22 genes were risk factors, which disclosed that necroptosis may play a tumor-promoting role in HCC. Interestingly, we found that a positive correlation between prognostic necroptotic genes occurred more frequently than negative connections. To further explore the clinical role of necroptosis genes in HCC, we collected 599 clinical samples from TCGA and ICGC datasets and the basic clinical parameters are shown in Table S2. Based on expression profiles of necroptosis genes, we implemented unsupervised clustering to analyze the HCC samples from two cohorts and classified patients into qualitatively different subgroups. Two distinct subgroups were ultimately identified, including 244 cases in subgroup A and 355 cases in subgroup B, which were referred to as necroptosis clusters (Ncluster) A and B, respectively (Figure 2B). Additionally, PCA results showed necroptosis genes could vividly distinguish one subgroup from another (Figure 2C). K-M curves for the two subgroups revealed the prominent survival advantage in Ncluster B compared to A (Figure 2D). Moreover, correlations between two clusters and other clinicopathological parameters, including age, gender, TNM stage, tumor stage, and tumor grade, are shown in the heatmap (Figure S2). Ncluster A was characterized by higher levels of necroptosis and more advanced clinical stage and grade, thus resulting in a poor prognosis.
Figure 2 Unsupervised clustering of necroptosis regulators. (A) Interaction network of necroptosis regulators in HCC. Necroptosis regulators were indicated by red circles. Risk factors were indicated by purple circles, and favorable factors were indicated by green circles. Different sizes of the circle represented the different P-values. Red lines showed the positive correlation, and blue lines showed the negative correlation. (B) Consensus clustering matrix for k =2. (C) The feature distribution between different two Nclusters was plotted via PCA. (D) Kaplan-Meier curves of the OS for the two Nclusters of HCC patients. Ncluster, necroptosis cluster; PCA, principal component analysis; OS, overall survival; HCC, hepatocellular carcinoma.
To further explore the potential biological behavior of necroptotic patterns, we performed GSVA analysis between two Nclusters. In KEGG signal pathways, Ncluster B was characterized by enhanced metabolism pathways, such as drugs, steroid hormone, bile acid, and fatty acid metabolism, while Ncluster A was characterized by upregulated cell endocytosis and phagocytosis pathways (Figure 3A). The tumor-infiltrating cells and immune functions play a critical role in the process of tumorigenesis and migration in HCC (54). In immune pathways of GSVA, we found Ncluster A possess more pathways for activated immune cell than Ncluster B (Figure S3A), suggesting N cluster A had a trend to immune-inflamed phenotype and Ncluster B was inclined to immune-excluded phenotype. Next, we investigated the difference between TICs and immune-related functions in two clusters by ssGSEA analysis. As vividly shown in Figure 3B, 10 immune functions were markedly downregulated and only IFN-response function upregulated in Ncluster B, which revealed that Ncluster B had less complicated TME. Consistent with results of immune pathways, there was significantly different TIC abundance in two subgroups, among which 19/23 TICs varied immensely, including 17 highly expressed types and only two downregulated kinds of TICs in Ncluster A (Figure 3C). In addition, we also determined 2988 Ncluster-related DEGs (Table S3) by limma R packages and conducted GO and KEGG enrichment analysis based on the Ncluster DEGs. Likely, DEGs were enriched in cell adhesion (EMT pathways) and T cell activation pathways (Figure S3B), PI3K-AKT pathways, and abundant immune cell regulation pathways (Figure S3C). Collectively, necroptosis could classify patients into two distinct subgroups and Ncluster A tends to have immune-inflamed phenotype and complex tumor-infiltrating patterns while Ncluster B was inclined to immune-excluded phenotype and metabolism phenotype.
Figure 3 Dramatic difference of biological features between Ncluster A and B. (A) Heatmap showed various KEGG pathways were enriched in Ncluster A and B. (B) Immune function and infiltrating immune cells (C) altered in two Nclusters shown by violin diagram. The heatmap was used to visualize these biological processes, and red represented activated pathways and blue represented inhibited pathways. Ncluster, necroptosis cluster.
Necroptosis Patterns Define Necroptosis-Related Gene Clusters and NRGscore
To further evaluate the underlying biological role of necroptosis in clinical outcomes and tumor-infiltrating traits, we screened 2988 Ncluster related DEGs and determined 1819 prognostic genes by univariate Cox analysis (Table S4). Similarly, an unsupervised clustering method was used to accurately classify HCC patients into three stable gene clusters, which was termed as necroptosis-related gene cluster (NRGcluster), including 203 cases in NRGcluster A, 321 cases in NRGcluster B, and 75 cases in NRGcluster C (Figure 4A). The three NRGclusters had dramatic prognostic differences in HCC patients: NRGcluster B was proven to possess better prognostic outcomes, while patients in NRGcluster C were associated with poorer outcomes (Figure 4B). Furthermore, the three NRGclusters exhibited the different Nclusters, necroptosis gene expression (Figure S4A), among which NRGclusters had nearly Ncluster B subgroups and NRGcluster C was included in Ncluster A. Besides, the three NRGclusters had distinct expression profiles of necroptosis genes (Figure S4A, B). Among them, NRGcluster A or C had more expression of 51 than NRGcluster B while only six necroptosis regulators were highly expressed in Ncluster B compared to A or C, which revealed that the expression level of necroptosis genes had a negative relationship with HCC clinical prognosis and NRGcluster could represent necroptosis patterns of each HCC sample.
Figure 4 Identification of necroptosis genomic classification and construction of NRGscore. (A) Consensus clustering matrix for k =3. (B) Kaplan-Meier curves indicated necroptosis genomic phenotypes were markedly related to the OS of 599 patients in TCGA and ICGC cohorts, of which 203 cases were in gene cluster A, 321 cases in gene cluster B, and 75 cases in gene cluster C. (C) The partial likelihood deviance plot. (D) The Lasso regression coefficient profiles. (E) Coefficient of three core necroptosis-related DEGs. (F) The alluvial diagram exhibited the correlation of Ncluster, NRGcluster, NRGscore, and survival status. (G) Differences in NRGscore between two Nclusters in TCGA and ICGC cohorts. The upper and lower ends of the boxes represented the interquartile range of values. The lines in the boxes represented median value, and dots showed outliers. (H) Differences in NRGscore among three NRGclusters in TCGC and ICGC cohorts. P < 0.05 *; P < 0.01**; P < 0.001***. NRGscore, necroptosis-related gene score; OS, overall survival; TCGA, the cancer genome atlas; ICGC, international cancer genome consortium; Lasso, least absolute shrinkage and selection operator; DEGs, differentially expressed genes; Ncluster, necroptosis cluster; NRGclusters, necroptosis-related genes clusters.
Integrating the evidence above, we found that necroptosis played an irreplaceable role in affecting prognosis, remodeling TME, and regulating immune reaction in HCC. However, the Ncluster or NRGcluster could not accurately clarify individual necroptosis patterns, and to further explore the heterogeneity and complexity of necroptosis, the Lasso algorithm was employed to quantify individual patients. Based on the optimal value of λ (λ=3), we constructed a necroptosis scoring system, termed as NRGscore (Figures 4C, D) by adopting multivariate COX regression analysis. Three hub genes were comprised, including GTPBP4 (Coefficient = 0.5632), CDC20 (Coefficient = 0.271), and SLC22A1 (Coefficient = -0.1136) (Figure 4E). The NRGscore was calculated by RNA expression multiplied by its corresponding coefficient.
The alluvial diagram visualized the quantification changes of patients and exhibited the interaction of Ncluster, NRGcluster, NRGscore, and survival state (Figure 4F). Results indicated that Ncluster A was mainly linked to a higher NRGscore, whereas Ncluster B exhibited a lower score (Figure 4G). Remarkably, NRGcluster C showed the highest NRGscore, followed by NRGcluster A, while NRGcluster B revealed the lowest scores (Figure 4H). Then, we investigated the relationship between NRGscore and necroptosis-driving regulators (RIPK1, RIPk2, and MLKL) and found NRGscore was positively linked to the expression of these three regulators (Figure S5). Additionally, we further stratified patients into high NRGscore and low NRGscore groups according to the median cut-off value derived from the R Survminer package. Next, we investigated the RNA expression levels of necroptosis genes, and the results exhibited that most (46/49) genes were highly expressed and only three (3/49) were downregulated in the high NRGscore subgroup (Figure S6A). Taken together, NRGcluster could accurately represent overall necroptosis levels and the NRGscore system could depict individual necroptosis patterns of HCC patients.
NRGscore Is Highly Predictive of Clinical Outcomes for HCC
Then, the prognostic value of the necroptosis scoring system in predicting patients’ survival outcomes was also estimated. We divided HCC patients into two training and testing cohorts randomly, including 300 cases in the training set, and 299 cases in the test set. K-M curves demonstrated that HCC patients in the high NRGscore group had lower OS rates than their counterparts in overall, training, and test sets (P < 0.001) (Figures 5A, C). Next, ROC curves were plotted to appraise the accuracy of the scoring system in predicting survival at 1, 3, and 5 years. (Overall set: AUC at 1-, 3-, and 5-year is 0.766, 0.718, and 0.685; Training set: 0.774, 0.756, and 0.824, respectively; Test set: 0.762, 0.691, and 0.530, separately; Figures 5D, F). Besides, we assessed the NRGscore of each HCC case amidst the overall, training, and test sets, which implied that HCC patients in the low NRGscore group had better survival status and less ratio of dead status than the high NRGscore group (Figures S6B, D). Therefore, a high NRGscore was considered a high NRGrisk while a low NRGscore was deemed a low NRGrisk. In addition, to further evaluate the clinical application of the NRGscore system, we downloaded an external HCC dataset (GSE54236) and calculated the NRGscore of each HCC patient. The survival curve demonstrated that the high NRGrisk had a poorer prognosis than the low risk (Figure S7A), and the ROC curve indicated the NRGscore system had a robust predictive ability (Figure S7B). We also found patients with high NRGrisk had advanced clinical stage and grade than low NRGrisk (Table S5). In general, the above results indicated that the necroptosis scoring system had a promising capacity to predict the survival of HCC patients.
Figure 5 Evaluation of predictive power of NRGscore in training and test cohorts. (A, C) Kaplan-Meier curves showed that the high NRGrisk group had a more inferior OS than the low NRGrisk group in TCGA and ICGC cohorts, of which 599 cases were in the overall set (A), 300 cases were in training set (B), and 299 cases were in the test set (C). (D-F) ROC of NRGscore scheme: Areas under the curve of 1-, 3-, and 5-year OS in the overall set (D), training set (E), and test set (F). (G) Nomogram predicting 1-, 3-, and 5- year OS of patients based on NRGscore and other clinical parameters, including gender, age, T stage, and N stage. The red dot presented the point of each parameter, and the length of the line segment reflected the contribution of factors to the outcome event. (H) Calibration plot of the nomogram for predicting the probability of 1-, 3-, and 5-year OS. The colored line was the fit line and represented the predicted value (the horizontal axis) corresponding to the actual value (the vertical axis). The gray diagonal was the ideal case. P < 0.01**; P < 0.001***. TCGA, the cancer genome atlas; ICGC, international cancer genome consortium; NRGscore, necroptosis-related gene score; NRGrisk, necroptosis-related gene risk; ROC, receiver operating characteristic curves; OS, overall survival.
To better build a convenient and applicable clinical prognosis evaluation method, we plotted a nomogram graph based on the NRGscore system and clinical features including gender, age, and pathologic stage (Figure 5G). The concordance index of the nomogram was 0.707 (95% CI 0.652–0.762). The calibration plot for the possibility of 1-, 3-, and 5-year survival exhibited good coincidence degrees of survival probability between the prediction and real observations (Figure 5H). These results demonstrated that the nomogram could be an effective approach to predicting the prognosis outcome of HCC patients for clinicians. All the above results disclosed that NRGscore could have a satisfying clinical prediction value in HCC.
NRGscore Accurately Predicts Immunotherapeutic Benefits for HCC
Based on the above results, Ncluster A was characterized by more abundant immune infiltrating cells and enhanced immune function than Ncluster B. Ncluster A could represent the immune-inflamed type and B represent the immune-excluded phenotype, which gave us a hint that the Ncluster A subgroup of patients could respond better to immunotherapy than B cluster. As shown in the alluvial diagram, high NRGscore had an overwhelming ratio of Ncluster A while low NRGscore were nearly derived from Ncluster B. Therefore, we speculated that patients with high NRGscore could be more sensitive to immunotherapy. First, we performed GSEA to investigate the biological of high and low NRGrisk, and the results (Figure 6A) revealed that high NRGrisk was characterized by cell proliferation and DNA repair-related pathways (cell cycle, DNA replication, nucleotide exclusion repair, and spliceosome, et al.) while low NRGrisk was characterized by enhanced metabolism pathways (drug metabolism of P450, fatty acid metabolism, and bile acid metabolism, et al.). Alternatively, the above results were also verified by other famous pathways, including hallmark pathways (Figure S8A) and reactome pathways (Figure S8B). Surprisingly, almost all immune pathways were enriched in the high NRGrisk group (Figure S8C), which further demonstrated high NRGrisk resembled Ncluster A group and could be more sensitive to clinical immunotherapy. Additionally, we further calculated the TME parts, including stromal, immune, and ESTIMATE parts. As expected, the stromal component was significantly downregulated in the high NRGrisk score group (Figure 6B). Previous research underlined the hub role of stromal activation in resistance to checkpoint immunotherapy and high NRGrisk patients with low stromal could get more clinical benefits from immunotherapy. To comprehensively evaluate the correlation between the scoring system and immunotherapy, we first estimated the TICs abundance based on the NRGscore. Results indicated that core genes in the scoring system had a close relation to various types of TICs. Among them, CDC20 was mainly related to T immune cells; GTPBP4 had a tight correlation with NK immune cells; SLC22A1 was mostly correlated to macrophages and DC immune cells (Figure 6C). Various immune functions altered with the NRGscore system; the heatmap indicated that a higher NRGscore always followed more inhibition of immune response (APC co-inhibition, MHC-I class, T cell co-inhibition and Treg cell, IFN-response inhibition; Figure S8D). Besides, three TICs were significantly positively correlated with the NRGscore, including Mast cells activated (R=0.22), T cells follicular helper (R=0.26), and T cells CD4 memory activated (R=0.32), while five TICs were significantly negatively correlated, including B cells naïve (R=-0.24), NK cells activated (R=-0.24), macrophages M1 (R=-0.25), Mast cells resting (R=-0.3), and T cells CD4 memory resting (R=-0.36) (Figure S9). To our knowledge, PD-L expression was considered as an important predictor of response to ICBs, and we found that NRGscore was positively associated with several known immune checkpoints (PD-1, PD-L1, PD-L2, CTLA4, TIM-3, IDO1, LAG3, TIGIT; Table S6). Integrating the results above, we found the NRGscore system could be dramatically correlated to the immune landscape of HCC, and its role in estimating the immunotherapy deserved further exploration.
Figure 6 NRGscore could accurately predict response to ICBs. (A) Different biological features in high and low NRGrisk groups by GSEA in KEGG pathways. (B) Violin diagrams exhibited the correlation of immune, stromal, and ESTIMATE scores with NRGrisk groups. (C) Heatmap indicated the relevance of TICs and three hub genes in NRGscore. Different colors represented the different degrees of correlation. (D) T cell exclusion score was positively associated with NRGscore in the GSE54236 HCC cohort. (E) Waterfall diagram of NRGscore with responder or non-responder to CTLA-4 cohort. (F) Patients who respond to ICB possess more NRGscore than non-responders in the TCGA cohort by the TIDE algorithm. (G) NRGscore of patients varied with different responses to immunotherapy in the advanced urothelial cancer cohort. CR/PR is short for complete response or partial response; SD/PD represented the stable disease or progressive disease. P < 0.05 *; P < 0.01**; P < 0.001***. ICBs, immune checkpoint blockades; NRGscore, necroptosis-related gene score; TICs, tumor-infiltrating cells; NRGrisk, necroptosis-related gene risk; TIDE, tumor immune dysfunction, and exclusion.
For the present limitation research, the clinical application of immunotherapy in HCC was still not widespread and there are few biomarkers to predict the response to ICBs (23). Therefore, we evaluated the response of checkpoint immunotherapy to patients with different NRGscore. In the cohort of GSE54236, we found NRGscore was positively linked to the level of T cell exclusion, which indicated a high NRGscore with strong immune inhibition and could respond to ICB therapy (Figure 6D). In another cohort of metastatic melanoma, patients with response to CTLA-4 therapy possessed a higher NRGscore (Figure 6E). In addition, we also found similar results in the TCGA HCC cohort (Figure 6F) and IMvigor210 cohort treated with atezolizumab (Figure 6G). IPS values were elevated in high NRGscore subgroups (Table S7), which indicated that the high NRGscore subgroup was relatively responsive to immunotherapy. In other cancer cohorts with ICB therapy, the NRGscore system could also accurately predict ICB response (Figure S10). Collectively, the above results elucidated that the NRGscore system could play a non-negligible role in predicting the immunotherapy response in HCC patients, or even in other tumors.
Patients With High NRGrisk Were More Sensitive to Common Chemotherapy Regimens
As a widely accepted biomarker for prognosis and immunotherapy, tumor mutation burden (TMB) was also evaluated in our study. Although there were few differences in TMB between high and low NRGrisk subgroups (Figure S11A), and there is no linear correlation of TMB with NRGscore (Figure S11B), the gene mutation type altered in two subgroups in the waterfall curve (Figures 7A, B). Results indicated that the high NRGscore subgroup presented more extensive TMB than the low score subgroup. Among them, TP53 was the first mutated gene with a 39% mutation alteration rate in the high score subgroup, while CTNNB1 was the top rank mutated gene with a 29% mutation alteration rate in the low score subgroup (Figures 7A,B). Consistent with a recent study, CTNNB1, TP53 mutation-associated pathways, and metabolism profiles were identified in HBV-related HCC (55). Besides, we found that the ratio of tumor stem cells in HCC was significantly positively related to NRGscore (R=0.34) (Figure 7C). Considering the frequent use of chemotherapy in the treatment of HCC, we further explored the response of patients with 138 different types of drugs. Analyses of consequences revealed that several prevalent clinical chemotherapies of HCC exhibited low IC50 in the high NRGscore group, including Cisplatin, Camptothecin, Doxorubicin, Paclitaxel, and Sorafenib. The finding suggested that patients with high NRGscore were more sensitive to the treatment of chemotherapy drugs than those with low scores in HCC (Figures 7D, H).
Figure 7 Characteristics of NRGscore with tumor somatic mutation, stem cell, and chemotherapy sensitivity. (A, B) The waterfall plot of tumor somatic mutation was drawn in those with high NRGrisk group (A) and low NRGrisk (B) respectively. Each column indicated individual HCC patients, and the upper bar diagram exhibited TMB. The right number represented the mutation frequency, and the bar diagram on the right exhibited the proportion of each variant type, including missense, nonsense, splice, frameshift, and multiple mutations. (C) Correlation between NRGscore and component of stem cell. (D-H) Box diagrams showed the chemotherapy response between high NRGrisk and low NRGrisk groups: (D) Cisplatin, (E) Camptothecin, (F) Doxorubicin, (G) Paclitaxel, (H) Sorafenib. NRGscore, necroptosis-related gene score; NRGrisk, necroptosis-related gene risk; HCC, hepatocellular carcinoma; TMB, tumor mutation burden; IC50, half-maximal inhibitory concentration.
Discussion
The NRGscore, a tool designed to evaluate necroptosis patterns of each HCC patient, is a robust biomarker for predicting clinical outcomes and for guiding rational and effective immunotherapy. Our findings revealed a survival benefit trend was observed in Ncluster A compared to B, suggesting necroptosis played a non-negligible role in the progression and prognosis in HCC. GSVA analysis disclosed that Ncluster A had more TME-infiltrating cells and more complex immune-response than Ncluster B, substantiating Ncluster could be a surrogate and favorable biomarker for TME. However, individual necroptosis characteristics could not be elucidated dependent on Ncluster, and we built a novel NRGscore scheme, determined by three core genes (GTPBP4, CDC20, and SLC22A1). NRGscore could accurately and robustly predict 1-, 3-, and 5-year OS of HCC. More importantly, the high NRGscore group had more sensitivity to clinical chemotherapy and ICB-based immunotherapy than the low score group, which demonstrated the NRGscore could be meaningful guidance for HCC treatment.
Mounting evidence demonstrates that necroptosis plays a pivotal role in predicting clinical outcomes of cancer, regulating cancer progression and metastasis, remodeling TME, and thus affecting immunotherapy and chemotherapy (56). Although necroptosis had been shown to perform an antitumor function in cancers, much evidence also demonstrated that necroptosis may play a tumor-promoting role and trigger cancer metastasis (11, 57). The paradox phenomena could be explained by high inflammation TME caused by necroptosis and following elevated ROS levels could further promote cancer progression and metastasis (56). Our results uncovered that higher necroptosis patterns were, worse outcomes did HCC patients have due to most necroptotic regulators being risky genes in HCC. In addition, high NRGscore patients always had more mRNA levels of three necroptosis-driving genes (MLKL, RIPK1, and RIPK3), further validating that necroptosis could play a tumor-promoting role in HCC. In another research on liver cancer, necroptosis directly determines the cancer subtype by driving cell releasing damage-associated molecular patterns (DAMPs) which could reshape the TME (58), causing the switch from HCC to intrahepatic cholangiocarcinoma (ICC) (15). Here, we found that the Ncluster A subgroup had a cline to the high NRGscore group while Ncluster B resembled low NRGscore patients in the alluvial diagram. Ncluster A and high NRGscore were characterized by more TME-infiltrating cells and enhanced APC-co-inhibitor and checkpoint gene expression, corresponding to immune-inflamed phenotype. Ncluster B and low NRGscore were characterized by enhanced metabolism and decreased stem cells and proliferation, corresponding to immune-exclusion or metabolism phenotype. TIDE results suggested that the high NRGscore group with lower stromal activation could achieve clinical benefits in immunotherapy, consistent with research that emphasized the irreplaceable role of stromal activation in resistance to PD-1/PD-L1 inhibitors (45, 59, 60).
Furthermore, we elucidated the predictive value in another three independent cohorts (HCC cohort with TIDE results; metastatic melanoma cohort treated with Nivolumab; advanced urothelial cancer cohort with the intervention of atezolizumab (IMvigor210)). In line with the results of TCGA HCC, we observed a significant positive trend between NRGscore and T cell exclusion in the GSE54236 HCC cohort. In the metastatic melanoma cohort, responders to anti-CTLA4 possessed a higher NRGscore than non-responders. Similarly, patients’ sensitive to anti-PD-L1 had more NRGscore than non-responders. Furthermore, we input our NRGscore system into the TIDE website to predict other cancer cohorts undertaking ICBs and the results revealed that NRGscore could behave well in predicting response. We also found different somatic gene mutations in high and low NRGrisk groups and that TP53 mutation was most common in high risk and CTNNB1 for low risk. The results may help us to gain a comprehensive understanding of precision immunotherapy and promoting the necroptosis could enhance the response to ICBs. Intriguingly, we observed that NRGscore could perform as a candidate biomarker for classic chemotherapy and target therapy (Sorafenib, Doxorubicin, Cisplatin, Camptothecin, and Paclitaxel) for HCC.
Additionally, we have validated the pivotal role of three necroptosis-related genes in HCC. It has been widely accepted that the anaphase-promoting complex (APC) drives and governs cell cycle progression (anaphase initiation and late mitosis exit (61, 62)) by interacting with two essential activators - CDC20 and CDH1 (63). However, the two proteins performed different roles in tumor development. CDC20 was considered as an oncogene while CDH1 was deemed as a suppressor in multicancer (64–66). CDC20 has been reported to be overexpressed in HCC tissues and positively related to the tumor, TNM stage, and ki-67 expression (67). In tumorigenesis of HCC, APCCDC20 could stabilize the HIF-1a by degrading oxygen-dependent prolyl hydroxylase enzymes3 (PHD3) (68). Suppressing CDC20 (depletion of endogenous or pharmacological inhibitors) in diverse cancer cell lines led to a mitotic arrest and apoptosis, suggesting targeting CDC20 might be a novel anti-cancer therapy, especially in cancer with high expression (69). GTPBP was a GTP-binding protein that is involved in 60S ribosome biogenesis (70). The previous report found that GTPBP4 could reduce TP53 accumulation and increased expression of GTPBP4 correlated with reduced survival (71). In colorectal cancer, GTPBP4 was demonstrated as an oncogene that could disrupt the actin cytoskeleton and promote metastasis of cancer (72). In addition, depletion of GTPBP4 could inhibit cell proliferation, and high expression correlated with a worse prognosis than the low expression of GTPBP4 in breast cancer (73). In a study of HCC, as an RNA binding protein, GTPBP4 could stratify patients into low- and high-risk, which could predict survival outcomes (74). The upregulated expression of GTPBP4 promotes the proliferation of liver cancer cells and promotes the growth of tumors in mice, while the downregulated expression of GTPBP4 inhibits the proliferation of liver cancer cells and inhibits the growth of tumors in mice (75). SLC22A1 (OCT1) is one member of organic cation transporters, which could uptake intracellular inactivation, a broad spectrum of endogenous and exogenous substrates as well as anticancer drugs (76–78). The OCT1 activity was reported to correlate with the sensitivity of tyrosine kinase inhibitors (TKI) in patients with chronic myeloid leukemia (CML) (79). The downregulation of OCT1 is associated with tumor progression and a worse patient survival (80).
The results of our study should be further validated in a prospective cohort of HCC patients undergoing immunotherapy or chemotherapy. Also, these three hub genes in NRGscore and their relationship with necroptosis should be further explored in basic experiments. Since not all high NRGscore patients would benefit from immunotherapy, more clinical parameters and appropriate evaluations of immune infiltration and TMB should be taken into consideration.
In the current study, we found that necroptosis patterns could become a novel index to predict HCC patients’ outcomes and high necroptosis levels might suggest a poor survival prognosis. Ncluster was built based on these necroptosis regulators, and two Nclusters had immense changes in prognosis; cell phagocytosis-, endocytosis-, and peroxisome-related pathways; TME cell infiltration characterization; and immune function. In addition, we designed a novel individual necroptosis pattern evaluation system: NRGscore scheme. NRGscore could accurately and robustly predict 1-, 3-, and 5-year OS, which was validated in another HCC cohort. The novel nomogram was comprised of NRGscore and other clinicopathological features to visualize the risk of HCC patients. In ICB-based immunotherapy of HCC, NRGscore had meaningful guidance for patients and the high NRGscore group could be more sensitive to ICB treatment, which was further validated in the metastatic melanoma cohort and advanced urothelial cancer cohort. The high NRGscore group had an increased ratio of cancer stem cells than its counterpart, which further demonstrates high NRGscore could be more sensitive to chemotherapy and targeted molecular therapy.
Conclusions
In the study, we found high expression of necroptosis could affect the TME of HCC and be closely related to the immune-inflamed phenotype. In predicting response to ICB, we concluded that high NRGrisk patients could be more sensitive than low NRGrisk. We systematically characterized the landscapes of necroptosis in HCC patients and suggested that NRGscore could be a prognostic marker and help to interpret the responses of HCC to chemotherapies and immunotherapies, providing new strategies for the treatment of cancers.
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 authors.
Author Contributions
LL and YX designed the study. JZ and TH performed most of the results and completed the manuscript together. SZ, YZ, and SM helped with all experiments and data analysis. FX, TB, and YT helped write the manuscript. All authors read and approved the final manuscript.
Funding
This work was supported by the National Natural Science Foundation of China (8217113876, 62171365).
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.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2022.933210/full#supplementary-material
References
1. Bray F, Ferlay J, Soerjomataram I, Siegel RL, Torre LA, Jemal A. Global Cancer Statistics 2018: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin (2018) 68:394–424. doi: 10.3322/caac.21492
2. Kulik L, El-Serag HB. Epidemiology and Management of Hepatocellular Carcinoma. Gastroenterology (2019) 156:477–491.e1. doi: 10.1053/j.gastro.2018.08.065
3. Singal AG, Lampertico P, Nahon P. Epidemiology and Surveillance for Hepatocellular Carcinoma: New Trends. J Hepatol (2020) 72:250–61. doi: 10.1016/j.jhep.2019.08.025
4. McGlynn KA, Petrick JL, El-Serag HB. Epidemiology of Hepatocellular Carcinoma. Hepatol (Baltimore Md.) (2021) 73(Suppl 1):4–13. doi: 10.1002/hep.31288
5. Amin MB, Greene FL, Edge SB, Compton CC, Gershenwald JE, Brookland RK, et al. The Eighth Edition AJCC Cancer Staging Manual: Continuing to Build a Bridge From a Population-Based to a More "Personalized" Approach to Cancer Staging. CA: Cancer J Clin (2017) 67:93–9. doi: 10.3322/caac.21388
6. Bruix J, Reig M, Sherman M. Evidence-Based Diagnosis, Staging, and Treatment of Patients With Hepatocellular Carcinoma. Gastroenterology (2016) 150:835–53. doi: 10.1053/j.gastro.2015.12.041
7. Llovet JM, Zucman-Rossi J, Pikarsky E, Sangro B, Schwartz M, Sherman M, et al. Hepatocellular Carcinoma. Nat Rev Dis Primers (2016) 2:16018. doi: 10.1038/nrdp.2016.18
8. Siegel RL, Miller KD, Jemal A. Cancer Statistics, 2019. CA: Cancer J Clin (2019) 69:7–34. doi: 10.3322/caac.21551
9. Degterev A, Hitomi J, Germscheid M, Ch'en IL, Korkina O, Teng X, et al. Identification of RIP1 Kinase as a Specific Cellular Target of Necrostatins. Nat Chem Biol (2008) 4:313–21. doi: 10.1038/nchembio.83
10. Su Z, Yang Z, Xu Y, Chen Y, Yu Q. Apoptosis, Autophagy, Necroptosis, and Cancer Metastasis. Mol Cancer (2015) 14:48. doi: 10.1186/s12943-015-0321-5
11. Liu X, Zhou M, Mei L, Ruan J, Hu Q, Peng J, et al. Key Roles of Necroptotic Factors in Promoting Tumor Growth. Oncotarget (2016) 7:22219–33. doi: 10.18632/oncotarget.7924
12. Najafov A, Chen H, Yuan J. Necroptosis and Cancer. Trends Cancer (2017) 3:294–301. doi: 10.1016/j.trecan.2017.03.002
13. Moriwaki K, Balaji S, McQuade T, Malhotra N, Kang J, Chan FK. The Necroptosis Adaptor RIPK3 Promotes Injury-Induced Cytokine Expression and Tissue Repair. Immunity (2014) 41:567–78. doi: 10.1016/j.immuni.2014.09.016
14. Stoll G, Ma Y, Yang H, Kepp O, Zitvogel L, Kroemer G. Pro-Necrotic Molecules Impact Local Immunosurveillance in Human Breast Cancer. Oncoimmunology (2017) 6:e1299302. doi: 10.1080/2162402X.2017.1299302
15. Seehawer M, Heinzmann F, D'Artista L, Harbig J, Roux PF, Hoenicke L, et al. Necroptosis Microenvironment Directs Lineage Commitment in Liver Cancer. Nature (2018) 562:69–75. doi: 10.1038/s41586-018-0519-y
16. Finn RS, Qin S, Ikeda M, Galle PR, Ducreux M, Kim TY, et al. Atezolizumab Plus Bevacizumab in Unresectable Hepatocellular Carcinoma. N Engl J Med (2020) 382:1894–905. doi: 10.1056/NEJMoa1915745
17. Yau T, Park JW, Finn RS, Cheng AL, Mathurin P, Edeline J, et al. CheckMate 459: A Randomized, Multi-Center Phase III Study of Nivolumab (NIVO) vs Sorafenib (SOR) as First-Line (1L) Treatment in Patients (Pts) With Advanced Hepatocellular Carcinoma (aHCC). Ann Oncol (2019) 30:874–+. doi: 10.1093/annonc/mdz394.029
18. Zhu AX, Finn RS, Cattan S, Edeline J, Ogasawara S, Palmer DH, et al. KEYNOTE-224: Pembrolizumab in Patients With Advanced Hepatocellular Carcinoma Previously Treated With Sorafenib. J Clin Oncol 36 (2018) 36:209–209. doi: 10.1200/JCO.2018.36.4_suppl.209
19. Van Allen EM, Miao D, Schilling B, Shukla SA, Blank C, Zimmer L, et al. Genomic Correlates of Response to CTLA-4 Blockade in Metastatic Melanoma. Science (2015) 350:207–11. doi: 10.1126/science.aad0095
20. Havel JJ, Chowell D, Chan TA. The Evolving Landscape of Biomarkers for Checkpoint Inhibitor Immunotherapy. Nat Rev Cancer (2019) 19:133–50. doi: 10.1038/s41568-019-0116-x
21. Le DT, Uram JN, Wang H, Bartlett BR, Kemberling H, Eyring AD, et al. PD-1 Blockade in Tumors With Mismatch-Repair Deficiency. New Engl J Med (2015) 372:2509–20. doi: 10.1056/NEJMoa1500596
22. Ruiz de Galarreta M, Bresnahan E, Molina-Sánchez P, Lindblad KE, Maier B, Sia D, et al. β-Catenin Activation Promotes Immune Escape and Resistance to Anti-PD-1 Therapy in Hepatocellular Carcinoma. Cancer Discovery (2019) 9:1124–41. doi: 10.1158/2159-8290.CD-19-0074
23. Pinter M, Jain RK, Duda DG. The Current Landscape of Immune Checkpoint Blockade in Hepatocellular Carcinoma: A Review. JAMA Oncol (2021) 7:113–23. doi: 10.1001/jamaoncol.2020.3381
24. Tang R, Xu J, Zhang B, Liu J, Liang C, Hua J, et al. Ferroptosis, Necroptosis, and Pyroptosis in Anticancer Immunity. J Hematol Oncol (2020) 13:110. doi: 10.1186/s13045-020-00946-7
25. Kaczmarek A, Vandenabeele P, Krysko DV. Necroptosis: The Release of Damage-Associated Molecular Patterns and its Physiological Relevance. Immunity (2013) 38:209–23. doi: 10.1016/j.immuni.2013.02.003
26. Aaes TL, Kaczmarek A, Delvaeye T, De Craene B, De Koker S, Heyndrickx L, et al. Vaccination With Necroptotic Cancer Cells Induces Efficient Anti-Tumor Immunity. Cell Rep (2016) 15:274–87. doi: 10.1016/j.celrep.2016.03.037
27. Villa E, Critelli R, Lei B, Marzocchi G, Cammà C, Giannelli G, et al. Neoangiogenesis-Related Genes are Hallmarks of Fast-Growing Hepatocellular Carcinomas and Worst Survival. Results prospective study Gut (2016) 65:861–9. doi: 10.1136/gutjnl-2014-308483
28. Zubiete-Franco I, García-Rodríguez JL, Lopitz-Otsoa F, Serrano-Macia M, Simon J, Fernández-Tussy P, et al. SUMOylation Regulates LKB1 Localization and its Oncogenic Activity in Liver Cancer. EBioMedicine (2019) 40:406–21. doi: 10.1016/j.ebiom.2018.12.031
29. Mayakonda A, Lin DC, Assenov Y, Plass C, Koeffler HP. Maftools: Efficient and Comprehensive Analysis of Somatic Variants in Cancer. Genome Res (2018) 28:1747–56. doi: 10.1101/gr.239244.118
30. Glickman ME, Rao SR, Schultz MR. False Discovery Rate Control is a Recommended Alternative to Bonferroni-Type Adjustments in Health Studies. J Clin Epidemiol (2014) 67:850–7. doi: 10.1016/j.jclinepi.2014.03.012
31. Wilkerson MD, Hayes DN. ConsensusClusterPlus: A Class Discovery Tool With Confidence Assessments and Item Tracking. Bioinformatics (2010) 26:1572–3. doi: 10.1093/bioinformatics/btq170
32. Hänzelmann S, Castelo R, Guinney J. GSVA: Gene Set Variation Analysis for Microarray and RNA-Seq Data. BMC Bioinf (2013) 14:7. doi: 10.1186/1471-2105-14-7
33. Chen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling Tumor Infiltrating Immune Cells With CIBERSORT. Methods Mol Biol (Clifton N.J.) (2018) 1711:243–59. doi: 10.1007/978-1-4939-7493-1_12
34. Zhu J, Zhou Y, Wang L, Hao J, Chen R, Liu L, et al. CXCL5/CXCL8 is a Promising Potential Prognostic and Tumor Microenvironment-Related Cluster in Hepatocellular Carcinoma. J Gastrointest Oncol (2020) 11:1364–80. doi: 10.21037/jgo-20-556
35. Zhu J, Wang L, Zhou Y, Hao J, Wang S, Liu L, et al. Comprehensive Analysis of the Relationship Between Competitive Endogenous RNA (ceRNA) Networks and Tumor Infiltrating-Cells in Hepatocellular Carcinoma. J gastrointestinal Oncol (2020) 11:1381–98. doi: 10.21037/jgo-20-555
36. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene Set Enrichment Analysis: A Knowledge-Based Approach for Interpreting Genome-Wide Expression Profiles. Proc Natl Acad Sci United States America (2005) 102:15545–50. doi: 10.1073/pnas.0506580102
37. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, et al. Limma Powers Differential Expression Analyses for RNA-Sequencing and Microarray Studies. Nucleic Acids Res (2015) 43:e47. doi: 10.1093/nar/gkv007
38. Jiang P, Gu S, Pan D, Fu J, Sahu A, Hu X, et al. Signatures of T Cell Dysfunction and Exclusion Predict Cancer Immunotherapy Response. Nat Med (2018) 24:1550–8. doi: 10.1038/s41591-018-0136-1
39. Charoentong P, Finotello F, Angelova M, Mayer C, Efremova M, Rieder D, et al. Pan-Cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep (2017) 18:248–62. doi: 10.1016/j.celrep.2016.12.019
40. Gajewski TF, Schreiber H, Fu YX. Innate and Adaptive Immune Cells in the Tumor Microenvironment. Nat Immunol (2013) 14:1014–22. doi: 10.1038/ni.2703
41. Joyce JA, Fearon DT. T Cell Exclusion, Immune Privilege, and the Tumor Microenvironment. Science (2015) 348:74–80. doi: 10.1126/science.aaa6204
42. García-Mulero S, Alonso MH, Pardo J, Santos C, Sanjuan X, Salazar R, et al. Lung Metastases Share Common Immune Features Regardless of Primary Tumor Origin. J Immunother Cancer (2020) 8(1). doi: 10.1136/jitc-2019-000491
43. Necchi A, Joseph RW, Loriot Y, Hoffman-Censits J, Perez-Gracia JL, Petrylak DP, et al. Atezolizumab in Platinum-Treated Locally Advanced or Metastatic Urothelial Carcinoma: Post-Progression Outcomes From the Phase II IMvigor210 Study. Ann Oncol (2017) 28:3044–50. doi: 10.1093/annonc/mdx518
44. Hugo W, Zaretsky JM, Sun L, Song C, Moreno BH, Hu-Lieskovan S, et al. Genomic and Transcriptomic Features of Response to Anti-PD-1 Therapy in Metastatic Melanoma. Cell (2017) 168:542. doi: 10.1016/j.cell.2017.01.010
45. Mariathasan S, Turley SJ, Nickles D, Castiglioni A, Yuen K, Wang Y, et al. Tgfβ Attenuates Tumour Response to PD-L1 Blockade by Contributing to Exclusion of T Cells. Nature (2018) 554:544–8. doi: 10.1038/nature25501
46. Geeleher P, Cox N, Huang RS. Prrophetic: An R Package for Prediction of Clinical Chemotherapeutic Response From Tumor Gene Expression Levels. PloS One (2014) 9:e107468. doi: 10.1371/journal.pone.0107468
47. Carafa V, Nebbioso A, Cuomo F, Rotili D, Cobellis G, Bontempo P, et al. RIP1-HAT1-SIRT Complex Identification and Targeting in Treatment and Prevention of Cancer. Clin Cancer Res an Off J Am Assoc Cancer Res (2018) 24:2886–900. doi: 10.1158/1078-0432.CCR-17-3081
48. Schneider AT, Gautheron J, Feoktistova M, Roderburg C, Loosen SH, Roy S, et al. RIPK1 Suppresses a TRAF2-Dependent Pathway to Liver Cancer. Cancer Cell (2017) 31:94–109. doi: 10.1016/j.ccell.2016.11.009
49. Xiong J, Wang L, Fei XC, Jiang XF, Zheng Z, Zhao Y, et al. MYC is a Positive Regulator of Choline Metabolism and Impedes Mitophagy-Dependent Necroptosis in Diffuse Large B-Cell Lymphoma. Blood Cancer J (2017) 7:e0. doi: 10.1038/bcj.2017.61
50. Nicolai S, Pieraccioli M, Peschiaroli A, Melino G, Raschellà G. Neuroblastoma: Oncogenic Mechanisms and Therapeutic Exploitation of Necroptosis. Cell Death Dis (2015) 6:e2010. doi: 10.1038/cddis.2015.354
51. Chen W, Lin W, Wu L, Xu A, Liu C, Huang P. A Novel Prognostic Predictor of Immune Microenvironment and Therapeutic Response in Kidney Renal Clear Cell Carcinoma Based on Necroptosis-Related Gene Signature. Int J Med Sci (2022) 19:377–92. doi: 10.7150/ijms.69060
52. Dey A, Mustafi SB, Saha S, Kumar Dhar Dwivedi S, Mukherjee P, Bhattacharya R. Inhibition of BMI1 Induces Autophagy-Mediated Necroptosis. Autophagy (2016) 12:659–70. doi: 10.1080/15548627.2016.1147670
53. Nicolè L, Sanavia T, Cappellesso R, Maffeis V, Akiba J, Kawahara A, et al. Necroptosis-Driving Genes RIPK1, RIPK3 and MLKL-P are Associated With Intratumoral CD3(+) and CD8(+) T Cell Density and Predict Prognosis in Hepatocellular Carcinoma. J Immunother Cancer (2022) 10(3). doi: 10.1136/jitc-2021-004031
54. Ruf B, Heinrich B, Greten TF. Immunobiology and Immunotherapy of HCC: Spotlight on Innate and Innate-Like Immune Cells. Cell Mol Immunol (2021) 18:112–27. doi: 10.1038/s41423-020-00572-w
55. Gao Q, Zhu H, Dong L, Shi W, Chen R, Song Z, et al. Integrated Proteogenomic Characterization of HBV-Related Hepatocellular Carcinoma. Cell (2019) 179:561–77.e22. doi: 10.1016/j.cell.2019.08.052
56. Gong Y, Fan Z, Luo G, Yang C, Huang Q, Fan K, et al. The Role of Necroptosis in Cancer Biology and Therapy. Mol Cancer (2019) 18:100. doi: 10.1186/s12943-019-1029-8
57. Seifert L, Werba G, Tiwari S, Giao Ly NN, Alothman S, Alqunaibit D, et al. The Necrosome Promotes Pancreatic Oncogenesis via CXCL1 and Mincle-Induced Immune Suppression. Nature (2016) 532:245–9. doi: 10.1038/nature17403
58. Pasparakis M, Vandenabeele P. Necroptosis and its Role in Inflammation. Nature (2015) 517:311–20. doi: 10.1038/nature14191
59. Ravi R, Noonan KA, Pham V, Bedi R, Zhavoronkov A, Ozerov IV, et al. Bifunctional Immune Checkpoint-Targeted Antibody-Ligand Traps That Simultaneously Disable Tgfβ Enhance the Efficacy of Cancer Immunotherapy. Nat Commun (2018) 9:741. doi: 10.1038/s41467-017-02696-6
60. Nagarsheth N, Wicha MS, Zou W. Chemokines in the Cancer Microenvironment and Their Relevance in Cancer Immunotherapy. Nat Rev Immunol (2017) 17:559–72. doi: 10.1038/nri.2017.49
61. Fung TK, Poon RY. A Roller Coaster Ride With the Mitotic Cyclins. Semin Cell Dev Biol (2005) 16:335–42. doi: 10.1016/j.semcdb.2005.02.014
62. Weinstein J, Jacobsen FW, Hsu-Chen J, Wu T, Baum LG. A Novel Mammalian Protein, P55cdc, Present in Dividing Cells is Associated With Protein Kinase Activity and has Homology to the Saccharomyces Cerevisiae Cell Division Cycle Proteins Cdc20 and Cdc4. Mol Cell Biol (1994) 14:3350–63. doi: 10.1128/mcb.14.5.3350-3363.1994
63. Wang Z, Wan L, Zhong J, Inuzuka H, Liu P, Sarkar FH, et al. Cdc20: A Potential Novel Therapeutic Target for Cancer Treatment. Curr Pharm Des (2013) 19:3210–4. doi: 10.2174/1381612811319180005
64. Zhou Z, He M, Shah AA, Wan Y. Insights Into APC/C: From Cellular Function to Diseases and Therapeutics. Cell Div (2016) 11:9. doi: 10.1186/s13008-016-0021-6
65. Chang DZ, Ma Y, Ji B, Liu Y, Hwu P, Abbruzzese JL, et al. Increased CDC20 Expression is Associated With Pancreatic Ductal Adenocarcinoma Differentiation and Progression. J Hematol Oncol (2012) 5:15. doi: 10.1186/1756-8722-5-15
66. Kim JM, Sohn HY, Yoon SY, Oh JH, Yang JO, Kim JH, et al. Identification of Gastric Cancer-Related Genes Using a cDNA Microarray Containing Novel Expressed Sequence Tags Expressed in Gastric Cancer Cells. Clin Cancer Res (2005) 11:473–82. doi: 10.1158/1078-0432.473.11.2
67. Li J, Gao JZ, Du JL, Huang ZX, Wei LX. Increased CDC20 Expression is Associated With Development and Progression of Hepatocellular Carcinoma. Int J Oncol (2014) 45:1547–55. doi: 10.3892/ijo.2014.2559
68. Shi M, Dai WQ, Jia RR, Zhang QH, Wei J, Wang YG, et al. APC(CDC20)-Mediated Degradation of PHD3 Stabilizes HIF-1a and Promotes Tumorigenesis in Hepatocellular Carcinoma. Cancer Lett (2021) 496:144–55. doi: 10.1016/j.canlet.2020.10.011
69. Wang L, Zhang J, Wan L, Zhou X, Wang Z, Wei W. Targeting Cdc20 as a Novel Cancer Therapeutic Strategy. Pharmacol Ther (2015) 151:141–51. doi: 10.1016/j.pharmthera.2015.04.002
70. Fuentes JL, Datta K, Sullivan SM, Walker A, Maddock JR. In Vivo Functional Characterization of the Saccharomyces Cerevisiae 60S Biogenesis GTPase Nog1. Mol Genet Genomics (2007) 278:105–23. doi: 10.1007/s00438-007-0233-1
71. Lunardi A, Di Minin G, Provero P, Dal Ferro M, Carotti M, Del Sal G, et al. A Genome-Scale Protein Interaction Profile of Drosophila P53 Uncovers Additional Nodes of the Human P53 Network. Proc Natl Acad Sci U.S.A. (2010) 107:6322–7. doi: 10.1073/pnas.1002447107
72. Yu H, Jin S, Zhang N, Xu Q. Up-Regulation of GTPBP4 in Colorectal Carcinoma is Responsible for Tumor Metastasis. Biochem Biophys Res Commun (2016) 480:48–54. doi: 10.1016/j.bbrc.2016.10.010
73. Pawitan Y, Bjöhle J, Amler L, Borg AL, Egyhazi S, Hall P, et al. Gene Expression Profiling Spares Early Breast Cancer Patients From Adjuvant Therapy: Derived and Validated in Two Population-Based Cohorts. Breast Cancer Res (2005) 7:R953–64. doi: 10.1186/bcr1325
74. Wang L, Zhang Z, Li Y, Wan Y, Xing B. Integrated Bioinformatic Analysis of RNA Binding Proteins in Hepatocellular Carcinoma. Aging (Albany NY) (2020) 13:2480–505. doi: 10.18632/aging.202281
75. Chen J, Zhang J, Zhang Z. Upregulation of GTPBP4 Promotes the Proliferation of Liver Cancer Cells. J Oncol (2021) 2021:1049104. doi: 10.1155/2021/1049104
76. Zhang S, Lovejoy KS, Shima JE, Lagpacan LL, Shu Y, Lapuk A, et al. Organic Cation Transporters are Determinants of Oxaliplatin Cytotoxicity. Cancer Res (2006) 66:8847–57. doi: 10.1158/0008-5472.CAN-06-0769
77. Jonker JW, Schinkel AH. Pharmacological and Physiological Functions of the Polyspecific Organic Cation Transporters: OCT1, 2, and 3 (SLC22A1-3). J Pharmacol Exp Ther (2004) 308:2–9. doi: 10.1124/jpet.103.053298
78. Thomas J, Wang L, Clark RE, Pirmohamed M. Active Transport of Imatinib Into and Out of Cells: Implications for Drug Resistance. Blood (2004) 104:3739–45. doi: 10.1182/blood-2003-12-4276
79. Crossman LC, Druker BJ, Deininger MW, Pirmohamed M, Wang L, Clark RE. hOCT 1 and Resistance to Imatinib. Blood (2005) 106:1133–4. doi: 10.1182/blood-2005-02-0694
Keywords: necroptosis, tumor microenvironment, immunotherapy, chemotherapy, tumor-infiltrating cells, hepatocellular carcinoma
Citation: Zhu J, Han T, Zhao S, Zhu Y, Ma S, Xu F, Bai T, Tang Y, Xu Y and Liu L (2022) Computational Characterizing Necroptosis Reveals Implications for Immune Infiltration and Immunotherapy of Hepatocellular Carcinoma. Front. Oncol. 12:933210. doi: 10.3389/fonc.2022.933210
Received: 30 April 2022; Accepted: 06 June 2022;
Published: 07 July 2022.
Edited by:
Liang Cheng, Harbin Medical University, ChinaCopyright © 2022 Zhu, Han, Zhao, Zhu, Ma, Xu, Bai, Tang, Xu and Liu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Lei Liu, liulei84207@tmmu.edu.cn; Yungang Xu, yungang.xu@uth.tmc.edu
†These authors have contributed equally to this work and share first authorship