The Role of GNMT and MMP12 Expression in Determining TACE Efficacy: Validation at Transcription and Protein Levels

Introduction

Primary liver cancer is the sixth most common cancer in the world, with the third highest mortality rate, and hepatocellular carcinoma (HCC) accounts for 75%-80% of primary liver cancers.¹ Currently, the 5-year survival rate of HCC is 18%, one of the key reasons may be the high heterogeneity of HCC.^2,3 Due to this high heterogeneity, even when receiving the same treatment, there are considerable differences in patient outcomes, with several patients still not deriving benefits from the therapy.⁴ Therefore, the realization of precision individualized medical care is crucial for improving the prognosis of HCC patients.

Transarterial chemoembolization (TACE) is a treatment that leads to localized ischemic necrosis and induces cell death by infusion of embolic agents and chemotherapeutic drugs into the hepatic artery.⁵ It is not only a first-line therapy for HCC patients with intermediate stage [Barcelona Clinic Liver Cancer (BCLC)-B stage] disease but also an important platform for downstaging, conversion therapy and combination therapy.^6,7 The overall survival (OS) rate of HCC patients treated with TACE was 70.3% at 1 year, 51.8% at 2 years, 40.4% at 3 years, and 32.4% at 5 years.⁸ However, only a subset of patients exhibited a survival benefit from TACE.^8,9 At present, there are several criteria for evaluating the efficacy of TACE. The modified Response Evaluation Criteria in Solid Tumors (mRECIST) criteria can evaluate the changes in tumor areas after TACE and then evaluate the efficacy of TACE.¹⁰ The concept of “TACE refractoriness/TACE failure” was initially put forward by the Japan Society of Hepatology (JSH), it refers to a condition where certain patients did not respond effectively to TACE treatment, represented by outcomes such as tumor progression after two or more consecutive TACE procedures, inadequate responses in the targeted lesions (viable lesions > 50%), or the emergence of new tumor lesions despite changing the chemotherapeutic agents.¹¹ However, there is still no consensus on this concept in the international community, and TACE refractoriness does not fully represent the poor prognosis.¹² More importantly, the above criteria cannot reliably reflect the long-term prognosis of patients directly and cannot accurately stratify patients before treatment. Although studies have used clinical information or extracted imaging features to develop a model to predict the efficacy of TACE. Niu et al¹³ constructed a computed tomography-based radiomics nomogram for the prediction of TACE refractoriness, Maesaka et al¹⁴ indicated that hypovascular hepatic nodules are a predictive factor for TACE refractoriness. Most of these studies focused on the concept of TACE refractoriness, paid more attention to the deposition effect of lipiodol, or analyzed related anatomical or clinical features rather than directly assessing long-term prognosis based on the molecular mechanism that affects TACE. Fako et al¹⁵ and He et al¹⁶ had previously developed gene models to predict the prognosis of patients undergoing TACE. However, these models were not validated at the protein level, and the included patients were affected by treatment modalities beyond TACE.

In this study, we explored the molecular mechanisms that affect the effectiveness of TACE and constructed a gene signature to precisely predict the effectiveness of TACE in patients, which was validated at the protein expression level. Furthermore, given the encouraging results of interventional therapy combined with targeted therapy and immunotherapy, we also investigated the correlation between TACE effectiveness and the immune microenvironment and immunotherapy of HCC.

Materials and Methods

Data Collection

A total of three datasets were included in this study. GSE104580 (Figure 1A) contained 81 tissue RNA samples from TACE responders and 66 tissue RNA samples from TACE non-responders. All samples were derived from tumor biopsies of patients with HCC before TACE treatment. Extramural reviewers utilized the modified Response Evaluation Criteria in Solid Tumors to assess the response to TACE within 3 months. Within this study, the primary endpoint for the final grouping was the tumor response. Patients achieving a complete response (CR) or partial response (PR) were categorized as TACE responders, while those with stable disease (SD) or progressive disease (PD) were classified as non-responders. Before reaching the primary endpoint, TACE was the only treatment.¹⁷ We randomly divided GSE104580 into two cohorts at a ratio of 3:2: cohort 1 (n = 89, 49 responders and 40 non-responders) was used as the training set, and cohort 2 (n = 58, 32 responders and 26 non-responders) was used as the verification set. GSE14520 (Figure 1B) contained tumor tissue mRNA transcriptional data of 247 Chinese patients with HCC at the initial stage of treatment. After excluding 7 patients with incomplete clinical information and 51 patients who received treatment other than TACE or surgical resection, the remaining 189 patients were divided into the resection group (n = 85) and TACE group (n = 104), according to whether patients received TACE treatment after resection. Among the TACE group, 74 patients received adjuvant TACE after surgical tumor resection (adjuvant set), and 30 patients received TACE therapy only after tumor relapse (relapse set). Notably, the purpose of adjuvant TACE is to reduce the risk of tumor recurrence.¹⁸ Patients who were considered to have a high probability of recurrence, such as those with vascular invasion, multiple tumor nodules, or tumor size > 10 cm, received adjuvant TACE therapy.^15,19,20 The expression profiles of GSE104580 and GSE14520 were downloaded from the GEO database (http://www.ncbi.nlm.nih.gov/GEO). Then, we included patients who were treated with TACE as the initial treatment at the Cancer Hospital Chinese Academy of Medical Sciences from February 2014 to February 2018. Inclusion criteria for the study were as follows: (1) diagnosis with HCC by biopsy, (2) availability of pathological HCC tissues, (3) age ≥ 18, (4) Child-Pugh stage A or B, (5) ECOG performance score ≤ 1, (6) initial treatment of TACE, and (7) life expectancy of over 3 months. The exclusion criteria were as follows: (1) use of another antitumor treatment after TACE, (2) occurrence of other malignant tumors, (3) extrahepatic metastasis and (4) incomplete clinical data. Finally, 29 samples of available archived formalin-fixed paraffin-embedded (FFPE) puncture pathological tissue of HCC were included in the external validation set (Figure 1C). The purpose of selecting puncture tissue before treatment was to avoid the influence of resection on the research as much as possible. The institutional review board approved the use of these FFPE tissues and waived the requirement for informed consent. The clinical baseline characteristics of patients in the GSE14520 and external verification sets are shown in Tables 1 and 2, respectively.

Table 1 The Clinical Baseline Characteristics of Patients in GSE14520

Table 2 The Clinical Baseline Characteristics of Patients in the External Validation Set

Figure 1 The grouping and inclusion/exclusion of datasets in this study. (A) The description and grouping of GSE104580. (B) After applying the exclusion criteria, group GSE14520 was based on the treatment approach. (C) The description and inclusion/exclusion of external validation set. Abbreviation: TACE, transarterial chemoembolization.

Development of the TACE Effectiveness Model

Identification of Differentially Expressed Genes

The “Limma” package was used to identify the differentially expressed genes between TACE responders and non-responders in cohort 1 (|Fold Change| > 1.3, p < 0.05). Studies have shown that the early recurrence of tumors after resection was related to the early dissemination of the primary tumor,^18,21 which may mean that in the adjuvant TACE set, patients with recurrence within 1 year may have poor reactivity to TACE. Similar to Wei et al²⁰ to investigate the genetic signature related to TACE responsiveness, we compared the differences in transcriptional characteristics between patients with recurrence within 1 year (n = 17) and without recurrence after 2 years (n = 33) in the adjuvant set (|fold change| > 1.3, p < 0.05), and we enrolled these 50 patients in cohort 3 (Figure 1B). Finally, the common differentially expressed genes of cohort 1 and cohort 3 were obtained and visualized by the “venndiag” package.

Calculation of TACE Effectiveness Score

The intersection of differential genes was used to construct a TACE effectiveness model in cohort 1. The “glmnet” R package was used to integrate the differential gene expression data and patient response to TACE, and logistic regression analysis was utilized to develop the prediction model in cohort 1. Subsequently, the regression coefficients of each signature gene and the formula for calculating the effectiveness score of TACE were determined:

To balance the false positive rate and false negative rate, we set the threshold of this model to 0.5. In general, the TACE effectiveness score of patients was calculated based on the expression of signature genes and the corresponding regression coefficients. When the score was greater than 0.5, we speculated that this patient may exhibit better responsiveness to TACE (high TACE effectiveness), indicating a more favorable prognosis associated with undergoing the TACE treatment. Conversely, if the score was below 0.5, the patient’s prognosis following TACE treatment may be less favorable (low TACE effectiveness). Eventually, all patients were divided into the TACE high effectiveness group (HE group) and the TACE low effectiveness group (LE group).

Verification of the TACE Effectiveness Model

First, receiver operating characteristic (ROC) curve analysis was used to illustrate the accuracy of our model in evaluating the effectiveness of TACE in HCC patients. In addition to the evaluation of cohort 1 and cohort 2, we also evaluated the prediction accuracy of our model in cohort 3. As mentioned earlier, we regard patients with recurrence within 1 year as TACE low effectiveness and patients without recurrence after 2 years as TACE high effectiveness. Subsequently, the “survival” and “survfit” packages were utilized to develop Kaplan‒Meier survival curves to analyze the difference in prognosis between the HE and LE groups in GSE14520. OS was defined as the time interval from initial TACE treatment to death from any cause, disease-free survival (DFS) was defined as the time interval between initial TACE treatment (adjuvant TACE) and the recurrence of disease or death from any cause, and progression-free survival (PFS) was defined as the time interval from initial TACE treatment to tumor progression or death for any reason. Furthermore, a survival heatmap was used to show the relationship between patient prognosis, TACE effectiveness score and signature gene expression. The prcomp function in the “stats” package was used for principal component analysis (PCA).

Immunohistochemistry Analysis

To further verify the ability of our signature to predict the effectiveness of TACE, immunohistochemical (IHC) staining of the paraffin-embedded HCC biopsy tissues of patients in the external verification was performed according to the manufacturer’s instructions. Anti-GNMT antibody (Abcam, #ab224753, 1:1000) and anti-MMP12 antibody (Abcam, #ab137444, 1:100) were used for IHC. The staining score indicative of relative protein expression levels was calculated according to the mean optical density (OD) value of the randomly selected visual field by ImageJ (2.9.0/1.53t). Following this, we employed the staining scores of the two proteins as relative protein expression levels to calculate the relative TACE effectiveness scores for each patient in the external validation set using the same formula and methodology. Subsequently, patients were categorized into HE and LE groups based on the median score of relative TACE effectiveness scores. The above operations were performed by two independent pathologists.

Immune Infiltration Analysis

CIBERSORTx²² (https://cibersortx.stanford.edu) was used to evaluate and quantify the infiltration levels of various immune cells in each sample. Subsequently, we further compared the difference in immune cell infiltration between the HE and LE groups, as well as TACE responders and non-responders. In addition, simulation studies were carried out using the Tumor Immune Dysfunction and Exclusion (TIDE) algorithm (http://tide.dfci.harvard.edu/), which can simulate the primary mechanisms of tumor immune escape, including T-cell dysfunction in tumors, cytotoxic T lymphocyte (CTL) invasion levels, and the exclusion properties of T cells in tumors with low levels of CTL invasion, and the final TIDE score can be given to predict the tumor response to immunotherapy.²³

Functional Enrichment Analysis

The potential molecular mechanisms and biological functions related to TACE effectiveness were explored by gene set enrichment analysis (GSEA) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis. The GSEA method was used to understand the enriched pathways in both groups. The molecular signature dataset was analyzed by Java GSEA v.4.2.2 and h. all. v7.5.1 symbols. gmt [Hallmarks] [threshold criteria: |NES| > 1.4, p < 0.05 and false discovery rate (FDR) < 0.21]. The differentially expressed genes (|fold change| > 1.2, p < 0.05) between TACE responders and non-responders in GSE104580, and the HE and LE groups in the TACE group were included in the KEGG enrichment analysis (https://www.kegg.jp/kegg/rest/keggapi.html, FDR < 0.1, p < 0.05).

Statistical Analysis

R version 3.6.1 and 4.1.1 were used to conduct statistical analyses. Pearson’s or Spearman’s rank correlation coefficients were used to investigate the relationship between variables. The continuous variables and normal distributions recorded for the two groups were compared using a t-test or Mann‒Whitney U-test. Categorical variables were compared using Fisher’s exact or chi-squared tests. The Kaplan‒Meier (KM) method was used to plot receiver operating characteristic (ROC) curves. The statistical significance was determined using Log rank tests (statistical significance: p < 0.05).

Results

Construction of a Gene Signature for Predicting the Effectiveness of TACE

As shown in Figure 2, we first identified the differentially expressed genes between TACE responders and non-responders in cohort 1 and cohort 3 (Figure 2A and B). Then, we obtained the intersection of the differentially expressed genes in the two cohorts (Figure 2C), namely, matrix metallopeptidase 12 (MMP12), glycine N-methyltransferase (GNMT) and anterior gradient 2 (AGR2). Among them, MMP12 and AGR2 exhibited a downregulated trend, while GNMT exhibited an upregulated trend in TACE responders compared with TACE non-responders (Figure 2A and B).

Figure 2 Screening of differential genes related to the effectiveness of TACE. Differentially expressed genes between responders and non-responders in cohort 1 (A) and cohort 3 (B). (C) The Venn diagram showed the intersection of differential genes related to TACE effectiveness in Cohort 1 and Cohort 3. Abbreviations: MMP12, matrix metalloproteinase-12; GNMT, glycine N-methyltransferase; AGR2, anterior gradient 2.

Afterward, using these three genes, we constructed a gene signature for calculating the TACE effectiveness score of patients by logistic regression. Ultimately, the TACE effectiveness score for each patient was calculated based on the following equation:

TACE Effectiveness Score = (expression level of GNMT * 0.393480404) - (expression level of MMP12 * 0.295095776) - 0.882549992

All patients were divided into HE and LE groups according to the TACE effectiveness score, the grouping threshold was established at 0.5. When the effectiveness score surpasses 0.5, the patient falls into the HE group and, conversely falls into the LE group.

The 2-Gene Signature Can Accurately Predict the Effectiveness of TACE in Patients

ROC analysis illustrated the accuracy of our 2-gene signature in predicting TACE effectiveness. As shown in Figure 3A-C, the area under the curve (AUC) was 0.83 in cohort 1, 0.83 in cohort 2 and 0.80 in cohort 3.

Figure 3 Evaluation of the predictive ability of the TACE effectiveness model by the receiver operating characteristic (ROC) analysis and Kaplan-Meier (KM) analysis. The ROC curves of our model for predicting patients’ response to TACE in cohort 1 (A), cohort 2 (B), and cohort 3 (C). KM analysis of overall survival and disease-free survival between low effectiveness and high effectiveness groups in relapse set (D), adjuvant set (E and F), and TACE group (G). Abbreviations: AUC, the area under the curve; OS, overall survival; DFS, disease-free survival; pval, p-value; HR, hazard ratio; CI, confidence interval.

Subsequently, we compared the OS and DFS between patients in the HE and LE groups in each subgroup of GSE14520, including the relapse set, adjuvant set, TACE group and resection group. For patients with a high risk of recurrence, adjuvant TACE was conducted, and other patients were treated with TACE at the time of recurrence. At this time, the long-term prognosis of patients depends largely on their responsiveness to TACE.^15,20 In the relapse set, the OS of patients in the HE group was significantly higher than that of patients in the LE group (Figure 3D, p < 0.05), which indicated that the HE group had a higher responsiveness to TACE than the LE group among patients who received only TACE treatment after tumor recurrence. In the adjuvant set, patients in the HE group had a significantly better OS and DFS than the LE group and resection group (Figure 3E and F, p < 0.05), and there was no significant difference in OS and DFS between the LE group and resection group (Figure 3E and F), which demonstrated that most of the patients in the LE group did not benefit significantly from TACE treatment. In the whole TACE group, patients in the HE group still had a better prognosis, and no significant difference was found between the LE and resection groups (Figure 3G). ROC curve analysis showed the accuracy of our model in predicting OS and PFS (Supplementary Figure 1A and B). Interestingly, when comparing the prognosis of HE and LE group patients in the resection group, we found no significant difference in OS or DFS (Figure 4A and B), and the univariate and multivariate Cox regression analyses revealed that the grouping generated by our model does not emerge as a significant factor influencing the prognosis of surgical patients, which indicated the specificity of our model (Supplementary Figure 1C-F).

Figure 4 Evaluation of the prediction performance of our model by Kaplan-Meier (KM) analysis, principal component analysis, and survival heat map. KM analysis of overall survival and disease-free survival between low and high effectiveness groups in resection group (A, B). The survival heatmap displayed the distribution of patients’ effectiveness score, prognostic performance, and the expression of signature genes in the TACE group (C), adjuvant set (D), and relapse set (E). The principal component analysis of patients in cohort 1 (F), cohort 2 (G), GSE104580 (H), and TACE group (I) is based on the TACE effectiveness model. Abbreviations: TACE, transarterial chemoembolization; MMP12, matrix metalloproteinase-12; GNMT, glycine N-methyltransferase; OS, overall survival; DFS, disease-free survival; pval, p-value; HR, hazard ratio.

The survival heatmap revealed that with the decrease in MMP12 expression and the increase in GNMT expression, the effectiveness score increased gradually; when the effectiveness score was less than 0.5, the OS of patients was lower, and more patients died, and when the score was higher than 0.5, the prognosis of patients was better (Figure 4C-E). Furthermore, the PCA yielded a two-way distribution for patients between the HE and LE groups in various cohorts (Figure 4F-I).

Finally, the statistically significant results of univariate regression were incorporated into multivariate Cox regression analysis, which demonstrated that the LE group was an independent risk factor for OS (HR: 3.890, 95% CI: 1.673–9.043, p < 0.05) and DFS (HR: 2.451, 95% CI: 1.119–5.369, p < 0.05) in patients treated with TACE (Figure 5A-D).

Figure 5 Univariate and multivariate analysis to evaluate the factors influencing the effectiveness of TACE. Univariate analysis to evaluate the factors influencing overall survival (OS) and disease-free survival (DFS) of patients treated with TACE in TACE group (A) and adjuvant set (B). The statistically significant factors identified through univariate analysis were included in the multivariate analysis to assess the independent risk factors influencing OS and DFS in patients treated with TACE in the TACE group (C) and adjuvant set (D). Bolded text indicates statistical significance (p < 0.05). Abbreviations: LE, low effectiveness; HE, high effectiveness; AFP, alpha-fetoprotein; BCLC, Barcelona Clinic Liver Cancer; TACE, transarterial chemoembolization; ALT, alanine aminotransferase; CLIP, Cancer of the Liver Italian Program; OS, overall survival; DFS, disease-free survival; CI, confidence interval.

Further Verification of IHC Analysis

A total of 29 patients were selected for inclusion in the external validation set. The baseline characteristics of the two groups are shown in Table 2, and there was no significant difference in baseline data including age, TACE cycle, number and dimension of tumor between the two groups (p > 0.05). After obtaining the IHC staining scores of GNMT and MMP12 in the external verification set, we scored the relative TACE effectiveness of each patient through our model and then similarly divided the patients into the HE and LE groups (Figure 6A). The relative expression of GNMT was significantly higher and MMP12 was lower in the HE group than in the LE group (Figure 6B and C, p < 0.05). Additionally, KM survival and univariate and multivariate Cox regression analyses were performed. The results of the KM survival curve illustrated that TACE was more effective for HE group patients (Figure 7A and B, p < 0.05). The median OS was 881 days (95% CI: 719–1043 days) in the HE group and 273 days (95% CI: 190–356 days) in the LE group. The median PFS of the HE group was 458 days (95% CI: 309–607 days), and that of the LE group was only 136 days (95% CI: 82–190 days). We included the statistically significant indicators in univariate regression (Supplementary Figure 2A and B) in multivariate regression analysis (Figure 7C and D) and finally verified that LE status was an independent risk factor for OS (HR: 2.661, 95% CI: 1.042–6.797, p < 0.05) and PFS (HR: 2.607, 95% CI: 1.054–6.448, p < 0.05). In addition, the TACE cycle and tumor dimension were also independent influencing factors for OS (Figure 7C).

Figure 6 Immunohistochemical (IHC) staining analysis on external validation set to validate the expression trends of the signature genes. (A) Hematoxylin eosin staining and the IHC staining analysis of GNMT and MMP12 in high and low effectiveness group patients. The difference of GNMT (B) and MMP12 (C) relative expression levels between high and low effectiveness group patients. (* p < 0.05; *** p < 0.001). Abbreviations: MMP12, matrix metalloproteinase-12; GNMT, glycine N-methyltransferase; HE, hematoxylin-eosin staining.

Figure 7 Validation for the predictive ability of the TACE effectiveness model by Kaplan-Meier (KM) analysis and multivariate analysis in external validation set. KM analysis of overall survival (A) and progression-free survival (B) in external validation set. The statistically significant factors identified through univariate analysis were included in the multivariate analysis to assess the independent risk factors influencing overall survival (C) and progression-free survival (D) in patients undergoing TACE in an external validation set. Bolded text indicates statistical significance (p < 0.05). Abbreviations: LE, low effectiveness; HE, high effectiveness; AFP, alpha-fetoprotein; BCLC, Barcelona Clinic Liver Cancer; TACE, transarterial chemoembolization; OS, overall survival; DFS, disease-free survival; HR, hazard ratio; CI, confidence interval.

Immune Infiltration Analysis

We investigated the relationship between the HCC immune microenvironment and the efficacy of TACE. As shown in Figure 8A, in terms of the TACE group, the infiltration level of M0 macrophages and activated mast cells in the LE group was significantly higher than that in the HE group (p < 0.05), while the level of resting mast cells was significantly lower than that in the HE group (p < 0.05). The same results were verified in GSE104580, whether it was non-responders vs responders or LE vs HE divided according to the model (Figure 8B and C). The TIDE algorithm was then applied to gain insights into the mechanisms involved in HCC immune evasion and forecast the potential response to immunotherapy for patients with different efficacies of TACE. Figure 8D showed that the TIDE scores of HE group patients were significantly lower than those of LE group patients (p < 0.05), which indicated that patients who exhibit better responsiveness to TACE may also demonstrate a favorable response to immunotherapy. This result was verified in GSE104580 (LE vs HE) (Figure 8F), and showed the same trend in non-responder vs responder, while it did not reach statistical significance (Figure 8E). Similarly, an analogous result was found in terms of T-cell exclusion; that is, more significant T-cell exclusion occurred in patients with poor TACE effectiveness. In addition, HE or responder group patients had higher microsatellite instability (MSI) expression and fewer myeloid-derived suppressor cells (MDSCs) in both datasets (p < 0.05, Figure 8D–F). Research has shown that CD24⁺ HCC cells possessed a significantly high survival ability after exposure to doxorubicin under hypoxia,²⁰ and there was a significant negative correlation between the expression of CD24 and our effectiveness score (r = −0.54 in GSE104580 and −0.40 in GSE14520, p < 0.05, Supplementary Figure 3A and B).

Figure 8 Immune infiltration and immune simulation analysis to evaluate the infiltration levels of various immune cells and the efficacy of immunotherapy between different groups. The immune infiltration analysis showed the differential infiltration of various immune cells between the low and high effectiveness group in TACE group (A), the non-responder and responder groups in GSE104580 (B), low and high effectiveness group in GSE104580 (C). The immune simulation analysis showed the difference in various immune scores between the low and high effectiveness groups in the TACE group (D), non-responder and responder groups in GSE104580 (E), and LE and HE group in GSE104580 (F). TIDE, tumor immune dysfunction, and exclusion; IFNG, interferon-gamma; MSI, microsatellite instability; MDSC, myeloid-derived suppressor cells; CAF, cancer-associated fibroblasts; TAM, tumor-associated macrophages; LE, low effectiveness; HE, high effectiveness; Re, responder. (* p < 0.05; ** p < 0.01; *** p < 0.001). Abbreviation: ns, no statistical significance.

Functional Enrichment Analysis Based on the Signature

To explore the mechanism of the poor response to TACE in the LE group, GSEA and KEGG enrichment analyses were conducted. First, GSEA showed that the hypoxia signaling pathway, unfolded protein response, G2M checkpoint, mechanistic target of rapamycin complex 1 (mTORC1) signaling and MYC targets version 1 pathway were enriched in the LE group (|NES| > 1.4, p < 0.05 and false discovery rate (FDR) < 0.21). These results were consistent in the TACE group and GSE104580 (Figure 9A and B). Additionally, KEGG enrichment analysis was carried out. Ninety-seven differentially expressed genes between responders and non-responders in GSE104580 and 98 differentially expressed genes between the HE and LE groups in the TACE group were screened out (|Fold Change| > 1.2, p < 0.05). The results revealed that the differentially expressed genes related to the efficacy of TACE were mainly enriched in metabolic pathways, the hypoxia-inducible factor-1 (HIF-1) signaling pathway, the glycolysis pathway and the bile secretion pathway (FDR < 0.1, p < 0.05, Figure 9C and D).

Figure 9 Exploration for related pathways influencing the effectiveness of TACE by functional enrichment analysis. Gene set enrichment analysis exhibited that five hallmark pathways were all enriched in the low effectiveness patients in the TACE group (A) and non-responders in GSE104580 (B). Kyoto Encyclopedia of Genes and Genomes enrichment analysis exhibited that differentially expressed genes between low and high effectiveness groups in TACE (C), non-responder, and responder groups in GSE104580 (D) were enriched in several same pathways. Abbreviation: TACE, transarterial chemoembolization.

Discussion

Currently, TACE holds a crucial position within the treatment paradigm of HCC. At least 30% of patients receive TACE treatment in the course of HCC treatment,²⁰ whether as a first-line therapeutic approach or as a pivotal component of combination or conversion therapy. However, even when the same TACE procedure is applied, patients may exhibit substantial differences in survival benefits. The majority of previous studies have focused on the performance of lipiodol deposition or the prediction of patient prognosis based on other imaging features or clinical characteristics, while this phenomenon has not been explored at the molecular level (which may be attributed to the practical clinical situation: the diagnosis of HCC does not rely on pathological examination, resulting in a limited availability of sufficient TACE treatment-related samples). Indeed, the true biological status of lesions filled with lipiodol remains unclear, and the precise molecular mechanisms underlying TACE treatment and its impact on efficacy have yet to be determined.

In this study, we collected a total of three datasets for model construction and validation. Among them, patients from the adjuvant TACE group were considered to have a high risk of recurrence, as they were likely to have residual tumor cells even after surgical resection. Therefore, we hypothesize that the OS and DFS of these patients primarily depend on the efficacy of adjuvant TACE, which is consistent with previous research. Additionally, FFPE samples from HCC patients who received TACE monotherapy were collected as an external validation set. The purpose of selecting patients receiving TACE monotherapy was to minimize the impact of other treatment interventions on the outcomes, thereby enhancing the reliability of the validation results. However, as mentioned earlier, the main disadvantage of this approach is the limited sample size.

Our 2-gene signature can accurately differentiate between HE and LE patients. The results indicated that in various cohorts, including the external validation set, patients in the HE group exhibited significantly higher OS, DFS and PFS than those in the LE group. The study findings revealed that the median OS for TACE treatment in HCC was 19.4 months, accompanied by a median PFS spanning from 3 to 9 months.⁸ Within this study, LE group patients displayed a notably lower median OS of only 237 days and a median PFS of 136 days, the HE group exhibited a substantially higher median OS of 881 days and a median PFS of 458 days. Furthermore, in the study conducted by Fako et al,¹⁵ TACE Non-Responders in the TACE group demonstrated a median OS of 39.5 months, whereas, in our study, the LE group of patients within the TACE group showed a lower median OS of only 31.9 months. Moreover, in comparison to the resection group, whether in the TACE group or the adjuvant set, patients in the LE group did not experience a significant improvement in prognosis, even receiving additional TACE treatment after resection. This suggests that patients in the LE group did not derive significant benefits from TACE therapy. Furthermore, our model did not demonstrate any efficacy in the resection group, indicating that it was effective only in patients who received TACE treatment. This observation suggests that the signature genes were not solely associated with prognosis but also have underlying connections to the molecular mechanisms of TACE therapy. Multiple studies have demonstrated that GNMT, as a tumor suppressor in HCC, can inhibit HCC proliferation through various pathways.^24–27 Loss of GNMT has been shown to result in spontaneous hepatic steatosis and HCC in mice.²⁴ What truly intrigues us, however, is that Curtis et al²⁸ discovered that the lack of GNMT can reduce glucose production from glycogen and gluconeogenic sources. This effect may occur through the redirection of glucose precursor flux toward pathways involved in regulating the liver redox state and tumorigenesis, including de novo glutathione synthesis and the pentose phosphate pathway. Theoretically, the reduction in hepatic glucose production not only contributes to the occurrence and development of HCC but also sustains the growth and proliferation of HCC by conserving nutrients and energy resources.^29–32 Curtis et al found an increase in glucose-6-phosphate dehydrogenase (G6PDH) and a decrease in the ratio of glucose-6-phosphate (G6P) to ribulose 5-phosphate (R5P) in GNMT-knockout (GNMT-KO) mice. The former is a key regulator of the flux in the oxidative pentose phosphate pathway, which is routinely observed to be higher in HCC models with elevated pentose phosphate pathway flux,³³ and the latter serves as an indicator of increased flux in the oxidative pentose phosphate pathway.³⁴ Importantly, the reducing potential generated through the pentose phosphate pathway can alleviate oxidative stress and provide protection to tumor cells.³⁵ This implied that the loss of GNMT may enable HCC cells to develop resistance against oxidative stress induced by TACE, potentially by increasing the flux in the oxidative pentose phosphate pathway. On the other hand, increased de novo synthesis of glutathione (GSH) was observed in the livers of GNMT-KO mice. GSH, as a crucial intracellular antioxidant, may also serve as a significant factor in combating oxidative stress responses induced by TACE therapy.³⁶ In addition, the loss of GNMT can stimulate dormant hepatic cancer stem/progenitor cells to replicate,³⁷ and a study found that hepatic cancer stem cells are highly involved in resistance to hypoxia and chemotherapy.²⁰ MMP12, a member of the matrix metalloproteinase (MMP) family, plays a crucial role in extracellular matrix degradation, and its overexpression is associated with vascular invasion in HCC.³⁸ Veronica et al³⁹ demonstrated that levels of tryptase stored in mast cell granules can predict patients’ response to TACE, whose mechanism may involve the activation of MMPs by tryptase, leading to the degradation of the extracellular matrix, promoting the release of angiogenic factors from the matrix, and providing space for neovascular growth. In addition, the study illustrated that mast cells can promote tumor angiogenesis, thereby reducing the efficacy of TACE, which is consistent with our immune infiltration analysis.³⁹

Ischemia and hypoxia caused by hepatic artery embolization are important factors contributing to the efficacy of TACE. Studies have shown that TACE can increase the expression of HIF-1α pathway components.^40,41 In our analysis, we found that patients in the LE group already exhibited high expression of hypoxia signaling pathway components before the procedure. Considering this and the significant enrichment of the glycolysis pathway (which is known to increase under hypoxic conditions).⁴² In the LE group, we reasonably speculate that the hypoxic microenvironment may influence the efficacy of TACE, and numerous studies have indicated that hypoxia can induce resistance to stress and apoptosis in HCC.^43,44 Additionally, MMPs and glucose metabolism-related enzymes are both involved in hypoxia-induced effects on HCC.⁴³ Considering the impact of GNMT on glucose metabolism mentioned earlier, we propose that the loss of GNMT and the overexpression of MMP12 may create a hypoxic tumor microenvironment, ultimately diminishing the effectiveness of TACE. This is a hypothesis we are currently dedicated to investigating and validating.

Our signature demonstrates superior accuracy in predicting the effectiveness of TACE, with AUC values exceeding 0.80, in contrast to other signatures. Additionally, the 2-gene signature incurs lower detection costs compared to multi-gene signatures. Most notably, our study includes biopsy pathology samples exclusively from patients treated solely with TACE as external validation, and the results are validated at the protein level, enhancing the overall reliability of the model. However, this study has several limitations. Firstly, the data in the external validation set were retrospectively collected, introducing potential confounding factors due to the non-randomized design. Nevertheless, in the comparison of OS and PFS, we employed the median of the effectiveness score as the grouping threshold to ensure better comparability of other variables between the two groups (Table 2). Additionally, to minimize the impact of other treatment modalities on study results, we exclusively selected patients who solely underwent TACE treatment as an external validation set. Nonetheless, because the diagnosis of HCC is not solely dependent on pathology, resulting in the limited samples of puncture biopsy pathology, we were constrained to analyze using a restricted number of samples. Moreover, in the section on immune infiltration analysis, insufficient tissue samples hindered further validation. Despite these constraints, the study incorporated three independent datasets and underwent multiple queues for mutual validation, thereby enhancing result reliability. Finally, cellular experiments were not included in this study, but this is an area we are actively exploring.

Conclusion

In summary, we have developed a 2-gene signature consisting of GNMT and MMP12 that can accurately predict the effectiveness of TACE. Additionally, we proposed that the hypoxic microenvironment may be a critical factor influencing TACE efficacy. Our study has identified novel targets and provided new insights for future research, advancing TACE treatment toward the era of precision medicine.