Next Article in Journal
Unsupervised Deep Learning Registration of Uterine Cervix Sequence Images
Previous Article in Journal
Chronic Lymphocytic Leukemia Progression Diagnosis with Intrinsic Cellular Patterns via Unsupervised Clustering
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

The Immunological Contribution of a Novel Metabolism-Related Signature to the Prognosis and Anti-Tumor Immunity in Cervical Cancer

1
Department of Obstetrics and Gynecology, Shanghai General Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai 201620, China
2
Reproductive Medicine Center, Department of Obstetrics and Gynecology, Shanghai General Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai 201620, China
*
Authors to whom correspondence should be addressed.
These authors contributed equally to this work.
Cancers 2022, 14(10), 2399; https://doi.org/10.3390/cancers14102399
Submission received: 4 April 2022 / Revised: 7 May 2022 / Accepted: 9 May 2022 / Published: 13 May 2022

Abstract

:

Simple Summary

Cervical cancer is the most commonly diagnosed gynecological malignant carcinoma worldwide. It is crucial to develop reliable prognostic models to predict clinical outcomes and identify patients who will benefit from different treatment strategies. In current study, we identified a reliable metabolism-related signature composed of ALOX12B, CA9, FAR2, F5 and TDO2 for the prognosis and anti-tumor immunity in cervical cancer. Patients with high-risk score underwent apparently worse prognosis and displayed lower infiltration of tumor infiltrating lymphocytes. Additionally, the metabolism-based risk score could also predict the prognosis of patients with cervical cancer based on the expression of immune checkpoints. Since this risk score signature achieves a good performance in predicting clinical outcome, we genuinely expect that our study could provide an effective prognostic tool for guidance of personalized treatment of cervical cancer patients.

Abstract

Cervical cancer is the most frequently diagnosed malignancy in the female reproductive system. Conventional stratification of patients based on clinicopathological characters has gradually been outpaced by a molecular profiling strategy. Our study aimed to identify a reliable metabolism-related predictive signature for the prognosis and anti-tumor immunity in cervical cancer. In this study, we extracted five metabolism-related hub genes, including ALOX12B, CA9, FAR2, F5 and TDO2, for the establishment of the risk score model. The Kaplan-Meier curve suggested that patients with a high-risk score apparently had a worse prognosis in the cervical cancer training cohort (TCGA, n = 304, p < 0.0001), validation cohort (GSE44001, n = 300, p = 0.0059) and pan-cancer cohorts (including nine TCGA tumors). Using a gene set enrichment analysis (GSEA), we observed that the model was correlated with various immune-regulation-related pathways. Furthermore, pan-cancer cohorts and immunohistochemical analysis showed that the infiltration of tumor infiltrating lymphocytes (TILs) was lower in the high-score group. Additionally, the model could also predict the prognosis of patients with cervical cancer based on the expression of immune checkpoints (ICPs) in both the discovery and validation cohorts. Our study established and validated a metabolism-related prognostic model, which might improve the accuracy of predicting the clinical outcome of patients with cervical cancer and provide guidance for personalized treatment.

1. Background

Cervical cancer is the most commonly diagnosed gynecological malignant carcinoma and accounts for an estimated 604,000 new cases and 342,000 deaths annually worldwide [1]. During the past few decades, although cervical cancer screening programs and comprehensive treatment strategies, including emerging anti-tumor immunotherapy, have reduced the incidence and mortality rates in most areas of the world, the prognosis of advanced cervical cancer patients is still not improved. Therefore, it is crucial to develop reliable prognostic models to predict clinical outcomes and identify patients who will benefit from different treatment strategies.
Metabolic reprogramming in cells and changes in energy metabolism levels have been identified as an emerging hallmark of cancer [2,3]. Increased aerobic glycolysis, fatty acid (FA) metabolism and glutamine decomposition contribute to malignant transformation, the invasion-metastasis cascade, tumor microenvironment (TME) stress and the treatment resistance of cancers [4,5,6,7]. Previous studies have shown that the Warburg effect and mitochondrial dysfunction favored the metabolic adaptation and survival of cervical cancer cells [8], and that reprogramming of fatty acid metabolism was associated with lymph node metastasis of cervical cancer [9]. Meanwhile, metabolic profiles could also distinguish cervical precancerous lesions from the normal cervical epithelium [10].
It is increasingly clear that crosstalk between abnormal metabolism and immune escape assumes a key role in the process of tumor progression [11,12]. For instance, cervical cancer cells can secrete lactate to convert the phenotype of tumor macrophages [13]. Recently, the application of bioinformatic analysis in predicting the prognosis and treatment response of patients with a malignant tumor has attracted rising attention. Within the context, some researchers have combined metabolomics with genomics to demonstrate the relationship between metabolism and immune infiltration [14,15,16]. However, to date, the association between the metabolism-based risk score model and TME landscape in cervical cancer remains uncharted territory.
In the current study, we established a novel prognostic metabolism-related risk score based on cancer genomics, bioinformatics and immunohistochemical analysis. The association between the risk score and infiltration of tumor infiltrating lymphocytes (TILs) was explored in a cervical cancer training cohort, a validation cohort and pan-cancer cohorts. Moreover, we also utilized the risk score to predict the prognosis of cervical cancer patients in the context of different expression of immune checkpoints.

2. Materials and Methods

2.1. Data Retrieval and Identification of Differentially Expressed Metabolism-Related Genes

The metabolism-related gene sets were downloaded from the Molecular Signature Database (MSigDB) via the Gene Set Enrichment Analysis tool (GSEA, http://software.broadinstitute.org/gsea/index.jsp (accessed on 3 May 2021)). The RNA sequencing transcriptomics data and corresponding clinical information were retrieved for a total of 304 cervical cancer tissues and another 32 cancer types from the TCGA database (https://tcga-data.nci.nih.gov/tcga/ (accessed on 3 November 2020)). We excluded samples whose overall survival (OS) or survival status were not available. We also obtained GTF files from Ensembl (http://asia.ensembl.org (accessed on 3 November 2020)) for annotation of the mRNA. Besides, 24 normal and 28 cervical cancer samples (GSE63514) based on the GPL570 (Affymetrix Human Genome U133 Plus 2.0 Array) platform, as well as 300 early cervical cancer tissues (GSE44001) based on the GPL14951 (Illumina HumanHT-12 WG-DASL V4.0 R2 expression beadchip) platform were, respectively, collected from the GEO database (http://www.ncbi.nlm.nih.gov/geo (accessed on 3 May 2021)).
The differentially expressed genes (DEGs) were identified using the R package limma. The volcano plot of DEGs was plotted and visualized with the R package ggplot2. Meanwhile, the overlapping genes of DEGs and metabolism-related gene sets were presented using the Venn diagram online tool (http://www.bioinformatics.com.cn/static/others/jvenn/example.html (accessed on 4 May 2021)).

2.2. Functional Enrichment Analysis

Gene ontology (GO) analysis was performed using the R package clusterProfiler based on the DEGs between high- and low-score groups to speculate the possible function terms. The p values were adjusted via the BH method and an adjusted p value < 0.05 was considered statistically significant. Herein, we displayed some significant enrichment outcomes in the aspect of the biological process (BP). Furthermore, we also evaluate the enrichment levels of 50 hallmark pathways using R packages GSVA and msigdbr. Based on the GSVA score, we conducted the differential analysis for these hallmark pathways with the R package limma between high- and low-score groups. In addition, we also resorted to the online tool EMTome (www.emtome.org (accessed on 19 May 2021)) to obtain the enriched networks of five metabolism-related hub genes.

2.3. Construction of the Metabolism-Related Risk Score Signature

The association of the 211 overlapped genes with clinical outcomes of cervical cancer patients was analyzed via univariate Cox regression analysis. Those genes with a p value less than 0.05 were selected for further multivariate Cox regression analysis. In this manner, we identified five hub genes that had independent prognostic values. In this step, the R packages plyr and survival were applied. Then we developed the risk score signature consisting of the five hub genes, and the total risk score of this biosignature was calculated as the following formula: Risk Score = i = 1 N β i G e n e i , where N is the number of hub metabolism-related genes, βi refers to the regression coefficient and Genei represents the expression level of each gene identified by multivariate Cox analysis. The cut-off for this risk model was defined as its median value. Subsequently, the patients were classified into high-risk and low-risk groups according to the median threshold. The ROC curve was plotted via the R package timeROC to assess the predictive potential of this signature for overall survival. The Kaplan-Meier survival curve was also conducted to evaluate the difference in patient prognosis of two subgroups using the R packages survival and survminer.

2.4. Predictive Power of This Signature in the Validation Cohort and Pan-Cancer Cohorts

The predictive performance of the risk score derived from the set of five metabolic genes was evaluated in the training cohort, as well as the validation cohorts. The multivariate Cox regression analysis demonstrated that the risk score served as an independent prognostic factor and, subsequently, we constructed an OS-related nomogram with this score using the R package rms. Furthermore, the C-index of this model was measured and the result suggested that the signature had reliable predictive power. Moreover, the calibration plots of this nomogram were presented to measure its predictive accuracy in comparison to the actual curve of the survival time.
The expression profile of early cervical cancer patients in GSE44001 was retrieved as a testing cohort to validate the model externally. The Kaplan-Meier analysis was used to explore the universality of this signature, while the ROC curve was also plotted via the R package timeROC. Moreover, we validated the predictive power of this model in other cancer types. We performed a pan-cancer analysis in 33 cancer types from the TCGA project by univariate Cox regression analysis, taking advantage of this risk score. The R packages survival and plyr were applied in this step.

2.5. Measurement of Tumor-Infiltrating Immune Cells and the Potential Response of Patients for Immunotherapy

An integrated list of representative marker genes of tumor-infiltrating immune cell types was acquired from Charoentong’s research, which involved a total of 366 microarrays of immune cells summarized from 37 previous studies [17]. To evaluate the infiltration of immune cells, we conducted the single-sample gene set enrichment analysis (ssGSEA) algorithm using these marker gene sets of different immune cells with the R package GSVA, which could measure the normalized enrichment score of multiple immune cell types. The tumor abundance of these immune cells between high-and low-risk groups was also plotted using the ggplot2 package. In the meantime, we downloaded the enriched results from ImmuCellAI, which estimated the abundance of 24 immune cells in the TCGA-CESC cohort [18]. The immune network of 24 ImmuCellAI cell types in the TCGA cohort was illustrated by the R packages reshape2, corrplot and igraph. In this Circos plot, we shed light on the correlation among these cell types and their survival impact by Spearman correlation analyses and univariate Cox regression analyses. Correlation between the established signature and these immune cells was also calculated by Spearman analyses using the R software.
Kaplan-Meier analysis was carried out to clarify the relationship between the prognosis of patients with similar expression levels of immune checkpoints and the established signature. The R packages survival and survminer were utilized to draw the plot and conduct multiple comparisons between different survival curves.

2.6. Tissue Microarray and Immunohistochemical (IHC) Staining

The study was approved by the Institutional Ethics Committee of Shanghai General Hospital. Cervical cancer (n = 32) and normal cervix (n = 32) tissue microarrays were purchased from Shanghai Zuocheng Biotech (Shanghai, China). We performed the IHC analysis as previously described [19]. The antibodies used in the study were listed as follows: anti-ALOX12B rabbit polyclonal antibody (NBP1-89409, Novus Biologicals, Colorado, CO, USA), anti-FAR2 rabbit polyclonal antibody (NBP1-90435, Novus Biologicals, Colorado, CO, USA), anti-F5 rabbit polyclonal antibody (20963-1-AP, Proteintech, Wuhan, China), anti-CA9 rabbit polyclonal antibody (11071-1-AP, Proteintech, Wuhan, China), anti-TDO2 rabbit polyclonal antibody (15880-1-AP, Proteintech, Wuhan, China), anti-CD4 mouse monoclonal antibody (67786-1-Ig, Proteintech, Wuhan, China), anti-CD8 mouse monoclonal antibody (66868-1-Ig, Proteintech, Wuhan, China), anti-CD57 rabbit polyclonal antibody (19401-1-AP, Proteintech, Wuhan, China) and anti-CD68 mouse monoclonal antibody (66231-2-Ig, Proteintech, Wuhan, China). The immunoreactivity score (IRS) was used to evaluate the expression level of each protein. The staining intensity was scored as: negative = 0, weak = 1, moderate = 2 and strong = 3. The staining extent was scored as: 0 (no positive cells), 1 (≤25% positive cells), 2 (26–49% positive cells), 3 (50–74% positive cells) and 4 (≥75% positive cells). IRS = extent score × intensity score. The numbers of CD4+ cells, CD8+ cells, CD57+ cells and CD68+ cells at the tumor site were counted under five randomly selected microscopic fields.

2.7. Scoring of Immune Cell Infiltration

We performed the scoring of immune cell infiltration (CD4+ cells, CD8+ cells, CD57+ cells and CD68+ cells) according to the methods used in the previous study [20]. Briefly, stained samples were assessed and scored on a five-point scale for the infiltration level of cells into epithelial or stromal areas, with regard to the range of infiltration, as the following scale: no positive events found on slide = 1, rare positive events observed = 2, low density of infiltration = 3, medium density of infiltration = 4 and high density of infiltration = 5.

2.8. Statistical Analysis

All of the statistical process was completed within the R software (version 4.0.4). The Wilcoxon test was applied to compare between two groups, while the Kruskal-Wallis test was used to compare among more than two groups. The Kaplan-Meier plot was performed to present survival curves for different groups, and the log-rank test was employed to evaluate the significance of statistical differences. Spearman analysis was carried out to determine the correlation coefficient. For all the analyses above, a two-tailed p < 0.05 was considered as statistically significant.

3. Results

3.1. Construction of a Metabolism-Related Risk Score Signature in Cervical Cancer

The whole flow diagram of this study is presented in Figure 1A. To screen differentially expressed metabolism-related genes in cervical cancer, we extracted the expression profiling data of a cohort from GSE63514 and collected a panel of 1378 metabolism-related genes from MSigDB. As demonstrated in Figure 1B,C, the volcano plot visualized the DEGs and a total of 211 metabolism-related DEGs showed significant dysregulation in the GEO dataset. After being intersected with the expression profile of the TCGA-CESC cohort, 206 overlapping genes were taken into account for further study.
The univariate Cox regression was implemented to examine the prognostic value of metabolism-related DEGs based on the transcriptome data from the TCGA-CESC project. Those DEGs whose p < 0.05 subsequently underwent a multivariate Cox regression. Eventually, five genes (ALOX12B, CA9, F5, FAR2 and TDO2) with significant regression coefficients were incorporated to set up a prognostic risk model (Figure 1D and Supplementary Table S1). The risk score based on the set of the five metabolic genes was calculated with the following formula: risk score = −0.1083206 × ALOX12B + 0.1375361 × CA9 − 0.3202668 × F5 + 0.2241991 × FAR2 + 0.2649977 × TDO2.
The relationship between the risk score and clinicopathological characteristics was explored. We found that patients with a more advanced clinical stage displayed higher risk scores (p = 0.0244), whereas there was no significant association between the signature and the histological grade of TCGA-CESC patients (Figure 1E). There were also statistically significant differences between the high- and low-risk groups in terms of the enrichment score of gene sets correlated with hypoxia, EMT and angiogenesis (p < 0.0001, Figure 1F). Among the five genes, the expression of ALOX12B and F5 was decreased, while that of CA9, FAR2 and TDO2 was relatively elevated in cervical cancer samples, as indicated by the IHC results (Figure 1G).

3.2. Verification of the Metabolism-Related Risk Score Signature

Then we calculated the metabolism-related risk score of all the CESC patients according to the formula and re-divided them into high-risk and low-risk groups based on the median cut-off value. Evidently, the high-risk group comprised more death cases than the low-risk group, and the differential expression levels of these genes conformed to their prognostic impact (Figure 2A). A Kaplan-Meier survival analysis revealed that patients in the high-risk group exhibited significantly worse clinical outcomes (p < 0.0001, Figure 2A). The time-dependent ROC curve indicated that the risk score had a relatively higher accuracy in predicting the 3-year (AUC = 0.767, NPV = 0.846, PPV = 0.583) and 5-year OS (AUC = 0.779, NPV = 0.901, PPV = 0.514) (Figure 2A). A prognostic nomogram was also established based on the predictive indicator of the risk score, which proved to be an independent prognostic factor in CESC patients (Figure S1A). The c-index of this model was 0.742 and the calibration curve demonstrated that the nomogram-predicted overall survival was close to the actual values of the 1-, 3- and 5-year survival rates (Figure S1B). The decision curve analysis (DCA) further indicated that the risk score benefited patients with CESC in clinical practice (Figure S1C).
An independent GEO dataset (GSE44001) was used to validate the aforementioned risk score. Three hundred early cervical cancer patients were categorized into low- and high-risk groups based on the median cut-off value, and the heatmap of expression levels of the five hub genes validated their differential distribution between the two subgroups (Figure 2B). Similarly, the Kaplan-Meier survival analysis revealed that patients in the high-risk group displayed significantly worse clinical outcomes (p = 0.0059, Figure 2B), and the time-dependent ROC curve also suggested that the risk score had a relatively higher reliability in predicting the 3-year (AUC = 0.582, NPV = 0.938, PPV = 0.13) and 5-year DFS (AUC = 0.642, NPV = 0.921, PPV = 0.238) (Figure 2B).
Meanwhile, the pan-cancer validation analysis in the other 32 cancer types of the TCGA project revealed that the risk score composed of the five metabolism-related genes could successfully discriminate between patients with better or worse clinical outcomes in the other eight cancer types (Figure 2C), such as breast invasive carcinoma (BRCA), head and neck squamous cell carcinoma (HNSC), brain lower grade glioma (LGG), liver hepatocellular carcinoma (LIHC), pancreatic adenocarcinoma (PAAD), etc. Therefore, these results suggested that the metabolism-based risk score might be a promising prognostic classifier, which could be applied to predict the clinical outcome of patients with various types of malignancy.

3.3. Functional Enrichment Analysis of the Metabolism-Related Risk Score Signature

According to the estimated infiltration levels of 24 immune cells from the ImmuCellAI database, we visualized the immune landscape of patients in the TCGA-CESC cohort via an immune network (Figure 3A). Among the 24 immune cells, B cells, Tfh cells, CD4 T cells and CD8 T cells were protective factors with significance, while neutrophils and monocytes were significant risk factors. Our established signature showed a positive correlation with these risk factors but it was negatively correlated with the protective immune cells (Figure 3B).
GO and GSEA analyses were conducted to shed light on the possible function of these DEGs and the underlying mechanism. As shown in Figure 3C, significantly enriched GO terms involved ‘T cell mediated immunity’, ‘humoral immune response’ and ‘antigen processing and presentation via MHC class Ib’. In addition, GSEA analysis suggested a series of enriched hallmark pathways that were activated in the high-score group compared with the low-score group (Figure 3D). For instance, the upregulated pathways involved ‘glycolysis’, ‘DNA repair’, ‘epithelial-mesenchymal transition’ and ‘TGFβ signaling’, while the suppressed pathways included ‘inflammatory response’, ‘IL6-JAK-STAT3 signaling’ and ‘IFNα response’. These results suggested that the risk score might play a crucial role in the regulation of the immune response.
Furthermore, we also obtained the enrichment terms of the five metabolism-related hub genes via the online tool EMTome and the enriched networks are shown in Supplemental Figure S2. Moreover, we quantified the enrichment score of the metabolic-related pathways, and plotted their abundance in high- and low-risk groups (Supplemental Figure S3A,B).

3.4. Correlation between the Metabolism-Related Risk Score and Immune Landscape

The normalized enrichment scores (NES) of different immune cells in the TCGA-CESC project were calculated based on the gene expression profiles via the ssGSEA algorithm (Figure 3E). We found that specimens in the high-risk group were conferred a significantly lower infiltrating density of various immune cells, such as activated CD4 T cells, activated CD8 T cells, macrophages and CD56dim natural killer cells. In addition, we performed Spearman correlation analyses and found that the risk score and five hub genes were closely linked to the immune landscape (Supplementary Figure S4).
To verify the difference in infiltration levels of immune cells, we stained histology sections for immune markers including CD4, CD8, CD57 and CD68, and scored for the extent of infiltration both in the epithelial and stromal tumor compartments based on the density of positive staining for immune cell populations (Figure 4A). The IHC results revealed that in the high-risk group, CD4 T cells, CD8 T cells, macrophages and NK cells exhibited consistently significantly lower densities of infiltration both in the stromal and epithelial content. The outcome suggested that detectable differences in immune cell recruitment was correlated with our established risk score.
In the pan-cancer cohorts, our conclusion was also validated in a number of cancer types, as shown in Figure 4B and Supplementary Figure S5. It was obvious that patients with higher levels of our established risk score tended to exhibit lower infiltration of diverse immune cells across several cancer types.

3.5. Prognostic Value of the Established Signature in Cervical Cancer Patients for Immunotherapies

Meanwhile, we retrieved the list of genes encoding immunostimulators, immunoinhibitors, chemokines and receptors from the TISIDB website [21]. Notably, the high- and low-risk groups displayed distinct expression patterns of these immune-related molecules, which also depicted their disparity in the immune landscape (Figure 5). For instance, the expression levels of most immunostimulators, involving ICOS and TNFRSF members, were augmented in the low-risk group rather than the high-risk group (Figure 5A).
Furthermore, we stratified patients in the training cohort, as well as the validation cohort, with our established signature based on the expression levels of different immune checkpoints (ICPs), and found that a higher risk score indicated a worse prognosis in patients with similar levels of ICPs (Figure 6). Consistently, we observed that patients with a lower risk score on the basis of similar expression of ICPs experienced a more favorable clinical outcome. For example, patients with a high risk score and high PD-1 displayed a shortened OS or DFS time compared to those with a low risk score and high PD-1 (p < 0.001). Similar survival patterns could be also observed using the risk score and PD-L1, CTLA-4, CD47, CD38, CD28 and BTLA (Figure 6 and Supplementary Figure S6). These observations suggested that the risk score might serve as a predictive biomarker of treatment response to immunotherapies.

4. Discussion

Previous studies have established metabolism-related signatures for the survival prediction of patients with several types of cancer [14,22]. Herein, we developed a novel risk score signature based on the expression of five metabolism-related hub genes and evaluated its prognostic value in cervical cancer patients. This prognostic score model was confirmed in an independent cervical cancer validation cohort and pan-cancer cohorts. Furthermore, we validated the correlation between this risk score and immunologic features, which might potentially improve the accuracy of predicting the clinical outcome of patients, in combination with the conventional clinical staging. To our knowledge, this is the first study to establish a metabolism-related risk score model in cervical cancer patients, which could also be testified in various kinds of cancers simultaneously.
The TME is a complicated system which undergoes dynamic changes [23]. The complex components of the TME nurture the surrounding environment that is essential for tumor growth. More and more studies have revealed that specific metabolic patterns of the TME potentiate tumor progression or treatment resistance [24,25,26]. Under specific conditions such as hypoxia, metabolic reprogramming occurs to enhance cellular proliferation, and cancer cells have been proven to facilitate the metabolism of reactive oxygen species, lactate, lipids, glutamine and glucose, as well as amino acids [27]. For instance, rather than the oxidative phosphorylation conducted by normal cells, cancer cells are inclined to adopt lactate metabolism and glycolysis [28]. Besides, atypical lipid metabolism has also been linked to tumor recurrence and CD8+ T cell exhaustion, giving rise to post-chemotherapy evasion of immune surveillance [29,30]. Metabolic reprogramming and redox imbalances were also revealed to mediate the development and maintenance of dormant cancer cells in various malignancies, which was caused by endoplasmic reticulum (ER) stress responses and oxidative stress [31]. As a result, exploring the specific metabolic disorders and determining several metabolism-related genes that were implicated in tumorigenesis would help to predict the prognosis and therapeutic responsiveness. To this end, we developed and validated a novel metabolism-related risk score signature consisting of five metabolism-related hub genes to predict the clinical outcomes of cervical cancer patients. The AUC of the ROC curves of the TCGA cohort, based on this risk score signature model, was higher than 0.76 at the 3- and 5-year OS. Importantly, our results pinpointed that the high-risk group was characterized by a lower extent of immune infiltration. Bioinformatic enrichment analyses also identified several immune-related signaling pathways as potentially relevant pathways between the high- and low-score groups. Considering that the immunosuppressive microenvironment could facilitate tumor progression, the high-risk score patients were presumed to bear a higher tumor burden.
Biomarker-based patient stratification has gained much attention, which calls for more accurate evaluation of these molecular properties. Especially for those patients who receive immunotherapy, the comprehensive analysis of the risk score, as well as immune checkpoint expression, could hopefully foretell the reactivity of these patients and therefore screen out the appropriate patients for immunotherapy [32,33,34]. Immune checkpoints expressed on cancer cells or cancer-associated immune cells have drawn substantial attention as promising treatment targets nowadays [35]. A mounting number of studies have attempted to develop immune-based biomarker signatures to depict the survival rate and tumor progression of patients [35,36,37,38,39]. In this study, we analyzed the expression levels of these immune checkpoints in the context of cervical cancer. As Figure 6 shows, patients in the high-risk group exhibited shorter overall survival or disease-free survival on the condition that they had similar expression levels of ICPs. This exploration may potentially guide a more personalized treatment for immunotherapies.
As a member of our risk score model, tryptophan 2,3-dioxygenase 2 (TDO2) catalyzes the commitment step of the KYN metabolic process, which subsequently activates the AhR and contributes to an immunosuppressive TME, and supports the cancer immune escape [40,41]. FAR2, namely the fatty acyl-CoA reductase 2, has been found to be localized in the peroxisome and participate in the first step of wax biosynthesis [42]. Previous studies have shown that overexpression of FAR2 induced the upregulation of platelet-activating factor (PAF) and profibrotic cytokine TGF-β in a mouse mesangial cell line. Another research revealed that FAR2 mediated the de novo synthesis of PAF, a potent inflammatory mediator activating platelets, eosinophils, neutrophils and macrophages, in vitro [42,43]. As for carbonic anhydrase 9 (CA9), a pH-regulating transmembrane protein which is overexpressed in solid tumors, it was proven to repress the mitochondrial biogenesis, favor the Warburg phenotype and activate glycolysis [44,45]. CA9 equilibrates among hypoxia, iron metabolism, and redox regulation in tumor cells [46]. Beyond that, CA9 also upregulated amino acid transporters to increase the intracellular content. Ectopic expression of CA9 is a biomarker of poor prognosis in breast cancer, tongue squamous cell carcinoma, and pancreatic and lung cancers [47,48,49].
On the other side, high expression of coagulation factor 5 (F5) was found to be associated with improved overall survival of patients with breast cancers [50]. Therefore, tumor-derived F5 appears to be beneficial to patient survival, which is compatible with the tumor suppressor function proposed by our study. Consistently, the ectopic expression of F5 in breast tumors could also represent a more infiltrated microenvironment with both lymphoid and myeloid cells, such as T cells, NK cells and macrophages, according to a research conducted by Tinholt et al. [51]. ALOX12B encodes an enzyme involved in the conversion of arachidonic acid to 12R-hydroxyeicosatetraenoic acid, which has been proposed by Song et al. as a potential protective gene for the overall survival of patients with esophageal squamous cell carcinoma [52]. Egolf et al. identified ALOX12B as a driver gene of ferroptosis, a form of programmed tumor suppressive cell death featured by lipid peroxidation [53]. In the current study, we established the relationship between the metabolism-related gene signature constructed by these genes and the immune landscape of patients with cervical cancer, whereas the underlying mechanism of these hub genes still awaits to be unraveled.
Our study also has its limitations due to objective reasons. We analyzed the expression levels of metabolism-related hub genes in tumor specimens of the TCGA database, whereas we could not distinguish the cell origin of these hub genes. This means that our conclusion only represents the immunity pattern on the macro level. Further, single-cell sequencing data might better illustrate the complex association between tumor cells and the surrounding microenvironment. In the meantime, it has been a controversial issue whether the relatively small piece of tissues in TMAs could fully represent the characteristics of original tissues due to tumor heterogeneity, which might give rise to disparities in diagnosis. Since we affirmed our conclusion with an independent TMA cohort, the experimental results might also show some technical limitations to a certain extent.

5. Conclusions

In summary, our constructed metabolism-related prognostic model allows for a more accurate categorization of patients at different risk levels of cervical cancer. We also determined the infiltration of TILs, the expression pattern of immune-related molecules and the prognostic value of our signature for immunotherapies. Since this risk score signature achieves a good performance in predicting clinical outcomes, we genuinely expect that our study could provide an effective prognostic tool for the guidance of personalized treatment for cervical cancer patients.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/cancers14102399/s1; Figure S1: Internal validation of the prognostic potential of the risk score signature; Figure S2: The results of functional enrichment analyses of the five metabolism-related hub genes; Figure S3: The association between the risk score and metabolism-related pathways; Figure S4: Spearman correlation analyses of the established risk score as well as five metabolism-related hub genes (ALOX12B, FAR2, F5, CA9 and TDO2) and tumor-infiltrating cells; Figure S5: Differences in the infiltration levels of 28 immune cells between high- and low-score groups in the pan-cancer validation cohorts (LUSC, UCEC and LAML); Figure S6: Kaplan-Meier curves for patients in the TCGA-CESC and GSE44001 cohorts stratified by the risk score and expressions of immune checkpoint BTLA; Table S1. Univariate and multivariate Cox regression proportional hazards model to analyze the correlation between the expression levels of 206 candidate genes and overall survival of patients with cervical cancers.

Author Contributions

J.Z. and S.W. conceived and designed the whole project and drafted the manuscript; S.Y., X.L., M.M. and R.Y. analyzed the data and wrote the manuscript; S.Y. and J.Z. carried out data interpretations and helped data discussion; J.Z. and S.W. revised the manuscript. All authors have read and agreed to the published version of the manuscript.

Funding

This study was supported by Grants from the National Natural Science Foundation of China (No. 81971340 and 81502230), Natural Science Foundation of Shanghai (22ZR1450700), Shanghai Aging and Women and Children’s Health Research Project (No. 2020YJZX0215).

Institutional Review Board Statement

The study was approved by the Institutional Ethics Committee of Shanghai General Hospital (2022AWS0035). Cervical cancer (n = 32) and normal cervix (n = 32) tissue microarrays were purchased from Shanghai Zuocheng Biotech (Shanghai, China). We performed the IHC analysis as previously described [19].

Informed Consent Statement

Not applicable.

Data Availability Statement

The data presented in this study are available in this article (and Supplementary Material).

Acknowledgments

We sincerely acknowledge the contributions from the TCGA project and the GEO project.

Conflicts of Interest

The authors declare that they have no conflict of interest.

Abbreviations

AUC, area under curve; DCA, decision curve analysis; DEG, differentially expressed gene; EMT, epithelial to mesenchymal transition; ICP, immune checkpoint; OS, overall survival; DFS, disease-free survival; NPV, negative predictive value; PPV, positive predictive value; ROC, receiver operating characteristic; TCGA, The Cancer Genome Atlas; TIL, tumor infiltrating lymphocyte; TME, tumor microenvironment.

References

  1. Sung, H.; Ferlay, J.; Siegel, R.L.; Laversanne, M.; Soerjomataram, I.; Jemal, A.; Bray, F. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J. Clin. 2021, 71, 209–249. [Google Scholar] [CrossRef] [PubMed]
  2. Pavlova, N.N.; Thompson, C.B. The Emerging Hallmarks of Cancer Metabolism. Cell Metab. 2016, 23, 27–47. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  3. Sun, H.; Zhou, Y.; Skaro, M.F.; Wu, Y.; Qu, Z.; Mao, F.; Zhao, S.; Xu, Y. Metabolic Reprogramming in Cancer Is Induced to Increase Proton Production. Cancer Res. 2020, 80, 1143–1155. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  4. Cairns, R.A.; Harris, I.S.; Mak, T.W. Regulation of cancer cell metabolism. Nat. Rev. Cancer. 2011, 11, 85–95. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  5. Cao, Y. Adipocyte and lipid metabolism in cancer drug resistance. J. Clin. Investig. 2019, 129, 3006–3017. [Google Scholar] [CrossRef] [Green Version]
  6. Tasdogan, A.; Faubert, B.; Ramesh, V.; Ubellacker, J.M.; Shen, B.; Solmonson, A.; Murphy, M.M.; Gu, Z.; Gu, W.; Martin, M.; et al. Metabolic heterogeneity confers differences in melanoma metastatic potential. Nature 2020, 577, 115–120. [Google Scholar] [CrossRef]
  7. Bader, J.E.; Voss, K.; Rathmell, J.C. Targeting Metabolism to Improve the Tumor Microenvironment for Cancer Immunotherapy. Mol. Cell 2020, 78, 1019–1033. [Google Scholar] [CrossRef]
  8. Riera Leal, A.; Ortiz-Lazareno, P.C.; Jave-Suárez, L.F.; Ramírez De Arellano, A.; Aguilar-Lemarroy, A.; Ortiz-García, Y.M.; Barrón-Gallardo, C.A.; Solís-Martínez, R.; Luquin De Anda, S.; Muñoz-Valle, J.F.; et al. 17β-estradiol-induced mitochondrial dysfunction and Warburg effect in cervical cancer cells allow cell survival under metabolic stress. Int. J. Oncol. 2020, 56, 33–46. [Google Scholar] [CrossRef]
  9. Shang, C.; Wang, W.; Liao, Y.; Chen, Y.; Liu, T.; Du, Q.; Huang, J.; Liang, Y.; Liu, J.; Zhao, Y.; et al. LNMICC Promotes Nodal Metastasis of Cervical Cancer by Reprogramming Fatty Acid Metabolism. Cancer Res. 2018, 78, 877–890. [Google Scholar] [CrossRef] [Green Version]
  10. Chen, X.; Yi, C.; Yang, M.-J.; Sun, X.; Liu, X.; Ma, H.; Li, Y.; Li, H.; Wang, C.; He, Y.; et al. Metabolomics study reveals the potential evidence of meta-bolic reprogramming towards the Warburg effect in precancerous lesions. J. Cancer 2021, 12, 1563–1574. [Google Scholar] [CrossRef]
  11. Dias, A.S.; Almeida, C.R.; Helguero, L.A.; Duarte, I.F. Metabolic crosstalk in the breast cancer microenvironment. Eur. J. Cancer 2019, 121, 154–171. [Google Scholar] [CrossRef] [PubMed]
  12. Ngwa, V.M.; Edwards, D.N.; Philip, M.; Chen, J. Microenvironmental Metabolism Regulates Antitumor Immunity. Cancer Res. 2019, 79, 4003–4008. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  13. Stone, S.C.; Rossetti, R.A.M.; Alvarez, K.L.F.; Carvalho, J.P.; Margarido, P.F.R.; Baracat, E.C.; Tacla, M.; Boccardo, E.; Yokochi, K.; Lorenzi, N.P.; et al. Lactate secreted by cervical cancer cells modulates macrophage phenotype. J. Leukoc. Biol. 2019, 105, 1041–1054. [Google Scholar] [CrossRef] [PubMed]
  14. Fan, Y.; Li, X.; Tian, L.; Wang, J. Identification of a Metabolism-Related Signature for the Prediction of Survival in Endometrial Cancer Patients. Front. Oncol. 2021, 11, 630905. [Google Scholar] [CrossRef]
  15. Huang, R.; Li, G.; Wang, Z.; Hu, H.; Zeng, F.; Zhang, K.; Wang, K.; Wu, F. Identification of an ATP metabolism-related signature as-sociated with prognosis and immune microenvironment in gliomas. Cancer Sci. 2020, 111, 2325–2335. [Google Scholar] [CrossRef]
  16. Yu, S.; Hu, C.; Cai, L.; Du, X.; Lin, F.; Yu, Q.; Liu, L.; Zhang, C.; Liu, X.; Li, W.; et al. Seven-Gene Signature Based on Glycolysis Is Closely Related to the Prognosis and Tumor Immune Infiltration of Patients with Gastric Cancer. Front. Oncol. 2020, 10, 1778. [Google Scholar] [CrossRef]
  17. Charoentong, P.; Finotello, F.; Angelova, M.; Mayer, C.; Efremova, M.; Rieder, D.; Hackl, H.; Trajanoski, Z. Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep. 2017, 18, 248–262. [Google Scholar] [CrossRef] [Green Version]
  18. Miao, Y.-R.; Zhang, Q.; Lei, Q.; Luo, M.; Xie, G.-Y.; Wang, H.; Guo, A.-Y. ImmuCellAI: A Unique Method for Comprehensive T-Cell Subsets Abundance Prediction and its Application in Cancer Immunotherapy. Adv. Sci. 2020, 7, 1902880. [Google Scholar] [CrossRef] [Green Version]
  19. Zhang, J.; Chen, Y.; Chen, X.; Zhang, W.; Zhao, L.; Weng, L.; Tian, H.; Wu, Z.; Tan, X.; Ge, X.; et al. Deubiquitinase USP35 restrains STING-mediated interferon signaling in ovarian cancer. Cell Death Differ. 2021, 28, 139–155. [Google Scholar] [CrossRef]
  20. MacGregor, H.L.; Sayad, A.; Elia, A.; Wang, B.X.; Katz, S.R.; Shaw, P.A.; Clarke, B.A.; Crome, S.Q.; Robert-Tissot, C.; Bernardini, M.Q.; et al. High expression of B7-H3 on stromal cells defines tumor and stromal compartments in epithelial ovarian cancer and is associated with limited immune activation. J. Immunother. Cancer 2019, 7, 357. [Google Scholar] [CrossRef]
  21. Ru, B.; Wong, C.N.; Tong, Y.; Zhong, J.Y.; Zhong, S.S.W.; Wu, W.C.; Chu, K.C.; Wong, C.Y.; Lau, C.Y.; Chen, I.; et al. TISIDB: An integrated repository portal for tumor-immune system interactions. Bioinformatics 2019, 35, 4200–4202. [Google Scholar] [CrossRef] [PubMed]
  22. Lv, J.; Wang, J.; Shen, X.; Liu, J.; Zhao, D.; Wei, M.; Li, X.; Fan, B.; Sun, Y.; Xue, F.; et al. A serum metabolomics analysis reveals a panel of screening metabolic biomarkers for esophageal squamous cell carcinoma. Clin. Transl. Med. 2021, 11, e419. [Google Scholar] [CrossRef] [PubMed]
  23. Khalaf, K.; Hana, D.; Chou, J.T.-T.; Singh, C.; Mackiewicz, A.; Kaczmarek, M. Aspects of the Tumor Microenvironment Involved in Immune Resistance and Drug Resistance. Front. Immunol. 2021, 12, 656364. [Google Scholar] [CrossRef] [PubMed]
  24. Traba, J.; Sack, M.N.; Waldmann, T.A.; Anton, O.M. Immunometabolism at the Nexus of Cancer Therapeutic Efficacy and Resistance. Front. Immunol. 2021, 12, 657293. [Google Scholar] [CrossRef]
  25. Alsheikh, H.A.M.; Metge, B.J.; Ha, C.-M.; Hinshaw, D.C.; Mota, M.S.V.; Kammerud, S.C.; Lama-Sherpa, T.; Sharafeldin, N.; Wende, A.R.; Samant, R.S.; et al. Normalizing glucose levels reconfigures the mammary tumor immune and metabolic microenvironment and decreases metastatic seeding. Cancer Lett. 2021, 517, 24–34. [Google Scholar] [CrossRef]
  26. Sazeides, C.; Le, A. Metabolic Relationship between Cancer-Associated Fibroblasts and Cancer Cells. Adv. Exp. Med. Biol. 2021, 1311, 189–204. [Google Scholar]
  27. Reina-Campos, M.; Moscat, J.; Diaz-Meco, M. Metabolism shapes the tumor microenvironment. Curr. Opin. Cell Biol. 2017, 48, 47–53. [Google Scholar] [CrossRef]
  28. Han, A.; Schug, Z.T.; Aplin, A.E. Metabolic Alterations and Therapeutic Opportunities in Rare Forms of Melanoma. Trends Cancer 2021, 7, 671–681. [Google Scholar] [CrossRef]
  29. Long, J.; Zhang, C.-J.; Zhu, N.; Du, K.; Yin, Y.-F.; Tan, X.; Liao, D.-F.; Qin, L. Lipid metabolism and carcinogenesis, cancer development. Am. J. Cancer Res. 2018, 8, 778–791. [Google Scholar]
  30. Ma, L.; Hernandez, M.O.; Zhao, Y.; Mehta, M.; Tran, B.; Kelly, M.; Rae, Z.; Hernandez, J.M.; Davis, J.L.; Martin, S.P.; et al. Tumor Cell Biodiversity Drives Microenvironmental Reprogramming in Liver Cancer. Cancer Cell 2019, 36, 418–430.e6. [Google Scholar] [CrossRef] [Green Version]
  31. Payne, K.K. Cellular stress responses and metabolic reprogramming in cancer progression and dormancy. Semin. Cancer Biol. 2022, 78, 45–48. [Google Scholar] [CrossRef] [PubMed]
  32. Wu, J.; Li, L.; Zhang, H.; Zhao, Y.; Zhang, H.; Wu, S.; Xu, B. A risk model developed based on tumor microenvironment predicts overall survival and associates with tumor immunity of patients with lung adenocarcinoma. Oncogene 2021, 40, 4413–4424. [Google Scholar] [CrossRef] [PubMed]
  33. Sun, Y.; Wu, J.; Yuan, Y.; Lu, Y.; Luo, M.; Lin, L.; Ma, S. Construction of a Promising Tumor-Infiltrating CD8+ T Cells Gene Signature to Improve Prediction of the Prognosis and Immune Response of Uveal Melanoma. Front. Cell Dev. Biol. 2021, 9, 673838. [Google Scholar] [CrossRef] [PubMed]
  34. Chen, Q.-Y.; Chen, Y.-X.; Han, Q.-Y.; Zhang, J.-G.; Zhou, W.-J.; Zhang, X.; Ye, Y.-H.; Yan, W.-H.; Lin, A. Prognostic Significance of Immune Checkpoints HLA-G/ILT-2/4 and PD-L1 in Colorectal Cancer. Front. Immunol. 2021, 12, 679090. [Google Scholar] [CrossRef] [PubMed]
  35. Wang, Y.-Q.; Zhang, Y.; Jiang, W.; Chen, Y.-P.; Xu, S.-Y.; Liu, N.; Zhao, Y.; Li, L.; Lei, Y.; Hong, X.-H.; et al. Development and validation of an immune check-point-based signature to predict prognosis in nasopharyngeal carcinoma using computational pathology analysis. J. Immunother. Cancer 2019, 7, 298. [Google Scholar] [CrossRef]
  36. Yao, Y.; Yan, Z.; Lian, S.; Wei, L.; Zhou, C.; Feng, D.; Zhang, Y.; Yang, J.; Li, M.; Chen, Y. Prognostic value of novel immune-related genomic biomarkers identified in head and neck squamous cell carcinoma. J. Immunother. Cancer 2020, 8, e000444. [Google Scholar] [CrossRef]
  37. Jiang, M.; Wu, C.; Zhang, L.; Sun, C.; Wang, H.; Xu, Y.; Sun, H.; Zhu, J.; Zhao, W.; Fang, Q.; et al. FOXP3-based immune risk model for recurrence prediction in small-cell lung cancer at stages I-III. J. Immunother. Cancer 2021, 9, e002339. [Google Scholar] [CrossRef]
  38. Marliot, F.; Chen, X.; Kirilovsky, A.; Sbarrato, T.; El Sissy, C.; Batista, L.; Van den Eynde, M.; Haicheur-Adjouri, N.; Anitei, M.-G.; Musina, A.-M.; et al. Analytical validation of the Immunoscore and its associated prognostic value in patients with colon cancer. J. Immunother. Cancer 2020, 8, e000272. [Google Scholar] [CrossRef]
  39. Wang, Y.; Liu, Y.; Guan, Y.; Li, H.; Liu, Y.; Zhang, M.; Cui, P.; Kong, D.; Chen, X.; Yin, H. Integrated analysis of immune-related genes in endometrial carcinoma. Cancer Cell Int. 2020, 20, 477. [Google Scholar] [CrossRef]
  40. Cheong, J.E.; Sun, L. Targeting the IDO1/TDO2-KYN-AhR Pathway for Cancer Immunotherapy—Challenges and Opportunities. Trends Pharmacol. Sci. 2018, 39, 307–325. [Google Scholar] [CrossRef]
  41. Chen, F.; Xu, G.; Tian, W.; Gou, S. Breakdown of chemo-immune resistance by a TDO2-targeted Pt(IV) prodrug via attenuating endogenous Kyn-AhR-AQP4 metabolic circuity and TLS-promoted genomic instability. Biochem. Pharmacol. 2021, 193, 114785. [Google Scholar] [CrossRef] [PubMed]
  42. Noordmans, G.A.; Caputo, C.R.; Huang, Y.; Sheehan, S.M.; Bulthuis, M.; Heeringa, P.; Hillebrands, J.-L.; van Goor, H.; Korstanje, R. Genetic analysis of mesangial matrix expansion in aging mice and identification of Far2 as a candidate gene. J. Am. Soc. Nephrol. 2013, 24, 1995–2001. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  43. Backer, G.; Eddy, S.; Sheehan, S.M.; Takemon, Y.; Reznichenko, A.; Savage, H.S.; Kretzler, M.; Korstanje, R. FAR2 is associated with kidney disease in mice and humans. Physiol. Genom. 2018, 50, 543–552. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  44. Xu, J.; Zhu, S.; Xu, L.; Liu, X.; Ding, W.; Wang, Q.; Chen, Y.; Deng, H. CA9 Silencing Promotes Mitochondrial Biogenesis, Increases Putrescine Toxicity and Decreases Cell Motility to Suppress ccRCC Progression. Int. J. Mol. Sci. 2020, 21, 5939. [Google Scholar] [CrossRef]
  45. McDonald, P.C.; Chafe, S.C.; Brown, W.S.; Saberi, S.; Swayampakula, M.; Venkateswaran, G.; Nemirovsky, O.; Gillespie, J.A.; Karasinska, J.M.; Kalloger, S.E.; et al. Regulation of pH by Carbonic Anhydrase 9 Mediates Survival of Pancreatic Cancer Cells with Activated KRAS in Response to Hypoxia. Gastroenterology 2019, 157, 823–837. [Google Scholar] [CrossRef] [Green Version]
  46. Li, Z.; Jiang, L.; Chew, S.H.; Hirayama, T.; Sekido, Y.; Toyokuni, S. Carbonic anhydrase 9 confers resistance to ferroptosis/apoptosis in malignant mesothelioma under hypoxia. Redox Biol. 2019, 26, 101297. [Google Scholar] [CrossRef]
  47. Da Motta, L.L.; Ledaki, I.; Purshouse, K.; Haider, S.; De Bastiani, M.A.; Baban, D.; Morotti, M.; Steers, G.; Wigfield, S.; Bridges, E.; et al. The BET inhibitor JQ1 selectively impairs tumour response to hypoxia and downregulates CA9 and angiogenesis in triple negative breast cancer. Oncogene 2017, 36, 122–132. [Google Scholar] [CrossRef] [Green Version]
  48. Guan, C.; Ouyang, D.; Qiao, Y.; Li, K.; Zheng, G.; Lao, X.; Zhang, S.; Liao, G.; Liang, Y. CA9 transcriptional expression determines prognosis and tumour grade in tongue squamous cell carcinoma patients. J. Cell. Mol. Med. 2020, 24, 5832–5841. [Google Scholar] [CrossRef] [Green Version]
  49. Song, X.; Zhu, S.; Xie, Y.; Liu, J.; Sun, L.; Zeng, D.; Wang, P.; Ma, X.; Kroemer, G.; Bartlett, D.L.; et al. JTC801 Induces pH-dependent Death Specifically in Cancer Cells and Slows Growth of Tumors in Mice. Gastroenterology 2018, 154, 1480–1493. [Google Scholar] [CrossRef]
  50. Tinholt, M.; Garred, Ø.; Borgen, E.; Beraki, E.; Schlichting, E.; Kristensen, V.; Sahlberg, K.K.; Iversen, N. Subtype-specific clinical and prognostic relevance of tumor-expressed F5 and regulatory F5 variants in breast cancer: The CoCaV study. J. Thromb. Haemost. 2018, 16, 1347–1356. [Google Scholar] [CrossRef] [Green Version]
  51. Tinholt, M.; Stavik, B.; Tekpli, X.; Garred, Ø.; Borgen, E.; Kristensen, V.; Sahlberg, K.K.; Sandset, P.M.; Iversen, N. Coagulation factor V is a marker of tumor-infiltrating immune cells in breast cancer. Oncoimmunology 2020, 9, 1824644. [Google Scholar] [CrossRef] [PubMed]
  52. Song, J.; Liu, Y.; Guan, X.; Zhang, X.; Yu, W.; Li, Q. A Novel Ferroptosis-Related Biomarker Signature to Predict Overall Survival of Esophageal Squamous Cell Carcinoma. Front. Mol. Biosci. 2021, 8, 675193. [Google Scholar] [CrossRef] [PubMed]
  53. Egolf, S.; Zou, J.; Anderson, A.; Simpson, C.L.; Aubert, Y.; Prouty, S.; Ge, K.; Seykora, J.T.; Capell, B.C. MLL4 mediates differentiation and tumor suppression through ferroptosis. Sci. Adv. 2021, 7, eabj9141. [Google Scholar] [CrossRef] [PubMed]
Figure 1. Construction of a metabolism-related risk score signature in cervical cancer. (A) Schematic diagram of this study design. (B,C) Identification of differentially expressed metabolism-related genes (DEGs) between tumor and normal samples using GEO dataset (GSE63514) and annotation of GPL570 platform. The volcano plot (B) and Venn diagram (C) are shown. (D) The uni- and multi-variate Cox regression analysis results of the five metabolism-related hub genes in TCGA-CESC cohort. (E,F) The association between risk score and clinicopathological characters as well as enrichment scores of specific gene sets, including EMT, angiogenesis and hypoxia, of patients in TCGA cohort. The statistical difference of two groups was compared through the Wilcoxon test. * p < 0.05; **** p < 0.0001. (G) Representative immunostaining pictures of the five hub genes (ALOX12B, FAR2, F5, CA9 and TDO2) in tumor and normal tissues. Scale bar = 50 μm. The protein levels were plotted as a boxplot. * p < 0.05. MRGs, metabolism-related genes; ANT, adjacent non-tumor tissue; SCC, squamous cell carcinoma.
Figure 1. Construction of a metabolism-related risk score signature in cervical cancer. (A) Schematic diagram of this study design. (B,C) Identification of differentially expressed metabolism-related genes (DEGs) between tumor and normal samples using GEO dataset (GSE63514) and annotation of GPL570 platform. The volcano plot (B) and Venn diagram (C) are shown. (D) The uni- and multi-variate Cox regression analysis results of the five metabolism-related hub genes in TCGA-CESC cohort. (E,F) The association between risk score and clinicopathological characters as well as enrichment scores of specific gene sets, including EMT, angiogenesis and hypoxia, of patients in TCGA cohort. The statistical difference of two groups was compared through the Wilcoxon test. * p < 0.05; **** p < 0.0001. (G) Representative immunostaining pictures of the five hub genes (ALOX12B, FAR2, F5, CA9 and TDO2) in tumor and normal tissues. Scale bar = 50 μm. The protein levels were plotted as a boxplot. * p < 0.05. MRGs, metabolism-related genes; ANT, adjacent non-tumor tissue; SCC, squamous cell carcinoma.
Cancers 14 02399 g001
Figure 2. Verification of the metabolism-related risk score signature. (A) The left panel shows the risk curve and scatter plot of each sample in TCGA-CESC project reordered by the risk score, and the heatmap of expression profiles of the five hub genes. The middle panel displays the results of Kaplan-Meier analysis. The 3- and 5-year ROC curves curves of this optimal model (the right panel) revealed the AUC values. (B) The risk curve, heatmap of expression profiles (the left panel), results of Kaplan-Meier analysis (the middle panel) and ROC curves at 3 and 5 years (the right panel) of patients in GSE44001 cohort. (C) Forest plot of the univariate Cox regression analyses results of this risk score signature in all 33 types of cancer from TCGA database (the left panel). The Kaplan-Meier survival analyses and ROC curves of TCGA-HNSC and TCGA-LGG cohorts were plotted on the right panel. AUC, area under curve; ROC, receiver operating characteristic. * p < 0.05; ** p < 0.01; *** p < 0.001, **** p < 0.0001.
Figure 2. Verification of the metabolism-related risk score signature. (A) The left panel shows the risk curve and scatter plot of each sample in TCGA-CESC project reordered by the risk score, and the heatmap of expression profiles of the five hub genes. The middle panel displays the results of Kaplan-Meier analysis. The 3- and 5-year ROC curves curves of this optimal model (the right panel) revealed the AUC values. (B) The risk curve, heatmap of expression profiles (the left panel), results of Kaplan-Meier analysis (the middle panel) and ROC curves at 3 and 5 years (the right panel) of patients in GSE44001 cohort. (C) Forest plot of the univariate Cox regression analyses results of this risk score signature in all 33 types of cancer from TCGA database (the left panel). The Kaplan-Meier survival analyses and ROC curves of TCGA-HNSC and TCGA-LGG cohorts were plotted on the right panel. AUC, area under curve; ROC, receiver operating characteristic. * p < 0.05; ** p < 0.01; *** p < 0.001, **** p < 0.0001.
Cancers 14 02399 g002
Figure 3. Functional enrichment analysis of the metabolism-related risk score signature. (A) Immune network of the 24 ImmuCellAI cell types in the TCGA cohort. The size of each cell was calculated by the formula log10 (p-values of univariate Cox regression analyses). The color of each cell was used to represent the different survival impact of these cell types. The thickness of the lines estimated by Spearman correlation analyses depicted the strength of correlation between diverse cell types. Red represents positive correlation whereas negative is in blue. (B) Correlation between these cell types and our established risk score in the TCGA cohort. Spearman analyses were applied to calculate the correlation coefficients and p-value < 0.05 was enrolled. (C) Gene ontology (GO) enrichment analysis of the differentially expressed genes between high- and low-risk groups in the TCGA-CESC cohort. Adjusted p-value < 0.05 was considered statistically significant. (D) GSVA analysis of hallmark pathways in the TCGA cohort was performed. Differential analysis of GSVA score between high- and low-risk groups was displayed. (E) Patients in the high-risk group were associated with lower infiltrating density of most cell types according to Charoentong’s research. * p < 0.05; ** p < 0.01; *** p < 0.001.
Figure 3. Functional enrichment analysis of the metabolism-related risk score signature. (A) Immune network of the 24 ImmuCellAI cell types in the TCGA cohort. The size of each cell was calculated by the formula log10 (p-values of univariate Cox regression analyses). The color of each cell was used to represent the different survival impact of these cell types. The thickness of the lines estimated by Spearman correlation analyses depicted the strength of correlation between diverse cell types. Red represents positive correlation whereas negative is in blue. (B) Correlation between these cell types and our established risk score in the TCGA cohort. Spearman analyses were applied to calculate the correlation coefficients and p-value < 0.05 was enrolled. (C) Gene ontology (GO) enrichment analysis of the differentially expressed genes between high- and low-risk groups in the TCGA-CESC cohort. Adjusted p-value < 0.05 was considered statistically significant. (D) GSVA analysis of hallmark pathways in the TCGA cohort was performed. Differential analysis of GSVA score between high- and low-risk groups was displayed. (E) Patients in the high-risk group were associated with lower infiltrating density of most cell types according to Charoentong’s research. * p < 0.05; ** p < 0.01; *** p < 0.001.
Cancers 14 02399 g003
Figure 4. Correlation between the metabolism-related risk score and immune landscape. (A) Representative immunostaining pictures of the five hub genes and four cell types (CD4+, CD8+, CD57+ and CD68+ cells). The upper panel comprises images of five hub genes, images of four cell types were in the middle panel. Scale bar: 50 μm. The lower panel illustrates the infiltration scores of tumor-infiltrating CD4+, CD8+, CD57+ and CD68+ cells in the epithelial or stromal cell compartments. * p < 0.05. (B) Differences in the infiltration levels of 28 immune cells between high- and low-score groups in the pan-cancer validation cohorts. * p < 0.05; ** p < 0.01; *** p < 0.001.
Figure 4. Correlation between the metabolism-related risk score and immune landscape. (A) Representative immunostaining pictures of the five hub genes and four cell types (CD4+, CD8+, CD57+ and CD68+ cells). The upper panel comprises images of five hub genes, images of four cell types were in the middle panel. Scale bar: 50 μm. The lower panel illustrates the infiltration scores of tumor-infiltrating CD4+, CD8+, CD57+ and CD68+ cells in the epithelial or stromal cell compartments. * p < 0.05. (B) Differences in the infiltration levels of 28 immune cells between high- and low-score groups in the pan-cancer validation cohorts. * p < 0.05; ** p < 0.01; *** p < 0.001.
Cancers 14 02399 g004
Figure 5. Different expression patterns of immunostimulators (A), immunoinhibitors (B), receptors (C) and chemokines (D) between high- and low-score groups in TCGA training cohort and GSE44001 validation cohort. * p < 0.05; ** p < 0.01; *** p < 0.001.
Figure 5. Different expression patterns of immunostimulators (A), immunoinhibitors (B), receptors (C) and chemokines (D) between high- and low-score groups in TCGA training cohort and GSE44001 validation cohort. * p < 0.05; ** p < 0.01; *** p < 0.001.
Cancers 14 02399 g005
Figure 6. The potential of the established signature in predicting the immunotherapeutic benefits of cervical cancer patients. (A) Kaplan-Meier curves for patients in the TCGA-CESC cohort stratified by the risk score and expressions of immune checkpoints, such as PD1, PD-L1, CTLA-4, CD28, CD38 and CD47. * p < 0.05; ** p < 0.01; *** p < 0.001; **** p < 0.0001. (B) Kaplan-Meier curves for patients in the GSE44001 validation cohort stratified by the risk score and expressions of immune checkpoints, such as PD1, PD-L1, CTLA-4, CD28, CD38 and CD47. * p < 0.05; ** p < 0.01.
Figure 6. The potential of the established signature in predicting the immunotherapeutic benefits of cervical cancer patients. (A) Kaplan-Meier curves for patients in the TCGA-CESC cohort stratified by the risk score and expressions of immune checkpoints, such as PD1, PD-L1, CTLA-4, CD28, CD38 and CD47. * p < 0.05; ** p < 0.01; *** p < 0.001; **** p < 0.0001. (B) Kaplan-Meier curves for patients in the GSE44001 validation cohort stratified by the risk score and expressions of immune checkpoints, such as PD1, PD-L1, CTLA-4, CD28, CD38 and CD47. * p < 0.05; ** p < 0.01.
Cancers 14 02399 g006
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Yu, S.; Li, X.; Ma, M.; Yang, R.; Zhang, J.; Wu, S. The Immunological Contribution of a Novel Metabolism-Related Signature to the Prognosis and Anti-Tumor Immunity in Cervical Cancer. Cancers 2022, 14, 2399. https://doi.org/10.3390/cancers14102399

AMA Style

Yu S, Li X, Ma M, Yang R, Zhang J, Wu S. The Immunological Contribution of a Novel Metabolism-Related Signature to the Prognosis and Anti-Tumor Immunity in Cervical Cancer. Cancers. 2022; 14(10):2399. https://doi.org/10.3390/cancers14102399

Chicago/Turabian Style

Yu, Sihui, Xi Li, Mingjun Ma, Rui Yang, Jiawen Zhang, and Sufang Wu. 2022. "The Immunological Contribution of a Novel Metabolism-Related Signature to the Prognosis and Anti-Tumor Immunity in Cervical Cancer" Cancers 14, no. 10: 2399. https://doi.org/10.3390/cancers14102399

APA Style

Yu, S., Li, X., Ma, M., Yang, R., Zhang, J., & Wu, S. (2022). The Immunological Contribution of a Novel Metabolism-Related Signature to the Prognosis and Anti-Tumor Immunity in Cervical Cancer. Cancers, 14(10), 2399. https://doi.org/10.3390/cancers14102399

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop