Matrix metalloproteinase 1 is a poor prognostic biomarker for patients with hepatocellular carcinoma

Hepatocellular carcinoma (HCC) remains an incurable malignancy despite the treatment methods being continually updated. Matrix metalloproteinases (MMPs) promote the progression of HCC; however, no consensus exists on which MMP plays the predominant role in HCCs. In the present study, we analyzed differentially expressed genes in HCCs, especially MMPs, compared with adjacent tissues using the Cancer Genome Atlas database. The KEGG enrichment pathway using differentially expressed genes included extracellular matrix–receptor interaction, which was correlated with MMPs. We found that among the MMP family, only MMP1, MMP3, MMP8, MMP9, MMP11, MMP12, MMP14, MMP15, MMP20, MMP21, and MMP24 significantly increased in HCCs compared with adjacent tissues. Crucially, survival and univariate analyses indicated that only MMPs 1, 9, 12, and 14 predict poor overall survival; however, multivariate Cox analysis and a nomogram demonstrated that only MMP1 is a poor prognostic biomarker for HCCs. In addition, we observed significant enrichment of uncharacterized cells but decreased macrophages in HCC tissues. Consistent with decreased macrophages in HCCs, MMP1 was negatively associated with macrophages but positively correlated with uncharacterized cells, indicating that the main producer of MMP1 is uncharacterized cells. Furthermore, MMP1 expression was negatively correlated with immune responses of HCCs. Taken together, our findings indicated that MMP1 is a poor and predominant prognostic biomarker for patients with HCC and that anti-MMP1 may be a novel therapy that is worth studying in depth.


Introduction
Liver cancer is the sixth most commonly diagnosed cancer and the fourth leading cause of cancer-related death worldwide [1,2].Hepatocellular carcinoma (HCC), the principal histologic type of liver cancer, represents approximately 75-85% of primary liver cancers and constitutes a major health problem worldwide [1,2].Despite relevant advances in its systemic treatment approaches, HCC remains a deadly disease that normally occurs against a background of chronic liver disease (e.g., cirrhosis) [3].Therefore, an urgent need exists to explore new biomarkers to have a clinically meaningful impact in the screening of HCC patients.
MMPs can promote cancer progression by increasing cancer cell growth, migration, invasion, metastasis, and angiogenesis [6].MMPs exert these effects by cleaving a diverse group of substrates, which include not only structural components of the extracellular matrix (ECM), but also growth-factor-binding proteins, growth factor precursors, receptor tyrosine kinases, cell adhesion molecules, and other proteinases [6].Direct evidence for the role of MMPs in tumor progression has been provided by xenograft experiments that have used cancer cells with decreased and increased expression levels of MMPs or tissue inhibitors of metalloproteinases (TIMPs), and also by carcinogenesis experiments with mice that either lacked a specific MMP or TIMP-1 or had organ-specific MMP or TIMP-1 overexpression [7][8][9][10][11][12][13]. MMPs are upregulated in almost every type of human cancer, including HCCs, and their expression is often associated with poor survival [6,[14][15][16][17][18][19].However, there is no consensus on which MMP plays the predominant role in HCCs, which limits their clinical predictive value for the prognosis of patients with HCC.Therefore, the prognostic role of MMP members in HCCs still needs to be fully uncovered.
In the present study, we characterized the gene expression profiles from the Cancer Genome Atlas (TCGA) database (https:// portal.gdc.com).We analyzed differentially expressed genes in HCCs, especially MMPs, compared with adjacent tissues using the TCGA database.Among the MMP family, only MMP1, -3, -8, -9, -11, -12, -14, -15, -20, -21, and -24 were significantly increased in HCCs compared with adjacent tissues.Crucially, survival analysis, univariate analysis, multivariate Cox analysis, and a nomogram demonstrated that only MMP1 is a poor prognostic biomarker for HCCs.In addition, we observed significant enrichment of uncharacterized cells but decreased macrophages in HCC tissues.Consistent with decreased macrophages in HCCs, MMP1 was negatively associated with macrophages but positively correlated with uncharacterized cells, indicating that the main producer of MMP1 is uncharacterized cells.Furthermore, MMP1 expression was negatively associated with immune responses in HCCs.Collectively, MMP1 is a poor and predominant prognostic biomarker for patients with HCC and anti-MMP1 may be a novel therapy that is worth studying in depth.

Gene expression data sets and functional enrichment
The RNA-seq transcriptome data of 371 HCC and 50 adjacent tissues were obtained from the TCGA database (https:// www.cancer.gov/ about-nci/ organ izati on/ ccg/ resea rch/ struc tural-genom ics/ tcga).The analysis method of RNA-seq is consistent with our previous described method [20][21][22].The Limma package (version: 3.40.2) for R was used to study the differential expression of mRNAs [23].The adjusted P value was analyzed to correct for false positive results in the database.The thresholds for screening the differential expression of mRNAs were defined as follows: adjusted P < 0.05 and Log2 (Fold Change) > 1.585 or Log2 (Fold Change) < − 1.585.To further confirm the underlying function of potential targets, the data were analyzed through functional enrichment.Kyoto Encyclopedia of Genes and Genomes (KEGG) Enrichment Analysis is a practical resource for the analytical study of gene functions and associated high-level genome functional information.To better understand the carcinogenesis of mRNA, the clusterProfiler package (version: 3.18.0)for R was employed to analyze the enrichment of the KEGG pathway [23].And GSEA enrichment was analyzed by R software clusterProfiler package and enrichplot package as our previous described method [20].

Immune infiltration estimations
To ensure reliable immune infiltration estimations, we used the EPIC score quantified by immunedeconv, an R package.The analysis method and R package were implemented using R (2020) version 4.0.3(R Foundation for Statistical Computing, Vienna, Austria) and the software packages ggplot2 and pheatmap [24][25][26].

Analysis of differential MMP family expression
For the 371 patients from the TCGA database, tumoral RNA-seq data were downloaded from the Genomic Data Commons data portal (TCGA), and 50 of the tumors also had mRNA expression data for paired normal tissue samples.All normal tissue sample data were obtained from GTEx V8 release version (https:// gtexp ortal.org/ home/ datas ets).Complete descriptions of the donor genders, multiple ethnicity groups, wide age range, biospecimen procurement methods, and sample fixation were provided in GTEx official annotations.Statistical analyses were performed using R v4.0.3 [27].A P value of < 0.05 was considered statistically significant.
Construction and validation of a risk model and Kaplan-Meier survival analysis.
In addition, we conducted Kaplan-Meier survival analysis with a log-rank test to compare the survival difference between the aforementioned two groups or more.timeROC (v 0.4) analysis was performed to compare the predictive accuracy of each gene and risk score.The least absolute shrinkage and selection operator (LASSO) regression algorithm for feature selection, using tenfold cross-validation, the above analysis used the R software package glmnet (v 4.1-1).For the Kaplan-Meier curves, P values and hazard ratios (HRs) with 95% confidence intervals (CIs) were generated through log-rank tests and univariate Cox proportional hazards regression.All of the aforementioned analytical methods and R packages were performed using R v4.0.3 [28].P < 0.05 was considered statistically significant.

Univariate analysis, multivariate Cox analysis, and the nomogram
Univariate and multivariate Cox regression analyses were performed to identify the proper terms for building the nomogram as our previous described [29].A forest was used to present the P value, HR, and 95% CI of each variable through the R package forestplot.The nomogram was developed based on the results of the multivariate Cox proportional hazards analysis to predict the X-year overall recurrence.The nomogram provided a graphical representation of the factors that can be used to calculate the risk of recurrence for an individual patient using the points associated with each risk factor through the R package rms [30][31][32].

Immunological correlation analysis
The data set comprised mRNA-seq data from TCGA tumors (see the TCGA Data Portal at https:// tcga-data.nci.nih.gov/ tcga/).A two-gene correlation map was realized using the R software package ggstatsplot, and a multi-gene correlation map was displayed using the R software package pheatmap.We used Spearman's correlation analysis to describe the correlation between quantitative variables without a normal distribution.A P value less than 0.05 was considered statistically significant.All of the aforementioned analysis methods and R packages were implemented using R v4.0.3 [26,33,34].Correlation between MMP1 expression and tumor mutation burden (TMB) and microsatellite instability (MSI) were also calculated using the same Spearman's correlation analysis.

Cancer stemness analysis by one-class logistic regression machine learning algorithm (OCLR)
The cancer stemness features of MMP1 hi and MMP1 low HCC patients were analyzed by OCLR as previously described [35].Briefly, for the Cancer Genome Atlas (TCGA) database, we download tumor RNA-seq (FPKM) from the Genomic Data Commons (GDC).Convert PFKM data to TPM and normalize the data log2 (TPM + 1), while keeping samples with clinical information recorded.We then used the OCLR algorithm constructed by Malta et al. [35] to calculate mRNAsi scores.Based on the characteristics of mRNA expression, the gene expression profile contains 11,774 genes.We use the same Spearman correlation (RNA expression data) and then subtract the minimum value and divide by the linear transformation of the maximum value mapping the dryness index to the range (0, 1).All the above analysis methods and R package were implemented by the R foundation for statistical computing (2020) version 4.0.3.

Immune checkpoint blockade (ICB) response analysis by tumor immune dysfunction and exclusion (TIDE)
Raw counts of RNA-sequencing data and corresponding clinical information from HCCs were obtained from The Cancer Genome Atlas (TCGA) dataset, in which the method of acquisition and application complied with the guidelines and policies.Potential ICB response was predicted with TIDE algorithm as previously described [36].TIDE is a computational method, which models two primary mechanisms of tumor immune evasion: the induction of T cell dysfunction in tumors with high infiltration of cytotoxic T lymphocytes (CTL) and the prevention of T cell infiltration in tumors with low CTL level.

Statistics analysis
R version 4.0.3 was used for all statistical analyses.The statistical details of all experiments are reported in the materials and methods, including statistical analysis performed and statistical significance.

Differentially expressed genes and enriched pathways in HCCs
We compared the differentially expressed genes to explore the molecular characteristics of HCCs and then retrieved the transcriptome profiling data of HCCs from the TCGA database.A total of 522 genes were distinguished as differentially expressed mRNAs, with 329 genes being upregulated and 193 genes being downregulated (Fig. 1A and B).First, we performed GSEA of differentially expressed genes between HCCs and adjacent tissues.Figure 1C indicates that the upregulated pathways include cell cycle, viral carcinogenesis, cellular senescence, MicroRNAs in cancer, and so on.Figure 1D demonstrates that the downregulated pathways include drug metabolism, metabolism of xenobiotics by cytochrome P450, bile secretion, and so on.The enriched KEGG signaling pathways were selected to demonstrate the primary biological actions of major potential mRNAs.The upregulated KEGG pathways included viral carcinogenesis, focal adhesion, and ECM-receptor interaction (Supplementary Fig. 1A).By contrast, the downregulated KEGG pathways included the PPAR signaling pathway, fatty acid degradation, and drug metabolism (Supplementary Fig. 1B).Then, we performed Gene Ontology (GO) analysis using differentially expressed genes.The upregulated GO included positive regulation of cell cycle, nuclear division, DNA replication, packaging, and DNA conformation change (Supplementary Fig. 1C).By contrast, the downregulated GO included epoxygenase P450 pathway, drug metabolic process, fatty acid metabolic process, and so on (Supplementary Fig. 1D).Collectively, these upregulated changes may play significant roles in the progression of HCCs and the downregulations indicate significantly decreased liver function.

Univariate and multivariate Cox analyses and the nomogram predicted MMP1 to be a poor prognostic factor for HCC
Subsequently, we conducted univariate and multivariate Cox analyses to evaluate the independent prognostic value of MMP1, -9, -12, and -14 in terms of the OS of patients with HCC.The univariate analysis indicated that the high groups of MMP 1, -9, -12, and -14 were significantly correlated with poor OS (Fig. 4A).Other variables related to poor survival included the p-TNM stage.Multivariate analysis revealed that the high group of MMP1 was independently associated with a significantly poorer OS of patients with HCC (Fig. 4B), which could serve as an Fig. 1 A Volcano plots of differentially expressed genes constructed using fold-change values and adjusted P values.The red point in the plot represents the over-expressed mRNAs, and the blue point indicates the down-expressed mRNAs with statistical significance.B Heat map of differentially expressed genes between HCC tissues (tumor) and adjacent tissues (control).C Upregulated GSEA enrichment pathways in HCC tissues.D Downregulated GSEA enrichment pathways in HCC tissues ◂ independent prognostic factor for HCC.To further develop a clinically applicable method that could predict the survival probability of a patient, we used a nomogram to construct a predictive model that considered clinicopathological covariates.On the basis of the univariate and multivariate analyses of OS rate, we generated a nomogram to predict the 1-, 3-, and 5-year OS rates in the discovery group using the Cox regression algorithm (Fig. 4C) and predicted the patients' odds of dying using generalized linear regression.The predictors included p-TNM stage and MMP1, which satisfied the criterion of P < 0.05 in the risk assessment.The calibration plots for the 1-, 3-, and 5-year OS rates were predicted well compared with an ideal model in the entire cohort (Fig. 4D).

Construction and validation of a risk model according to MMP1, -9, -12, and -14 in HCCs
We screened MMP-associated prognostics from the TCGA training set using univariate Cox regression analysis.Four MMPs in the TCGA data set were significantly correlated with OS.To investigate the prognosis of these four genes on HCC samples, we constructed a prognostic model based on the least absolute shrinkage and selection operator (LASSO) for dimensionality reduction.The dashed perpendicular line illustrates the first-rank value of log λ with minimum segment likelihood bias.Hence, MMP1 and MMP9 were selected for the subsequent multivariate analysis (Fig. 5A and B).HCC samples were categorized into low-and high-risk groups based on the median value of the MMP1 and MMP9 grades.The distribution of risk grades between the low-and high-risk groups, the survival status and survival time of patients in the two different risk groups, as well as the relative expression standards of MMP1 and MMP9 are presented in Fig. 5C.The survival analysis demonstrated that the OS of the low-risk group was longer than that of the high-risk group (P = 0.00111; Fig. 5D).The AUC was 0.693 at 1 year, 0.676 at 3 years, and 0.635 at 5 years, indicating a high predictive value of MMP1 in combination with MMP9 (Fig. 5E).

MMP1 is negatively associated with macrophages but positively correlated with uncharacterized cells in HCCs
Immune cells are a major component of the tumor microenvironment [40].To better understand the characterization of the TIICs signature in HCCs, we conducted immune infiltration estimations using EPIC scores from the TCGA database (available online: http:// tcgad ata.nci.nih.gov/ tcga/ tcgaA bout.jsp).In total, 50 adjacent and 371 HCC samples were included in this study.Figure 6A and B demonstrates that among the immune cells, macrophages and uncharacterized cells in HCC tissues were significantly different compared with adjacent tissues.Furthermore, among HCCs, macrophage infiltration in each sample was distinct, highlighting the heterogeneity in patients with HCC.We further analyzed the effects of MMP1, -9, -12, and -14 on the immune microenvironment in patients with HCC.The four characteristic MMPs exhibited strong correlations with macrophages and uncharacterized cells (Fig. 6C).However, the correlations of the four MMPs with B cells, T cells, NK cells, and endothelial cells were not statistically significant (Fig. 6C).Moreover, we further observed that the four MMPs were negatively associated with macrophages but positively correlated with uncharacterized cells in HCCs (Fig. 6C).

High expression of MMP1 is associated with low responses to immune checkpoint blockade (ICB) via increased stemness of HCCs
In recent years, immune checkpoint inhibitors have become indispensable immunotherapy for HCC patients [41].We  then analyzed the ICB responses of MMP1 hi and MMP1 low HCC patients by tumor immune dysfunction and exclusion (TIDE) scores.Figure 7A indicates that MMP1hi patients had a higher TIDE score, which represents lower responses to ICBs.Immune checkpoint molecules are highly associated with ICB responses, we next analyzed the expression of immune checkpoint molecules in the MMP1 hi and MMP1 low HCC patients.Figure 7B demonstrates that the expression of PDCD1LG1, CTLA4, TIM3, LAG3, PDCD1, and TIGIT but not PDCD1LG2 and SIGLEC15 significantly increased.Tumor mutation burden (TMB) and microsatellite instability (MSI) are biomarkers that can help predict a patient's response to immunotherapy [42,43].We further tested the association of TMB and MSI with the expression of MMP1 in HCCs. Figure 7C and D reveals that both TMB and MSI are negatively associated with the expression of MMP1 in HCCs.Previous studies have revealed that cancer progression involves the gradual loss of a differentiated phenotype and acquisition of progenitor and stem cell-like features [35].Significantly, cancer stem cells are not sensitive to many treatment approaches including immunotherapy [44].We then tested the stemness features of MMP1hi and MMP1low HCC patients by one-class logistic regression (OCLR) machine learning algorithm as previously described [35].Figure 7E reveals that MMP1hi HCC patients had higher OCLR scores, suggesting higher stemness features of MMP1hi HCCs.Together, these data indicated that MMP1 hi HCC patients are not sensitive to ICBs.
High expression of MMP1 is associated with poor prognosis of HCC patients.Having shown that MMP1 expression can be a predictive marker in HCCs.We then analyzed the expression of MMP1 during the progression of HCCs. Figure 8A shows that MMP1 expression is higher in stage III patients and in stage II patients than in stage I patients.Although there was no statistically significant difference between the expression of MMP1 in stage III patients and the expression in stage II patients, the expression of MMP1 in patients with stage III had a tendency to increase compared with patients with stage II.Next, we compared the clinicopathological characteristics of patients with high and low expression of MMP1.Table 1 indicates that MMP1hi group of patients exhibited more dead patients, higher TNM stage and Grade stage.We then compared the transcriptomes of patients with high and low MMP1 expression.Compared with the MMP1 low expression group, the MMP1 high expression group upregulated 48 genes and downregulated 66 genes (Fig. 8B and C).The upregulated KEGG pathways included platinum drug resistance, cell cycle, and transcription misregulation in cancer (Fig. 8D), and the downregulated KEGG pathway included PPAR signaling pathway, metabolism of xenobiotics by cytochrome P450, fatty acid degradation, and bile secretion (Fig. 8E).Similarly, the upregulated GO terms included nuclear division, chromosome segregation, and mitotic spindle assembly (Fig. 8F) and the downregulated GO terms included epoxygenase P450 pathway, fatty acid metabolic process, and alcohol metabolic process (Fig. 8G).Collectively, these upregulated changes may play significant roles in the progression of HCC patients and the downregulations indicate significantly decreased liver function.Finally, we selected some upregulated genes in the MMP1hi group, which is significant for the progression of HCC.These genes included CDK1, MKI67, MMP9, MMP12, CD24, and AFP (Fig. 8H).Taken together, high expression of MMP1 is associated with a poor prognosis of HCC patients.

Discussion
HCC, an aggressive cancer, has a poor long-term prognosis and ranks as the third leading cause of cancer-related death [17].Despite advances in targeted therapy and immunotherapy, HCC remains an incurable malignant solid tumor.Mounting evidence exists for the functions of MMPs, and furthermore, they are the pivotal mediators in microenvironment alteration-determined tumorigenesis [45,46].However, the relationships between MMPs and HCCs remain unclear.In this study, we aimed to systematically analyze the relationship between MMPs and prognosis in patients with HCC and then to establish a more precise prognostic biomarker associated with the MMP family.In the liver, all hepatic cells, such as hepatocytes, hepatic stellate cells (HSCs), hepatic macrophages (including resident KCs and infiltrated monocyte-derived macrophages [MoMFs]), and infiltrated leukocytes, are able to produce MMPs; however, HSCs are the major producers among them [47,48].In addition, MMPs are also produced by tumor cells in HCC [49,50].Through analyzing immune cell infiltration, we found that macrophages significantly decreased in HCC; however, uncharacterized cells significantly increased.The upregulated KEGG pathways in tumor tissues include ECM receptor interaction and focal adhesion, which are correlated to MMPs.Previous studies have suggested that MMPs may play significant roles in HCCs [49].However, no consensus exists regarding which one plays the predominant prognostic role in HCCs or which cell type is the main producer that plays a significant role in HCCs.Because macrophages decreased, while uncharacterized cells increased in HCCs.We speculated that uncharacterized cells, which mainly included hepatocytes, HSCs, and tumor cells in patients with HCC are the main producers of MMPs.This led us to hypothesize that MMPs expressed by tumor cells might be a superior biomarker for HCC prognosis.
A previous study reported that MMP9 is a significant biomarker in patients with HCC.[51] While the survival and univariate analyses revealed that high expression of MMP9 was associated with a poorer OS than low expression of MMP9, the multivariate Cox analysis indicated no significant difference between high and low expression of MMP9.Consistent with the previous study, overexpression of MMP12 mRNA was significantly associated with a poor OS of patients with HCC after hepatectomy according to the survival and univariate analyses [51].By contrast, multivariate Cox proportional hazard regression analysis revealed that MMP12 mRNA was not an independent factor for OS [51].Furthermore, immunohistochemical staining revealed that MMP14 expression in HCC tissues was significantly higher than in adjacent tissues [52].In terms of survival, MMP14 expression was negatively associated with both OS and DFS in univariate analyses, while it remained statistically significant for DFS but not OS in a multivariate Cox regression test [52].In the present study, MMP14 was also found to be increased in patients with HCC, exhibiting a poor OS for HCCs in the survival and univariate analyses.By contrast, the multivariate Cox regression test indicated no significant difference between high and low expression of MMP14.
Previous studies have demonstrated that MMP1 is considered a putative and predictive marker for many cancer types, including breast, oral, ovarian, and colorectal cancer [53][54][55][56].In addition, several studies have indicated that inhibiting MMP1 significantly suppresses HCC progression [57][58][59].However, the prognostic role of MMP1 in human HCCs has not been previously reported.In the present study, we found that MMP1 increased in patients with HCC.Within the HCCs, the expression of MMP1 displayed a heterogenous pattern.Using survival analysis, univariate analysis, multivariate Cox regression, and a nomogram, we found that high expression of MMP1 was associated with a poor OS in patients with HCC, which suggested that MMP1 may play more crucial roles than MMP 9, -12, and -14 in HCC progression.This led us to believe that MMP1 might be a superior biomarker for HCC prognosis.Considering this, we further examined the prognostic value of MMP1 expression in multiple cancers.Noteworthily, MMP1 expression also predicted survival outcomes in other cancers, including PAAD, LIHC, LUAD, CESC, KICH, KIRP, SARC, and UVM.
We also constructed and validated a risk model according to MMP1, -9, -12, and -14 in HCCs.Using the first-rank value of log λ with minimum segment likelihood bias analysis, MMP1 and MMP9 were selected for the subsequent multivariate analysis.Our results demonstrated that the OS of the low-risk group (low expression of MMP1 and MMP9) was longer than that of the high-risk group (high expression of MMP1 and MMP9).Thus, the multivariate Cox analysis indicated no significant difference between high and low expression of MMP9; however, MMP9 in combination with MMP1 can be a poor prognostic biomarker for HCCs.
MMPs can be produced by many cell types in HCCs, including tumor cells, macrophages, hepatocytes, and HSCs [47,48,57,58].We used Spearman's correlation analysis to describe the correlation between MMPs and immune cells.Remarkably, MMP1, -9, -12, and -14 were all negatively correlated with macrophages, whereas they were positively associated with uncharacterized cells in tumor tissues.The significantly decreased macrophages and increased uncharacterized cells in HCCs also confirmed these data.Uncharacterized cells may include tumor cells, hepatocytes, and HSCs.Given that HCCs have predominant tumor cells, we can speculate that MMP1, -9, -12, and -14 may be produced by tumor cells in HCCs, but further studies are required to clarify the producers of the four MMPs.
Immunotherapy has become an important approach in the treatment of HCC patients.However, only a minority of HCCs are sensitive to immunotherapy [60].Research on biomarkers of a response or primary resistance to immunotherapy is advancing.TMB, MSI, immune checkpoint molecules such as PD-L1, PD-1, and stemness features are previously used to predict the clinical responses of HCCs to immunotherapy [35,41,61].In the present study, we found that MMP1 hi HCC patients had a poorer response to immunotherapy compared with MMP1 low HCC patients.Consistent with this, TMB, and MSI are higher in MMP1 low HCC patients, which can partially explain the higher immune responses to immunotherapy.Previous studies have reported that the expression of PD-1 and PD-L1 are positively correlated with anti-PD1 or PD-L1 immunotherapy [62,63].However, in the present study, we found that MMP1 hi HCC patients had a higher expression of PD-1 and PD-L1.Given that the expression of PD-1 and PD-L1 can be induced by many factors, such as inflammatory status [64], therefore, we speculate the expression of PD-1 and PD-L1 may not a good indicator to predict the clinical responses to immunotherapy in HCCs.Cancer cells with higher stemness features are not sensitive to immunotherapy [35].Consistent with this, MMP1hi HCC patients showed higher stemness features, suggesting lower sensitivity to immunotherapy.Additionally, it has been demonstrated recently that MMPs modulate the microenvironment of tumors, thus promoting tumorigenesis.The dysregulation of MMPs/TIMPs leads to serious pathological conditions related to inflammation, uncontrolled cell growth, degradation of the extracellular matrix, increased cell migration, cell death resistance, replicative immortality, and the establishment of metastatic niche at secondary sites [65].Therefore, the MMPs/TIMPs imbalance in the tumor microenvironment in HCC may also affect the sensitivity to immunotherapy.Further studies should be performed to examine the roles of MMPs/TIMPs in HCC patients.Collectively, MMP1 expression can be a predictor for HCC patients with immunotherapy.
Our study has several limitations that require consideration.Although all analyses demonstrated that MMP1 is a superior biomarker for HCC progression, more studies are required to identify the roles and mechanisms of MMP1 during HCC progression.In addition, whether anti-MMP1 treatment alone or in combination with existing immunotherapies can yield significantly enhanced anti-tumor activity in vitro and in vivo must also be confirmed.Nevertheless, the expression of MMP1 is still a good predictor of the prognosis of patients with HCC, especially the predictive role for HCCs with immunotherapy.

Fig. 3 A
Fig. 3 A Prognostic analysis of gene signature in the TCGA set.The dotted line represents the median risk score and divides the patients into low-risk and high-risk groups.The curve of risk score.Survival status of the patients.More dead patients correspond to a higher risk score.Heat map of the expression profiles of MMP1, -9, -12, and -14

Fig. 4
Fig. 4 Hazard ratio and P value of constituents involved in univariate A and multivariate Cox regression B and some parameters of the MMP1, -9, -12, and -14 genes.C Nomogram for predicting the 1-, 2-, and 3-year overall survival of patients with HCC.D Calibration curve

Fig. 5 A
Fig. 5 A Coefficients of MMP1, -9, -12, and -14 depicted by the Lambda parameter.B Partial likelihood deviance versus log (λ) was drawn using the LASSO Cox regression model.C The dotted line represents the median risk score and divides the patients into low-risk and high-risk groups.The curve of risk score.Survival status of the

Fig. 6 A
Fig. 6 A Immune cell EPIC score heat map between HCC tissues (tumor) and adjacent tissues (control), where different colors represent the expression trend in different samples.B The EPIC score of

Fig. 7 A
Fig. 7 A The distribution of TIDE scores in MMP1 hi and MMP1 low groups in the prediction results.B The expression of immune checkpoint molecules in MMP1 hi and MMP1 low groups.C The correlation between TMB scores and MMP1 expression in HCC patients.D

Fig. 8 A
Fig. 8 A The MMP1 expression during HCC progression (Normal, Stage I, Stage II, and Stage III).B Volcano plots of differentially expressed genes between MMP1hi and MMP1low groups constructed using fold-change values and adjusted P values.The red point in the plot represents the over-expressed mRNAs, and the blue point indicates the down-expressed mRNAs with statistical significance.C Heat map of differentially expressed genes between MMP1hi HCC tissues andMMP1low HCC tissues.D Upregulated KEGG pathways in MMP1hi HCC tissues.E Downregulated KEGG pathways in MMP1hi HCC tissues.F Upregulated GO terms in MMP1hi HCC tissues.G Downregulated GO terms in MMP1hi HCC tissues.(H) The mRNA expression of CDK1, MKI67, MMP9, MMP12, CD24, and AFP between MMP1hi HCC tissues andMMP1low HCC tissues ◂

Fig. 9 A
Fig. 9 A The mRNA expression of MMP1 in multiple cancers.B Kaplan-Meier survival analysis of MMP1 in the low-and high-risk groups in multiple cancers