Exploration of Immunogenic Cell Death-Associated Genes for Predicting Prognosis and Immunological Characteristics in Cervical Squamous Cell Carcinoma

Wang, Kunyu; Chen, Yiran; Zhang, Yanan; Liu, Shuanghuan; Ao, Miao; Yang, Huansong; Li, Bin

doi:https://doi.org/10.1155/2023/1405635

Journal of Oncology

On this page

Abstract Introduction Materials and Methods Results Discussion Conclusion Data Availability Disclosure Conflicts of Interest Authors’ Contributions Acknowledgments Supplementary Materials References Copyright Related Articles

Special Issue

Novel Biomarkers for Predicting Response to Cancer Immunotherapy

View this Special Issue

Research Article | Open Access

Volume 2023 | Article ID 1405635 | https://doi.org/10.1155/2023/1405635

Exploration of Immunogenic Cell Death-Associated Genes for Predicting Prognosis and Immunological Characteristics in Cervical Squamous Cell Carcinoma

Kunyu Wang,¹Yiran Chen,¹Yanan Zhang,¹Shuanghuan Liu,¹Miao Ao,¹Huansong Yang,¹and Bin Li¹

Academic Editor: Jinghua Pan

Received13 Aug 2022

Revised28 Aug 2022

Accepted30 Aug 2022

Published30 Jan 2023

Abstract

Background. The tumor microenvironment (TME) has gradually entered the vision of researchers and is becoming a vital part of the occurrence of cervical squamous cell carcinoma (CSCC). However, understanding the specific composition of TME still confront enormous challenges, particularly immune and stromal components. Methods. In this study, we performed an unsupervised cluster analysis to determine the immunogenic cell death-associated subtype of CSCC patients. The differences in immune status, genomic alteration, and clinical outcomes between each subtype were compared. Subsequently, we screened vital prognostic factors. The HPA database was employed to verify the protein localization and the expression level between cancer and adjacent tissues. Results. CSCC patients were divided into three subtypes according to the expression of immunogenic cell death-associated genes. Cluster C has the highest survival rate because of the lower activation of tumor-related pathways. The immune score and stromal score of patients with Cluster B were the highest, so it may be considered that stromal tissue inhibits the anti-tumor effect of immunocytes. In addition, we constructed a risk score based on immunogenic cell death-associated genes to screen for vital markers. We systematically revealed the genomic alteration of vital markers. Conclusions. We have established a novel immunogenic cell death-associated risk scoring system in CSCC, and the expression of immunogenic cell death-associated genes may be a valuable biomarker for immunotherapy strategies. Our work may contribute to the development of new immunomodulators and develop new precision immunotherapies for CSCC.

1. Introduction

Cervical squamous cell carcinoma (CSCC) is a common gynecologic tumor worldwide that is regulated by many mechanisms and progresses slowly. According to the global cancer statistics in 2021, the incidence and mortality of CESC ranked second, and it is still one of the main causes of female death [1]. Nevertheless, the specific tumor pathogenesis has not yet been fully explained; existing studies generally believe that cervical squamous cell carcinoma is mainly caused by HPV infection and DNA methylation of involved genes [2, 3]. In most countries, common preventive measures for CESC patients include HPV vaccination and screening. Nevertheless, these measures help to reduce the incidence of cancer through early cancer prevention but have no therapeutic effect on CESC. The incidence of postoperative recurrence or metastasis in CESC patients is high, and the 5-year survival rate is as low as 10–20% [4, 5]. In short, there is no better treatment for cervical squamous cell carcinoma at present, so it is urgent to find new prognosis and treatment indicators.

Tumor tissue includes not only malignant cells but also tumor-related stromal cells and immunocytes, which together with cancer cells constitute the tumor microenvironment (TME) [6, 7]. TME has recently received extensive attention, which provides new clues for tumor diagnosis, treatment, and prognosis [8, 9]. Stromal cells in the tumor microenvironment have the potential to maintain cell stemness and promote the spread of tumor cells at the site of metastasis [10–12]. The tumor microenvironment can be estimated and quantified by using the “ESTIMATE” algorithm for calculating stromal and immune cells of malignant tumors based on gene expression values [13]. Many studies have shown that stromal score, immune score, and tumor purity measurement, which are calculated by the “ESTIMATE” algorithm, can be used as prognostic tumor biomarkers. For CESC patients, the tumor microenvironment immunomodulators can be used to predict responses to immunotherapy against PD-1/PD-L1 molecules [14, 15]. In addition, the risk score based on TME was proved to be an effective prognostic predictor of CESC patients [16–18]. However, the composition of TME is complex, with multiple innate and adaptive immunocyte populations that affect tumor immune escape and response to immunotherapy [19]. In addition, TME is characterized by multiple immunosuppressive factors. Therefore, it is necessary to explore and clarify the internal and external factors involved in the malignant transformation of tumors.

Immunogenic cell death is regulated by a series of immune-stimulatingdamage-associated molecules, referred to as immunogenic cell death-associated genes (ICDGs), such as protein disulfide isomerase family A member 3 (PDIA3), tumor necrosis factor (TNF), and high mobility group Box 1 (HMGB1) [20–22]. Previous studies have shown that HMGB1 derived from CESC cells enhanced irradiation resistance by activated CD8+ T cells, which reveals the importance of immunogenic cell death in clinical therapy for CESC [23]. The sensitivity of ICDGs to ICI strongly depends on the immune response time, in which immunocytes play an important role in the clinical effect prediction of CESC [24]. Oltean et al. obtained 38 patients with locally advanced cervical cancer who received cisplatin neoadjuvant chemoradiotherapy and detected the expression of cell death markers. Although these markers cannot predict the outcome of patients in terms of recurrence or survival, many markers were significantly associated with immune cell infiltration. Measuring ICD markers can reflect the effect of treatment on the tumor microenvironment through immunocytes’ recruitment and infiltration [25]. In this regard, evaluating different tumor immune microenvironments is beneficial to understand the immune characteristics of different subtypes, optimize immunotherapy, and improve the prognosis of CESC patients [26]. Although a large number of immune-related biomarkers and prognostic models have been developed [27, 28], few articles use immunogenic cell death as an actionable target.

In this study, we comprehensively analyzed the characteristics of 31 ICDGs genes to explore the predictive performance of immunogenic cell death on tumor immune microenvironment and clinical outcomes of CESC patients. In addition, we constructed a risk-scoring model based on 31 ICDGs to evaluate the prognostic value of immunogenic cell death in CESC. In a word, HMGB1 as the subject of our investigation may provide new clues for exploring the potential molecular mechanisms of different responses to immunotherapy in CESC patients, thus providing new insights into immunotherapy strategies for CESC.

2. Materials and Methods

2.1. Acquisition and Preprocessing Datasets

In this study, the samples with no survival status or PFI < 30 days were eliminated to ensure the accuracy of the datasets. First, the RNA-Seq data of CSCC patients were downloaded from the TCGA database (formatted as TPM and log (x + 1) conversion), and finally, a total of 243 tumor samples were included. The GSE44001 dataset was downloaded from the GEO database and grouped with the same admittance standard as the verification queue. A total of 299 samples were included. Using the combat function in the “sva” package to remove batch effects on GSE44001 and TCGA-CSCC, the modified queues were named as a meta-cohort. The GSE168652 dataset contains CSCC single-cell RNA transcriptome data and clinical data; its annotation information was obtained from the TISCH database. Cellular localization immunofluorescence staining and IHC images were from the HPA database. Mutation data (copy number variations, single-nucleotide variants) are processed based on the GSCA database. In addition, 31 immunogenic cell death-associated genes (ICDGs) were identified from the relevant literature [29].

2.2. Unsupervised Clustering Algorithm for Identifying Molecular Subtypes

In the meta-cohort, unsupervised consistent cluster analysis was performed to classify samples into different subtypes based on 31 ICDGs. The R package “consensus PlusterPlus” was employed to determine the number of molecular subtypes; 100 replicates were performed, and pltem = 0.8. Principal component analysis (PCA) was performed to verify whether each molecular subtype was relatively independent from each other. To verify the stability of the subtypes, GSVA was used to evaluate the differences in biological pathways between subtypes. C2.cp.kegg.v7.0. symbols as a reference gene set and FDR < 0.05 as a screening threshold. In this paper, the log-rank test was employed for analysis, and the Kaplan–Meier curve displayed the clinical prognosis of patients in the dataset.

2.3. Identification of the Immune Microenvironment among Molecular Subtypes

In the analysis of tumor microenvironment immunocytes’ infiltration, we used several algorithms, such as TIMER, CIBERSORT, QUANTISEQ, MCP-counter, XCELL, and EPIC, to evaluate the abundance of immunocytes in different samples. In addition, the ESTIMATE algorithm was used to calculate the immune score and stromal score to reflect the tumor microenvironment status. Thorsson et al. defined six immune characteristic subtypes based on transcriptome profiles of 33 solid tumors, including: wound healing (Immune C1), IFN-gamma dominant (Immune C2), inflammatory (Immune C3), lymphocyte depleted (Immune C4), immunologically quiet (Immune C5), and TGF-beta dominant (Immune C6). We compared the known immune subtypes with our predicted risk score, trying to reveal the subtype types represented by the score.

2.4. Construction and Validation of Risk Scores

To quantify the immune microenvironment status of each sample, we constructed the ICDGs-base risk score. In this paper, according to the expression of ICDGs, we employed the TCGA database for the training cohort, and the GEO database as the validation cohort for external verification. The Least Absolute Shrinkage and Selection Operator (LASSO) model was employed to remove redundant genes of ICDGs in the training cohort, and then the integration coefficient and gene expression value were calculated by multivariate Cox regression, and the risk scoring formula was established. According to the median risk score of each cohort, patients were divided into high-risk and low-risk subtypes. The time-dependent receiver operating characteristic curve was plotted to compare the prediction accuracy of the risk score and traditional clinicopathological parameters. The “survivalROC” package was used to draw the ROC curve and calculate the area under the curve (AUC). A prognostic-related nomogram was constructed using the “replot” package and validated using calibration curves.

2.5. Statistics

Kruskal–Wallis analysis was performed for differences between groups; χ² test was used for associations between covariables; the Kaplan–Meier method was used to compare survival differences between different groups; a value <0.05 was considered statistically significant.

3. Results

3.1. Construction of Immunogenic Cell Death-Associated Molecular Subtype

PCA analysis found that the sample information of the GSE44001 dataset and the TCGA-CSCC dataset was not consistent, so we used the combat function in the “sva” package to effectively reduce the batch effects on GSE44001 and TCGA-CSCC (Figures 1(a) and 1(b)). Then, we merged the TCGA-CSCC and GSE44001 datasets into a meta-cohort. In this study, a total of 31 immunogenic cell death-associated genes were collected for further analysis. The network diagram displayed the correlation between these genes and their prognostic and predictive performance (Figure 1(c)). Since the clinical outcomes of CSCC patients with different levels of immunogenic cell death were diverse, we then conducted an unsupervised clustering analysis to identify different immunogenic cell death patterns and stratified the patients into three clusters (Figures 1(d) and 1(e)). The results of the PCA analysis suggested that the three molecular subgroups were independent of each other (Figure 1(f)). The abovementioned results suggest that the patients can be effectively divided into different molecular subgroups according to the expression level of immunogenic cell death-related genes.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 1

Landscape of immunogenic cell death-associated genes in cervical squamous cell carcinoma. (a–b) Principal component analysis verified the sample distribution of two databases. (c) In the correlation among immunogenic cell death-associated genes, the size of the circle represented the effect of each gene on clinical outcomes. (d–e) Unsupervised clustering of CSCC samples based on their immunogenic cell death patterns to classify patients into three subgroups. (f) Principal component analysis verified the sample distribution of molecular subtypes.

3.2. The Transcriptome Data of Immunogenic Cell Death-Associated Molecular Subtypes

In order to further verify the accuracy of subtype classification based on ICDG expression, the heatmap was employed to display the expression levels of 31 ICDGs in different subtypes and clinical subtypes. We found that the expression of genes was different among different subgroups. For example, the TNF gene was higher in Cluster A than in the other two subgroups, and the FOXP3 gene was the highest in Cluster B. In addition, 31 ICDGs were found to be generally higher in clusters A and B than in cluster C (Figure 2(a)). Among these molecular subtypes, cluster C had a significant survival advantage, while cluster A represented the worst clinical results in the meta-cohort (Figure 2(b)).

(a)

(b)

(c)

(d)

Figure 2

Immunological characteristics of different molecular subtypes. (a) The heatmap was performed to display the expression level of 31 immunogenic cell death-associated genes in distinct subtypes. The cluster, menopause status, lymph vascular invasion, FIGO status, tumor grade, age, and project were utilized as the patients’ annotations. Red/blue represented high/low expression of genes. (b) Kaplan–Meier survival analysis for CSCC patients with distinct subtypes, . (c) The content of each tumor microenvironment immunocytes in three molecular subtypes. (d) The expression of each immunomodulator in three molecular subtypes. The median value: black lines in boxes, and the outliers: black dots out boxes.

3.3. ICDGs-Based Molecular Subtypes Reflect Immunological Characteristics

In order to explore the alteration of immunocytes’ infiltration of three molecular subtypes, ssGSEA was performed to draw a box diagram of the relative content of immunocytes in different subtypes (Figure 2(c)). The results displayed that almost all immunocytes increased in cluster B, and the infiltration of immunocytes in cluster A was lower. However, the K-M analysis showed that the survival rate of cluster A was higher than that of cluster B. We speculated that immunocytes’ infiltration, such as M2 macrophage infiltration, promoted the tumor’s malignant degree. In addition, the tumor microenvironment consists of immune cells and stromal tissues. The level of immunocyte infiltration is not the only prognostic factor. Therefore, we further detected the expression levels of immunomodulators in different molecular subtypes. Immunomodulators, such as CD274, PDCD1, IDO1, etc., have significant differences among the three molecular subgroups (Figure 2(d)). Second, the HLA recognition function refers to the unique synergy in the immune response that is used to transmit information. The HLA family also alters (Figure 3(a)). Therefore, the abovementioned results demonstrated that ICDGs participate in the regulation of the tumor immune response.

(a)

(b)

(c)

(d)

(e)

The tumor environment was also assessed using the “ESTIMATE” algorithm; the immune score was the lowest in cluster C. Then, the stromal purity (StromalScore) and tumor cell purity (ESTIMATEScore) in three clusters were also evaluated. The results displayed that the StromalScore and ESTIMATEScore of cluster B increased, while the StromalScore and ESTIMATEScore of cluster C decreased (Figure 3(b)). Therefore, as mentioned above, it is possible that the stroma tissue inhibits the immunocyte effect, and thus the prognosis of cluster B is poor.

A GSVA enrichment analysis was performed to identify the biological process of immunogenic cell death pathway. The results displayed that cluster A was enriched in tumor activation pathways, such as the MAPK signaling pathway and the JAK-STAT signaling pathway. On the other hand, cluster B is significantly associated with immune exclusion, such as autoimmune thyroid disease, antigen processing and presentation, and the toll-like receptor signaling pathway. Cluster C showed lower carcinogenicity (Figure 3(c)–3(e)). It is concluded that cluster C has a better clinical prognosis than other groups because of its insensitive carcinogenic pathway.

3.4. Identification and Validation of an ICDGs-Based Risk Score

We constructed risk scores to quantify individual levels of immunogenic cell death. The LASSO algorithm finally screened seven candidate prognostic genes from ICDGs. Then, the multivariate Cox regression algorithm was used to construct the risk-scoring model (Figure 4(a)). The integration coefficient of each gene is shown in Figure 4(b).

(a)

(b)

(c)

(d)

(e)

According to the median risk score, patients were classified as either a high-risk group or a low-risk group. The scatter plot showed that the mortality rate is proportional to the risk score (Figures 4(c) and 4(d)). The heatmap shows the expression levels of five prognostic genes in the risk group. HMGB1, PDIA3, and TNF were highly expressed in the high-risk group, while TLR4 and FOXP3 were highly expressed in the low-risk group (Figure 4(e)). According to Kaplan–Meier analysis, the OS rate of high-risk patients was significantly reduced ( value <0.001, Figures 5(c) and 5(d)). In addition, the 1-year, 3-year, and 5-year OSs of AUC values based on the TCGA-CSCC and GSE44001 cohorts were shown in Figures 5(a) and 5(b). We further studied whether the risk model had similar or better predictive validity with other clinical parameters of CSCC. We established a nomogram to predict the OS of patients, including FIGO_stage and risk score (Figure 5(e)). The calibration chart shows that the mortality rate can be accurately estimated by the nomogram (Figures 5(f) and 5(g)).

(a)

(b)

(c)

(d)

(e)

(f)

(g)

3.5. ICDGs-Based Risk Score Reflects Immune Microenvironment Status

Next, we explored whether the constructed risk score could be used to evaluate the tumor’s microenvironment status. We found that cluster B had the lowest risk score, followed by cluster C, and cluster A was the highest-scoring group (Figure 6(a)). Six immune characteristic subtypes were defined by Thorsson et al. We compared known subtypes with our constructed risk score model. C1 (wound healing) scored the highest, while C4 (lymphocyte depletion) scored the lowest (Figure 6(b)). Second, in the high-risk group, the proportion of the C2 subtype (IFN-gamma dominant) was 87%, and the proportion of the C1 subtype (wound healing) was 12%. In low-risk patients, the proportion of the C2 subtype (IFN-gamma dominant) was 69%, and the proportion of the C1 subtype (wound healing) was 30% (Figure 6(c)). Our results showed that the risk score was negatively correlated with the immune score, stromal score, and ESTIMATE score while positively correlated with the tumor score, which also revealed the cause of poor prognosis in high-risk patients (Figures 6(d)–6(g)). In addition, the TIMER, CIBERSORT, QUANTISEQ, McP-counter, XCELL, and EPIC algorithms were performed to estimate the abundance of immunocytes in different groups. Figure 6(h) shows high levels of immune cell infiltration in the low-risk group and the opposite in the high-risk group. Figure 6(i) shows that risk scores were negatively correlated with the majority of immune cell infiltration. Here, immunocytes mostly refer to B cells and T cells.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

Figure 6

ICDGs-based risk score reflects immune microenvironment status. (a–b) The correlation between the risk score and molecular subtypes/known immune subtypes was evaluated by the Kruskal–Wallis test. (c) A bar diagram displayed percentages of known immune subtypes in risk groups. (d–g) Correlation between risk score and ImmuneScore, StromalScore, TumorScore, and ESTIMATEScore. (h) The heatmap displayed the content of each tumor microenvironment’s immunocytes in risk score groups. (g) A correlation between risk score and tumor microenvironment immunocytes.

3.6. 5-ICDGs Emphasize Genomic Heterogeneity in Gynecological Oncology

The abovementioned results indicate that the 5-ICDGs gene can be employed to predict the clinical outcomes of CSCC patients. Nevertheless, its upstream mutation information has not been clearly interpreted. Therefore, this section introduces the 5-ICDGs of genomic heterogeneity in gynecological oncology. Figure 7(a) shows the landscape of mutations in gynecological tumors. Variant classification in gynecological tumors is mainly due to missense mutations. The variant type is mainly SNP, and the SNV class is mainly C > T and C>A. The most common mutation gene in 5-ICDGs is TLR4. Figure 7(b) shows that 5-ICDGs mutation frequency is the highest in uterine corpus endometrial carcinoma (UCEC), followed by OV. The 5-ICDGs mutation frequency was low in cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC). Figure 7(c) displayed copy number variation of 5-ICDGS in gynecological tumors. The light red represents heterozygous amplification; the light green indicates heterozygous deletion; the dark red represents homozygous amplification; the dark green represents homozygous deletion; the gray represents no CNV occurrence. HMGB1 and PDIA3 mainly occur by a heterozygous deletion in CESC, while TNF occurs by a heterozygous amplification.

(a)

(b)

(c)

A heatmap displayed the functional (inhibited or activated) genes in at least 5 cancer types. Results show that 5-ICDGs can be employed to activate apoptosis, the MT pathway, the PAS-MAPK pathway, etc., and inhibit cell cycle and DNA damage (Figure S1A). Figure S1B is a list of targeted drugs for 5-ICDGS based on the GDSC database. In addition, we explored differences in the survival of gynecologic tumors between gene mutations and wild-type and found that UCEC tumors were more susceptible to mutations (Figure S1C). The abovementioned analysis results suggest that patients with CESC may be less affected by gene mutation because there are other regulatory pathways.

3.7. HMGB1 is a Hub ICDGs with Tumor Promotion

Analyzing the relationship between 5-ICDGS and the clinical stage of gynecological tumors, we found that HMGB1 has a more profound influence on the clinical stage of CESC (Figure S1D). In addition, the expression level of 5-ICDGS in tumors and adjacent tissues of cervical squamous cell carcinoma is shown in Figure 8(a), and the expression level of HMGB1 in tumor tissues was increased, suggesting that HMGB1 plays an indispensable role in tumors as a cancer-promoting factor. The protein-protein interaction (PPI) network suggested that HMGB1 had the highest weight in 5-ICDGS, so we performed an analysis on HMGB1 (Figure 8(b)). Figure 8(c) shows the localization of HMGB1 protein in A231 cells; the HMGB1 protein is targeted for nucleoplasm. HE staining showed that glandular cells accounted for 10%, squamous epithelia cells accounted for 5%, smooth muscle cells accounted for 65%, and other cells accounted for 20% (Figure 8(d)). Consistent with the predicted results, immunohistochemistry showed that the HMGB1 protein was highly expressed in tumor tissues and low expressed in Paracancer tissues (Figure 8(e)). After preprocessing scRNA-seq data based on strict quality control indicators, tSNE technology was performed to visualize high-dimensionalscRNA-seq data and successfully classify cells into seven subtypes, which were then annotated as identifiable cell types according to TISCH database information (Figure 8(f)). The following main cell types were characterized: CD8+ T cells, endometrial stromal cells, endothelial, fibroblasts, malignant, mono/macro, and SMC. The high expression of the HMGB1 gene in stromal cells is fully consistent with our previous speculation that stromal tissue limits the anti-tumor effect of immunocytes.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 8

Predictive power of 5-ICDGS in gynecological tumors. (a) The expression of each 5-ICDGS in tumor and paracancer tissues; the median value: black lines in boxes; the outliers: black dots out boxes. (b) Protein-protein interaction schematic diagram of 5-ICDGS. (c–e) HMGB1 protein localization (c), an expression level (e), and tumor tissue composition (d), data were obtained from the HPA database. (f) Cells were clustered into 7 types via the tSNE dimensionality reduction algorithm, where each color represented the annotated phenotype of each cluster. The expression level of HMGB1 in each cell subpopulation.

4. Discussion

In view of the serious menace of cervical cancer to women’s health and the lack of public understanding of the genetic basis of cervical cancer, continuous work should be carried out to find remarks that have a causal relationship with cervical cancer [30, 31]. This study is one of such efforts aimed at detecting new possible pathogenic genes for cervical cancer using integrated genomics information. With the growth of high-throughput omics datasets in cancer research in recent years [32], it is well known that the application of a single level of genome measurement alone is not sufficient to completely solve the etiology of cancer prognosis [33, 34]. Based on the TCGA omics data measured on multiple platforms, we carried out in-depth research on gene expression and clinical outcomes.

Immunocyte infiltration in TME is considered to be a vital factor in the response to immunotherapy [35]. In recent decades, combinatorial evaluation of biomarkers for predictive responses has been better investigated. A well-developed biomarker tool that stratifies patients with the same type of cancer into different molecular subtypes may allow for more reasonable and personalized treatment. In our investigation, we noted the discriminative power of risk scores in patients with similar expression levels of immune checkpoint genes. This phenomenon also reveals a complex crosstalk of immunogenic cell death between immune cell infiltration, immune checkpoint genes, and patient clinical outcomes.

On the other hand, we learned that preexisting studies have developed some immune-related prognostic models to prove the relationship between the immune landscape and cervical cancer progression. For example, Chen et al. constructed six lncRNA immunoprognostic features to predict the prognosis of cervical cancer, revealing the association between the occurrence of cervical cancer and an anti-tumor immune effect [36]. Qi et al. employed the ESTIMATE algorithm to calculate an immune risk score and construct a nomogram model combining clinicopathological features for subsequent analysis. This nomogram model can effectively predict the clinical outcomes of cervical cancer [37]. However, the correlation between their models and immune infiltration was not clearly identified. In a recently published paper, researchers reclassified cervical cancer cohorts and established genetic characteristics derived from the ESTIMATE algorithm or constructed a specific TMEscore based on immunocytes’ infiltration mode to evaluate the relationship between immunocytes’ infiltration and prognosis [38, 39]. Considering the relatively enormous difference between tumors and normal tissues, we aim to explore the heterogeneity of the tumor immune landscape. Notably, we want to find vital genes with prognostic value, which may not necessarily be included in the list of known immune-related genomes used in previous studies [29]. Therefore, we analyzed the prognostic value of 28 immunocytes and developed a new risk score signature by multivariate Cox analysis. We comprehensively expounded the characteristics of this feature from the perspectives of immunocyte infiltration, immunomodulator gene expression, tumor microenvironment score, GSVA algorism, and so on. In addition, we verified the influence of key genes on the progression of cervical cancer.

The advantage of this study was to highlight the carcinogenic effect of cell death-related genes in cervical cancer and to predict the effect of immunotherapy and the degree of immune cell infiltration, which is unique in previous studies. Our research had some limitations. For example, patients in the TCGA cohort lack detailed support information, which hinders a more systematic investigation of the clinical and pathological features of CSCC patients. In addition, on account of the relative paucity of cervical cancer public datasets, the risk pattern can only be externally verified by one separate cohort. We assume that, although there is insufficient external verification, the signature may be reliable. In such a difficult environment, we verify the modeling algorithm. In the future, we will attempt to collect a multicenter, large sample dataset to better verify our conclusions.

5. Conclusion

In summary, our study proposes an immunogenic cell death-associated signature that predicts the immune landscape, clinical outcomes, and response to immunotherapy, regardless of the divergence of assays and platforms.

Data Availability

All datasets generated for this study are included within the article material, including the TCGA database (https://portal.gdc.cancer.gov/), the GEO dataset (https://www.ncbi.nlm.nih.gov/gds/): GSE44001, GSE168652, HPA database (https://www.proteinatlas.org/), and the GSCA database (https://bioinfo.life.hust.edu.cn/GSCA/#/).

Disclosure

The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Authors’ Contributions

Kunyu Wang and Yiran Chen contributed equally to this work. Kunyu Wang and Yanan Zhang downloaded the dataset, analyzed the data, and Kunyu Wang was a major contributor in writing the manuscript. Kunyu Wang and Bin Li reviewed the manuscript. Yiran Chen reanalyzed data in the revised manuscript. All authors read and approved the final manuscript.

Acknowledgments

This study was supported by the Special Research Fund for Central Universities, Peking Union Medical College (3332021030). The authors thank all members who participated in the data collection and analysis.

Supplementary Materials

Figure S1: Genome landscape of 5-ICDGS in gynecological tumors. (A) The heatmap displayed the functional (inhibited or activated) genes in at least 5 cancer types. (B) This is a list of targeted drugs for 5-ICDGS based on the GDSC database. (C) Effects of gene mutations on the clinical prognosis of gynecological tumors, including overall survival (OS), progression-free-survival (PFS), disease-specific survival (DSS), and disease-free interval (DFI). (D) Relationship between 5-ICDGs expression level and the clinical stage of a gynecological tumor. (Supplementary Materials)

References

R. L. Siegel, K. D. Miller, H. E. Fuchs, and A. Jemal, “Cancer statistics, 2021,” CA: A Cancer Journal for Clinicians, vol. 71, no. 1, pp. 7–33, 2021.
View at: Publisher Site | Google Scholar
T. A. Berman and J. T. Schiller, “Human papillomavirus in cervical cancer and oropharyngeal cancer: one cause, two diseases,” Cancer, vol. 123, no. 12, pp. 2219–2229, 2017.
View at: Publisher Site | Google Scholar
V. Gregoire and L. Bauwens, “Reply to sanguineti regarding the manuscript untitled Prevalence and distribution of cervical lymph node metastases in HPV-positive and HPV-negative oropharyngeal squamous cell carcinoma,” Radiotherapy & Oncology, vol. 161, p. 253, 2021.
View at: Publisher Site | Google Scholar
M. Rotman, A. Sedlis, M. R. Piedmonte et al., “A phase III randomized trial of postoperative pelvic irradiation in Stage IB cervical carcinoma with poor prognostic features: follow-up of a gynecologic oncology group study,” International Journal of Radiation Oncology, Biology, Physics, vol. 65, no. 1, pp. 169–176, 2006.
View at: Publisher Site | Google Scholar
W. Small Jr., W. R. Bosch, M. M. Harkenrider et al., “NRG oncology/RTOG consensus guidelines for delineation of clinical target volume for intensity modulated pelvic radiation therapy in postoperative treatment of endometrial and cervical cancer: an update,” International Journal of Radiation Oncology, Biology, Physics, vol. 109, no. 2, pp. 413–424, 2021.
View at: Publisher Site | Google Scholar
F. Ye, Y. Hu, J. Gao et al., “Radiogenomics map reveals the landscape of m6A methylation modification pattern in bladder cancer,” Frontiers in Immunology, vol. 12, Article ID 722642, 2021.
View at: Publisher Site | Google Scholar
F. Ye, Y. Liang, J. Hu et al., “DNA methylation modification map to predict tumor molecular subtypes and efficacy of immunotherapy in bladder cancer,” Frontiers in Cell and Developmental Biology, vol. 9, Article ID 760369, 2021.
View at: Publisher Site | Google Scholar
L. Peng, J. Li, J. Wu et al., “A pan-cancer analysis of SMARCA4 alterations in human cancers,” Frontiers in Immunology, vol. 12, Article ID 762598, 2021.
View at: Publisher Site | Google Scholar
S. Zuo, M. Wei, S. Wang, J. Dong, and J. Wei, “Pan-cancer analysis of immune cell infiltration identifies a prognostic immune-cell characteristic score (ICCS) in lung adenocarcinoma,” Frontiers in Immunology, vol. 11, p. 1218, 2020.
View at: Publisher Site | Google Scholar
P. R. Prasetyanti and J. P. Medema, “Intra-tumor heterogeneity from a cancer stem cell perspective,” Molecular Cancer, vol. 16, no. 1, p. 41, 2017.
View at: Publisher Site | Google Scholar
L. Deng, T. Meng, L. Chen, W. Wei, and P. Wang, “The role of ubiquitination in tumorigenesis and targeted drug discovery,” Signal Transduction and Targeted Therapy, vol. 5, no. 1, p. 11, 2020.
View at: Publisher Site | Google Scholar
S. Feng, T. Xia, Y. Ge et al., “Computed tomography imaging-based radiogenomics analysis reveals hypoxia patterns and immunological characteristics in ovarian cancer,” Frontiers in Immunology, vol. 13, Article ID 868067, 2022.
View at: Publisher Site | Google Scholar
M. Xu, Y. Li, W. Li et al., “Immune and stroma related genes in breast cancer: a comprehensive analysis of tumor microenvironment based on the cancer genome atlas (TCGA) database,” Frontiers of Medicine, vol. 7, p. 64, 2020.
View at: Publisher Site | Google Scholar
A. Marabelle, M. Fakih, J. Lopez et al., “Association of tumour mutational burden with outcomes in patients with advanced solid tumours treated with pembrolizumab: prospective biomarker analysis of the multicohort, open-label, phase 2 KEYNOTE-158 study,” The Lancet Oncology, vol. 21, no. 10, pp. 1353–1365, 2020.
View at: Publisher Site | Google Scholar
C. C. Balança, A. Salvioni, C. M. Scarlata et al., “PD-1 blockade restores helper activity of tumor-infiltrating, exhausted PD-1hiCD39+ CD4 T cells,” JCI insight, vol. 6, no. 2, Article ID e142513, 2021.
View at: Publisher Site | Google Scholar
J. Budczies, M. Bockmayr, C. Denkert et al., “Pan-cancer analysis of copy number changes in programmed death-ligand 1 (PD-L1, CD274) - associations with gene expression, mutational load, and survival,” Genes, Chromosomes and Cancer, vol. 55, no. 8, pp. 626–639, 2016.
View at: Publisher Site | Google Scholar
Y. P. Chen, Y. Zhang, J. W. Lv et al., “Genomic analysis of tumor microenvironment immune types across 14 solid cancer types: immunotherapeutic implications,” Theranostics, vol. 7, pp. 3585–3594, 2017.
View at: Publisher Site | Google Scholar
S. Feng, Y. Ge, X. Ma et al., “The genomic landscape of invasive stratified mucin-producing carcinoma of the uterine cervix: the first description based on whole-exome sequencing,” Journal of Translational Medicine, vol. 20, no. 1, p. 187, 2022.
View at: Publisher Site | Google Scholar
D. C. Hinshaw and L. A. Shevde, “The tumor microenvironment innately modulates cancer progression,” Cancer Research, vol. 79, no. 18, pp. 4557–4566, 2019.
View at: Publisher Site | Google Scholar
A. Annibaldi and P. Meier, “Checkpoints in TNF-induced cell death: implications in inflammation and cancer,” Trends in Molecular Medicine, vol. 24, no. 1, pp. 49–65, 2018.
View at: Publisher Site | Google Scholar
Y. Li, Z. Song, Q. Han et al., “Targeted inhibition of STAT3 induces immunogenic cell death of hepatocellular carcinoma cells via glycolysis,” Molecular oncology, vol. 16, no. 15, pp. 2861–2880, 2022.
View at: Publisher Site | Google Scholar
A. Ahmed and S. W. G. Tait, “Targeting immunogenic cell death in cancer,” Molecular oncology, vol. 14, no. 12, pp. 2994–3006, 2020.
View at: Publisher Site | Google Scholar
X. Jiao, S. Zhang, J. Jiao et al., “Promoter methylation of SEPT9 as a potential biomarker for early detection of cervical cancer and its overexpression predicts radioresistance,” Clinical Epigenetics, vol. 11, no. 1, p. 120, 2019.
View at: Publisher Site | Google Scholar
X. Zhang, Y. Wang, C. Qu, and J. Chen, “Pan-cancer analysis of PARP1 alterations as biomarkers in the prediction of immunotherapeutic effects and the association of its expression levels and immunotherapy signatures,” Frontiers in Immunology, vol. 12, Article ID 721030, 2021.
View at: Publisher Site | Google Scholar
T. Oltean, L. Lippens, K. Lemeire et al., “Association of cell death markers with tumor immune cell infiltrates after chemo-radiation in cervical cancer,” Frontiers in Oncology, vol. 12, Article ID 892813, 2022.
View at: Publisher Site | Google Scholar
X. Qu, Z. Shi, J. Guo, C. Guo, J. Qiu, and K. Hua, “Identification of a novel six-gene signature with potential prognostic and therapeutic value in cervical cancer,” Cancer Medicine, vol. 10, no. 19, pp. 6881–6896, 2021.
View at: Publisher Site | Google Scholar
F. Xu, J. Shen, and S. Xu, “Integrated bioinformatical analysis identifies GIMAP4 as an immune-related prognostic biomarker associated with remodeling in cervical cancer tumor microenvironment,” Frontiers in Cell and Developmental Biology, vol. 9, Article ID 637400, 2021.
View at: Publisher Site | Google Scholar
S. Yu, X. Li, J. Zhang, and S. Wu, “Development of a novel immune infiltration-based gene signature to predict prognosis and immunotherapy response of patients with cervical cancer,” Frontiers in Immunology, vol. 12, Article ID 709493, 2021.
View at: Publisher Site | Google Scholar
M. Xu, J. H. Lu, Y. Z. Zhong, and J. Jiang, “Immunogenic cell death-relevantdamage-associated molecular patterns and sensing receptors in triple-negative breast cancer molecular subtypes and implications for immunotherapy,” Frontiers in Oncology, vol. 12, Article ID 870914, 2022.
View at: Publisher Site | Google Scholar
P. Olusola, H. N. Banerjee, J. V. Philley, and S. Dasgupta, “Human papilloma virus-associated cervical cancer and health disparities,” Cells, vol. 8, no. 6, p. 622, 2019.
View at: Publisher Site | Google Scholar
A. Buskwofie, G. David-West, and C. A. Clare, “A review of cervical cancer: incidence and disparities,” Journal of the National Medical Association, vol. 112, no. 2, pp. 229–232, 2020.
View at: Publisher Site | Google Scholar
A. Blum, P. Wang, and J. C. Zenklusen, “SnapShot: TCGA-analyzed tumors,” Cell, vol. 173, no. 2, p. 530, 2018.
View at: Publisher Site | Google Scholar
K. Chaudhary, O. B. Poirion, L. Lu, and L. X. Garmire, “Deep learning-basedmulti-omics integration robustly predicts survival in liver cancer,” Clinical Cancer Research, vol. 24, no. 6, pp. 1248–1259, 2018.
View at: Publisher Site | Google Scholar
S. V. Vasaikar, P. Straub, J. Wang, and B. Zhang, “LinkedOmics: analyzing multi-omics data within and across 32 cancer types,” Nucleic Acids Research, vol. 46, no. D1, pp. D956–d963, 2018.
View at: Publisher Site | Google Scholar
J. E. Bader, K. Voss, and J. C. Rathmell, “Targeting metabolism to improve the tumor microenvironment for cancer immunotherapy,” Molecular Cell, vol. 78, no. 6, pp. 1019–1033, 2020.
View at: Publisher Site | Google Scholar
Q. Chen, L. Hu, D. Huang, K. Chen, X. Qiu, and B. Qiu, “Six-lncRNA immune prognostic signature for cervical cancer,” Frontiers in Genetics, vol. 11, Article ID 533628, 2020.
View at: Publisher Site | Google Scholar
X. Qi, Y. Fu, M. Zhang et al., “An innovative immune score-based prognostic nomogram for patients with cervical cancer,” BioMed Research International, vol. 2020, pp. 1–9, 2020.
View at: Publisher Site | Google Scholar
L. Peng, G. Hayatullah, H. Zhou et al., “Tumor microenvironment characterization in cervical cancer identifies prognostic relevant gene signatures,” PLoS One, vol. 16, no. 4, Article ID e0249374, 2021.
View at: Publisher Site | Google Scholar
Z. Zhao, J. Li, H. Li et al., “Integrative bioinformatics approaches to screen potential prognostic immune-related genes and drugs in the cervical cancer microenvironment,” Frontiers in Genetics, vol. 11, p. 727, 2020.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2023 Kunyu Wang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

341

Downloads

416

Citations