- 1The Affiliated Hospital of Medical School of Ningbo University, Ningbo, Zhejiang, China
- 2First Affiliated Hospital of Wenzhou Medical University, Wenzhou, Zhejiang, China
Reactive oxygen species play a crucial role in the prognosis and tumor microenvironment (TME) of malignant tumors. An ROS-related signature was constructed in gastric cancer (GC) samples from TCGA database. ROS-related genes were obtained from the Molecular Signatures Database. Consensus clustering was used to establish distinct ROS-related subtypes related to different survival and immune cell infiltration patterns. Sequentially, prognostic genes were identified in the ROS-related subtypes, which were used to identify a stable ROS-related signature that predicted the prognosis of GC. Correlation analysis revealed the significance of immune cell iniltration, immunotherapy, and drug sensitivity in gastric cancers with different risks. The putative molecular mechanisms of the different gastric cancer risks were revealed by functional enrichment analysis. A robust nomogram was established to predict the outcome of each gastric cancer. Finally, we verified the expression of the genes involved in the model using RT-qPCR. In conclusion, the ROS-related signature in this study is a novel and stable biomarker associated with TME and immunotherapy responses.
1 Introduction
Gastric cancer (GC) is associated with high morbidity and mortality rates (Smyth et al., 2020). Surgical resection of advanced GC is the only truly effective treatment (Berretta et al., 2016; van Boxel et al., 2019). However, GC patients have a higher rate of recurrence and metastasis after surgery (Catalano et al., 2005). Therefore, investigating the mechanisms that govern the occurrence and metastasis of gastric cancer, as well as identifying diagnostic markers and therapeutic targets, has become a hotspot.
Reactive oxygen species (ROS) are defined as oxygen-containing reactive species, such as superoxide anions (O2-) and hydroxyl radicals (OH-) (Collin, 2019). Normal levels of ROS participate in multiple signal transduction pathways and control cell proliferation, differentiation, and growth (Mittler, 2017). It appears that a redox imbalance in which ROS production exceeds the cellular capacity to scavenge ROS leads to oxidative stress (Pizzino et al., 2017). This phenomenon has been implicated in the development of GC. In previous studies, ROS have been found to influence GC progression through autophagy (Liu J. Z. et al., 2020; Cao et al., 2020). However, ROS-related genes have not yet been comprehensively identified as playing a significant role in GC.
NOS3 (endothelial NOS, eNOS) is an enzyme belonging to the family of nitric oxide synthases that mainly generates nitric oxide (NO) (Tenopoulou and Doulias, 2020). NO is involved in ROS-mediated malignancy and normal tissue damage (Nissanka and Moraes, 2018). Previous research has shown that NOS3 inhibits apoptosis and promotes tumor proliferation, invasiveness, and immunosuppression (Sun et al., 2020). Therefore, we proposed that the ROS-related genes can be used to predict tumor survival and immunotherapy response.
In our study, based on a systematic investigation of the expression level and clinical characteristics of ROS-related genes, we constructed a prognostic model for these genes in GC patients. In addition, the correlation between these genes and immune cell infiltration was also demonstrated. The prognostic value of these genes can predict tumor immune microenvironment in GC. Furthermore, clinical and immunological characterization of ROS-related genes will be beneficial for the optimization of immunotherapy in GC.
2 Methods and materials
2.1 Extraction and processing of gastric cancer data
The whole steps performed in the analysis were displayed in the Supplementary Figure S1 The transcription profile and clinical characteristics of gastric cancer were extracted from The Cancer Genome Atlas (TCGA) (https://portal.gdc.cancer.gov/) and the National Center for Biotechnology Information Gene Expression Omnibus (GEO) databases (https://www.ncbi.nlm.nih.gov/geo/). The GEO datasets that contain more than 200 human samples and have complete expression data and clinical information were selected as independent test cohorts. Finally, the GSE84433 as well as the GSE84437 which is a SuperSeries that composed of the GSE84426 and GSE84433 Subseries were chosen. The TCGA-STAD project containing 375 tumor and 32 normal tissues, GSE84437 containing 433 human samples, and GSE84433 containing 357 human tissues were selected for further analysis (Yoon et al., 2020). The transcription profile was converted into fragments per kilobase million (FPKM), and the batch effect of the data in GSE84437 was eliminated. The clinical information of TCGA and GEO datasets were presented in Supplementary File S2. The transcriptomic data of TCGA-STAD was downloaded using the “TCGAbiolinks” R package. The expression difference between tumor and normal tissues was identified using the “limma” algorithm and presented with a heatmap and a volcano plot (Ritchie et al., 2015). The log2 fold-change (log2FC) > 1 and false discovery rate (FDR) < 0.05 were set as the cut-off to screen the differently expressed genes (DEGs). ROS-related genes (Supplementary File S1) were extracted from the Molecular Signatures Database (MSigDB) database (http://www.gsea-msigdb.org/gsea/index.jsp).
2.2 Inner correlation between ROS-related genes
A protein-protein interaction (PPI) network of ROS-related genes was constructed using the STRING database (http://string-db.org/), revealing their putative connections. Then, the correlation between these genes and their mRNA expression was further explored by R packages “igraph” and “reshape2 (Csardi and Nepusz, 2006)”. Red lines show positive correlations, whereas blue lines indicate negative correlations.
2.3 Identification and verification of ROS-related clusters
First, the consensus clustering was conducted with the expression of prognostic ROS-related genes from the “limma” as well as “ConsensusClusterPlus” R packages (Wilkerson and Hayes, 2010; Ritchie et al., 2015). Through consistent cluster analysis, we obtained a cluster Max K value = 9 and cluster Num = 2. Then, K-M analysis of ROS-related clusters was conducted using the “survival” and “survminer” R packages (Rich et al., 2010). Additionally, a heatmap was applied to show the relationship between clusters and age, sex, tumor grade, clinical and T, N, and M stages, as well as the expression profiles of these genes in different clusters, demonstrating the success of this clustering. Sequentially, the immune cell infiltrations between different clusters were evaluated utilizing the “MCPcounter” R package, including B linage, CD8+ T cell, cytotoxic lymphocyte, fibroblasts, T cell, NK cell, monocytic linage, myeloid dendritic cell, and neutrophils (Becht et al., 2016).
2.4 Establishment of ROS-related signature
TTCGA-STAD was used as the training set. The GSE84437 and GSE84433 were used as the validation sets. DEGs were identified in ROS-related clusters. Then, least absolute shrinkage and selection operator (LASSO) Cox regression was applied with the DEGs in the TCGA-STAD project to obtain a ROS-related signature (Tibshirani, 1997). The risk score was derived using the formula =
2.5 Development and validation of a ROS-related nomogram
Uni- and multivariable Cox regression analyses were performed on TCGA-STAD and GSE84437 datasets to identify independent prognostic factors. Sequentially, the relationship between the risk score and clinical characteristics was presented in the heatmap. Then we constructed a prognostic nomogram using the “rms” R package (Iasonos et al., 2008). Calibration plots were used to assess the discriminative power of the nomogram (Van Calster et al., 2016). The accuracy of this nomogram was verified using ROC and decision curves (Mandrekar, 2010; Kerr et al., 2016).
2.6 Subgroup analyses and immune cell infiltration pattern
The limma algorithm was used to reveal the distinct risk of GC in different subgroups, including age, sex, tumor grade, clinical and T, N, and M stages, as well as immune subtypes (Ritchie et al., 2015). Moreover, Survival differences in the low- and high-risk subgroups were established using the K-M analysis. The ImmuneScore, StromalScore, and EstimateScore were assessed by “estimation of stromal and immune cells in malignant tumor tissues using expression data (ESTIMATE)” algorithm, as well as immune cell infiltrations were evaluated using the “tumor immune estimation resource (TIMER),” “cell-type identification by estimating relative subsets of RNA transcripts (CIBERSORT),” “CIBERSORT-absolute mode (ABS),” “QUANTISEQ,” “microenvironment cell populations-counter (MCPCOUNTER),” “XCELL,” and “estimating the proportion of immune and cancer cells (EPIC)” algorithms (Becht et al., 2016; Aran et al., 2017; Li et al., 2017; Chen et al., 2018; Plattner et al., 2020). In addition, single-sample gene set enrichment analysis (ssGSEA) and immune function analyses have been conducted for different risk patients (Zhu L. et al., 2021).
2.7 Prediction of immunotherapy response
The expression of immune-related genes was assessed using the “limma” algorithm. The tumor mutation burden (TMB), microsatellite instability (MSI), and tumor immune dysfunction and exclusion (TIDE) scores were then calculated. TMB score was conducted by R package “maftools” (Mayakonda et al., 2018), MSI score was obtained from previous research (Bonneville et al., 2017), and TIDE score was conducted by online database (http://tide.dfci.harvard.edu/). Finally, the immunophenoscore (IPS) of the different risks was evaluated in different subgroups. IPS refers to four main parts (effector cells, immunosuppressive cells, MHC molecules, and immunomodulators) determining the immunogenicity, and is calculated without bias using machine learning methods. The IPS of STAD patients were downloaded from The Cancer Immunome Atlas (TCIA) (https://tcia.at/home).
2.8 Correlation between ROS-related signature and gene mutation and drug sensitivity
To further explore the relevant factors of the ROS-related signature, we established the risk score in distinct wild and mutated GC types, including ARID1A, CSMD3, FAT4, FLG, LRP1B, MUC16, SYNE1, TP53, and TTN. Afterward, the relationship between signature and 5-fuorouracil, gemcitabine, cytarabine, dasatinib, etoposide, GSK690693, masitinib, and tipifarnib were evaluated according to inhibitory concentration (IC50) by the “pRRophetic” R package (Geeleher et al., 2014) based on the Genomics of Drug Sensitivity in Cancer database (GDSC, https://www.cancerrxgene.org/).
2.9 Functional enrichment analysis
Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed to evaluate the functions of DEGs in low- and high-risk GC (Kanehisa and Goto, 2000; The Gene Ontology Consortium, 2019). Gene set enrichment analysis (GSEA) was applied between different risk groups. The threshold of GSEA analysis was log2FC > 1 and p < 0.05, and the top five enriched pathways were identified (Subramanian et al., 2005). The “c2. cp.kegg.v7.4. symbols. gmt” file downloaded from the GSEA database (https://www.gseamsigdb.org/gsea/index.jsp). Moreover, gene set variation analysis (GSVA) was performed with reference to the methods of previous research, and the enriched functions of different risk groups were visualized using a heat map (Hänzelmann et al., 2013). The enrichment analyses were performed by the “limma” and “GSVA” R packages (Yu et al., 2012; Ritchie et al., 2015).
2.10 Cell culture
HGC-27 and GES-1 cells were obtained from the Cell Bank of Shanghai Institute of Biochemistry and Cell Biology (Shanghai, China). The cell lines were cultured in DMEM or RPMI-1640 supplemented with 10% fetal bovine serum (FBS), 100 μg/mL streptomycin, and 100 U/mL penicillin (Gibco).
2.11 Clinical sample
All samples were obtained with the approval of the Ethics Committee of Ningbo First Hospital. A total of 4 GC patients who signed the informed consents were recruited from the Ningbo First Hospital. The tumor samples and normal tissues were collected from the surgical specimen of the patients.
2.12 RT-qPCR
Total RNA was extracted using TRIzol reagent (Invitrogen) and reverse-transcribed into cDNA templates. β-Actin was used as an endogenous reference. The primer sequences used for amplification are listed in Supplementary File S3.
2.13 Statistical analysis
All statistical analyses were performed using the R software (version 4.1.3) (Supplementary File S4). Group comparisons were performed using the Wilcoxon test or Student’s t-test. Correlations between two continuous variables were assessed using Spearman’s correlation analysis. The Kruskal-Wallis test was used to compare the three groups. Adjusted p-value was calculated using Benjamini-Hochberg FDR. Statistical significance was set at p < 0.05.
3 Results
3.1 Expression and PPI network of ROS-related gene
First, 87 ROS-related genes were obtained from the GSEA database. Using the “limma” package, 31 differentially expressed ROS-related genes were identified (Figures 1A, B). The PPI network was constructed to identify ROS-related DEG interactions (Figure 1C). Meanwhile, Figure 1D shows the correlation network of DEGs, in which different colors indicate different correlation coefficients.
FIGURE 1. ROS-related genes in GC. (A) A heatmap of differentially expressed ROS-related genes in GC. (B) A volcano plot of ROS-related genes in GC. (C) A protein-protein interaction network. (D) A correlation network of ROS-related genes.
3.2 Identification of ROS-related molecular subtypes
The effect of ROS-related genes on GC samples was assessed using consensus clustering analysis. When k = 2, the cumulative distribution curve was most horizontal in the middle section and the heatmap of clustering showed a clear edge. Therefore, GC samples in TCGA-STAD were grouped into two clusters based on ROS gene expression (Figure 2A). The overall survival (OS) rates of the C1 and C2 clusters were significantly different, and the prognosis of the C1 cluster was worse than that of the C2 cluster (Figure 2B). Besides, the K-M survival curve of progression free survival (PFS) also demonstrated the poorer outcome of C1 cluster (Supplementary Figure S2). Figure 2C shows the gene expression profiles between the C1 and C2 clusters. We also analyzed differences in immune cell infiltration between the two clusters and found more immune cell infiltration in the C1 cluster, including B lineage, CD8 T cells, cytotoxic lymphocytes, T cells, NK cells, neutrophils, fibroblasts, monocytic lineage, and myeloid dendritic cells, than in the C2 cluster (Figure 3).
FIGURE 2. Consensus clustering of ROS-related genes in GC. (A) Clustering heatmap for k = 2. (B) K-M curve of C1 and C2. (C) A heatmap for correlation of clusters with clinical characteristics.
3.3 Establishment and validation of the signature
A 4-gene prognosis model including GPX3, DUSP1, NOS3, and TCIRG1 was constructed using LASSO Cox regression (Figures 4A, B). All of the signature genes were identified prognostic genes by univariable Cox analysis, including GPX3, DUSP1, NOS3, and TCIRG1 (Figure 4C). We also analyzed mutations in these genes in GC (Figure 4D). GC patients were divided into high- and low-risk groups according to the median risk score in the training cohort (Figures 4E, F) and test cohort (Figures 4G, H). Patients in the high-risk group had a higher risk of death and poorer prognosis than those in the low-risk group. As shown in Figures 5A, B, PCA demonstrated clear distinctions between the two risk groups.
FIGURE 4. Development of a ROS-related signature in GC. (A) A univariable Cox regression analysis of genes constructing the signature. (B) A waterfall of mutation profile of genes constructing the ROS-related signature. (C) LASSO regression of prognostic genes of cluster 1 and 2. (D) Cross-validation for tuning parameter selection. (E) Risk score of GC in TCGA-STAD dataset. (F) Correlation of survival status with risk score in TCGA-STAD dataset. (G) Risk score of GC in GSE84437 dataset. (H) Correlation of survival status in GSE84437 dataset.
FIGURE 5. Verification of the stability of ROS-related signature. (A) PCA analysis of low- and high-risk GC for TCGA-STAD and GSE84437 datasets. (B) t-SNE analysis for TCGA-STAD and GSE84437 datasets. (C) K-M curve for different risk GC in TCGA-STAD dataset. (D) K-M curve for patients in GSE84437 dataset. (E) ROC analysis of the ROS-related signature across TCGA-STAD project. (F) ROC analysis of ROS-related signature across GSE84437 dataset.
Low-risk patients were linked to better prognoses based on survival curves (Figure 5C) in the TCGA cohort. Similarly, in the test cohort, patients with lower-risk scores showed better survival (Figure 5D and Supplementary Figure S3). Additionally, the area under the curve (AUC) value indicated that our model had a good predictive power (Figures 5E, F). According to Cox regression analyses, the risk score was regarded as an independent prognostic factor for the training and testing cohorts (Figures 6A–D). The landscape of the four signature genes in the training group is shown in the heatmap (Figure 6E). The comparation analysis indicated that the stability of our signature is better than signatures of Mak et al., Zhou et al., and Liu et al. (C-index, 0.614 to 0.568, 0.597, and 0.607) (Supplementary Figure S4). To explore the genomic alterations in low- and high-risk GC, gene mutation analysis was conducted. Result showed that the TMB of low-risk GC is significantly higher than high-risk GC. Besides, some vital cancer-related genes such as TTN, TP53, MUC16, and ARID1A more frequently mutated in low-risk GC (Supplementary Figure S5).
FIGURE 6. Independently predictive ability of the signature. (A, B) Univariable and Multivariable analysis across TCGA dataset. (C) Univariable analysis across GSE84437 dataset. (D) Multivariable analysis of the signature across GSE84437 dataset. (E) A heatmap for the correlation of the signature with clinical characteristics.
A nomogram was developed to predict survival rates in GC (Figure 7A). The predicted results were in good agreement with the actual results according to the calibration plots (Figure 7B). Moreover, the nomogram was regarded as an independent prognostic factor (Figures 7C,D). The AUC of the nomogram was 0.747 (Figure 7E).
FIGURE 7. Establishment of a nomogram with ROS-related signature and clinical characteristics. (A) A nomogram contains ROS-related signature, age, gender, tumor grade, clinical and T, N and M stage. (B) A calibration plot for assessing the predictive power of nomograms at 1-, 3-, and 5-year. (C) Univariable analysis of the nomogram across TCGA dataset. (D) Multivariable analysis across TCGA dataset. (E) DCA analysis of the nomogram model.
3.4 Correlation of prognostic model with clinical features
Next, the correlation between the risk scores and clinical features was studied. Different subgroups, including age, grade, T and TNM stage, as well as immune subtype, had significantly distinct risk scores (Figure 8). In addition, subgroup analysis showed differences in survival between different risk groups in different GC subgroups, including age >65 years, female, male, M0, age ≤65 years, Grade III, stage III-IV, N1-3 and T3+4, further verifying the reliability of this model (Figure 9).
3.5 Correlation of the model with immune activity
The immune microenvironment plays a critical role in the development of GC. In our study, we found that patients in the high-risk group had higher immune, stromal, and ESTIMATE scores (Figure 10A). Next, the differences in immune cell infiltration between the different risk groups were assessed (Figure 10B). Infiltration of Monocytes, Macrophages M2, Mast cells resting, naïve B cells, and eosinophils was more abundant in the high-risk group, while infiltration of CD4 memory activated T cells, resting NK cells resting and Macrophages M1 was more abundant in the low-risk group (Figure 10C). In addition, APC co-stimulation, CCR, and Type II IFN responses were usually more significant in the high-risk group (Figure 10D).
FIGURE 10. Immune cell infiltration pattern in different risk GC. (A) Tumor purity analysis. (B) Immune cell infiltration of different risk GC. (C) ssGSEA analysis. (D) Immune function analysis of different risk GC.
Immune checkpoint molecules are also involved in the occurrence and development of GC. We found that the expression of most checkpoint genes was higher in the high-risk group (Figures 11A, B). TIDE, TMB, and MSI were used to predict the response to tumor immunotherapy. As shown in Figure 11C, the TIDE, exclusion, and dysfunction scores were all higher in the high-risk group than in the low-risk group. Moreover, the TMB and MSI scores were higher in the low-risk group (Figures 11D, E). Furthermore, we found that patients in the low-risk group responded better to immunotherapy (Figure 11F).
FIGURE 11. Capacity of immunotherapy response of the ROS-related signature. (A, B) Expression difference of immune checkpoint in different risk GC. (C) TIDE score. (D) TMB score. (E) MSI of low- and high-risk GC. (F) IPS of low- and high-risk GC in different subgroups.
3.6 Correlation of prognostic model with genetic mutations and predict chemotherapy drug sensitive
Gene mutations are crucial factors in tumor development. In the current study, the correlation between prognostic models and genetic mutations was analyzed. We found that the risk score was lower in the mutation group than in the wild-type group (Figure 12). Next, we predicted the chemoresponses of the subgroups to common chemotherapeutics (Figure 13). The results revealed that high-risk patients were more sensitive to 5-Fluorouracil, Gemcitabine, cytarabine, dasatinib, etoposide, GSK690693, masitinib, and tipipifarnib, suggesting that these patients could benefit from these chemotherapeutic agents.
3.7 Functional enrichment analysis
To further explore the functional enrichment of the signature, we performed an enrichment analysis of DEGs between the low- and high-risk groups. GO enrichment indicated that DEGs were enriched in “muscle system process,” “collagen−containing extracellular matrix,” and “extracellular matrix structural constituent” (Figures 14A, C). KEGG enrichment analysis showed significant enrichment of the “PI3K-AKT” and “cGMP-PKG” signaling pathways (Figures 14B, D). Additionally, from the heatmap of GSVA, significant difference of enriched functions between different risk groups were presented (Figure 14E).
FIGURE 14. Functional enrichment analysis. (A. C) GO analysis of DEGs between low- and high-risk GC. (B, D) KEGG analysis. (E) GSVA analysis of low- and high-risk GC.
3.8 RT-qPCR in GC
As shown in Figure 15A, the expression of GPX3, DUSP1, NOS3, and TCIRG1 was explored using the GEPIA2 site (http://gepia2.cancer-pku.cn/#index). Low expression of GPX3 and DUSP1 was observed in GC tissues as compared to that in normal tissues, whereas NOS3 and TCIRG1 were highly expressed in GC tissues. RT-qPCR was performed to verify the expression levels of these genes. The results demonstrated that NOS3 and TCIRG1 were highly expressed in tumor cells and samples compared to normal cells and tissues, while GPX3 and DUSP1 were expressed at low levels in tumor cells and samples (Figures 15B, C).
FIGURE 15. Expression of GPX3, DUSP1, NOS3, and TCIRG1. (A) Expression of DUSP1, GPX3, NOS3, and TCIRG1 in GEPIA2 online set. (B) RT-qPCR detected the expression of DUSP1, GPX3, NOS3, and TCIRG1 in gastric cancer cell and normal cell. (C) RT-qPCR detected the expression of DUSP1, GPX3, NOS3, and TCIRG1 in gastric cancer samples and normal tissues.
4 Discussion
With advances in high-throughput sequencing technologies, an increasing number of disease targets and prognostic biomarkers have been identified. However, reactive oxygen species-related prognostic biomarkers for GC remain limited. More and more researches have shown that elevated ROS levels are closely linked to the occurrence and development of cancer (Wang et al., 2019). A prognostic model based on these genes was constructed to elucidate the roles of ROS-related genes in GC. Next, the correlation of the prognostic model with immune-infiltrating cells was demonstrated. In addition, the model showed excellent predictive performance for the immune microenvironment of patients with GC. Furthermore, the clinical and immunological characteristics of the constructed prognostic model for GC may help optimize tumor immunotherapy.
ROS refers to a class of oxygen-containing compounds converted from molecular oxygen, whose chemical properties are more active than molecular oxygen (Wu and Cederbaum, 2003). Living organisms have a complex oxidative-antioxidant system (Birben et al., 2012). Under normal circumstances, ROS maintained within a stable range can play a positive role in anti-inflammatory and antibacterial activities (Zhu et al., 2021). When this balance fails to be maintained, ROS increase, promoting cell transformation and lead to the occurrence of malignant tumors (Kim et al., 2016). Studies have shown that ROS can disrupt mitochondrial function and promote oxidative stress, leading to gastric carcinogenesis (Lee and Yu, 2016).
Our 4-gene (GPX3, DUSP1, NOS3, TCIRG1) signature based on ROS-related genes was used for GC. Previous studies have reported that these genes play crucial roles in various cancers, including GC. For example, GPX3 inhibits the growth and spread of GC by regulating methylation and ROS (Xie et al., 2022). Cheng et al. found that DUSP1 activates the MAPK pathway in GC, leading to apatinib resistance (Teng et al., 2018). Moreover, Zhang et al. revealed that the expression level of OS3 was increased and significantly associated with poor patient prognosis in GC (Zou et al., 2021). In this study, we found that high-risk patients had a poorer prognosis.
ROS are not only associated with tumorigenesis and development but also with immune checkpoint inhibitors. The study has shown that ROS can induce the expression of PD-L1 by regulating JAK/STAT3 pathway (Liao et al., 2021). In our study, we also found that different risk groups had different immune checkpoint expression, and the survival time of patients with high expression of TGFB1, TGFBR1, CD27, BTLA, CD40, and IL-10 in the high-risk group was shorter than that in the low-risk group. In addition, we found that the constructed prognostic model could predict the response to immunotherapy, with patients in the high-risk group having a poorer response to immunotherapy.
The occurrence of tumors is an extremely complicated biological process, and there are various non-tumor cells in its microenvironment, such as cancer-associated fibroblasts (CAF), tumor-associated macrophages (TAM), immune cells and various lymphatic vessels, tissue fluid, and cytokines, which can maintain tumor growth and increase tumor heterogeneity, adaptability, and metastasis, eventually leading to the development of malignant tumors (Chen et al., 2021; Laha et al., 2021; Raskov et al., 2021; Zhou et al., 2021). Our study found a higher degree of macrophage M2 infiltration in the high-risk group, which may contribute to its poor prognosis. M2 macrophages are involved in tumorigenesis, growth, invasion, and metastasis, and often promote tumor cell growth, angiogenesis, and migration (Zhao et al., 2020). Numerous studies have shown that TAM present in tumor tissues usually exhibit an M2-like phenotype, which is closely associated with cancer treatment and prognosis (Liu et al., 2020; Pan et al., 2020).
GO and KEGG analyses based on DEGs between different risk groups were performed to further elucidate the underlying mechanisms of ROS-related genes. GO analysis revealed that the DEGs were enriched in “collagen−containing extracellular matrix,” “muscle system process,” and “extracellular matrix structural constituent.” KEGG showed significant enrichment of “PI3K-AKT” and “cGMP-PKG” signaling pathway. The PI3K-AKT signaling pathway is involved in the occurrence and development of various cancers, and can promote the proliferation and invasion of cancer cells, inhibit apoptosis, and promote tumor angiogenesis, thereby leading to the progression of malignant tumors (Alzahrani, 2019; Ediriweera et al., 2019; Noorolyai et al., 2019).
In addition, we evaluated the effectiveness of some common chemotherapeutic agents in different risk groups and found that patients in the high-risk group were more sensitive to eight chemotherapeutic agents (5-fluorouracil, gemcitabine, cytarabine, dasatinib, etoposide, GSK690693, masitinib, and tipifarnib), indicating that high-risk patients may benefit from these eight chemotherapy drugs.
In general, we have provided the most comprehensive elucidation of ROS-related genes in GC. We first constructed an ROS-related prognostic model for GC, confirmed the relationship between the model and immune checkpoint genes, as well as immune infiltration, and predicted the response of different risk patients to immunotherapy and chemotherapy. It provides a screening tool for the diagnosis and prognosis of gastric cancer, and a new way to dissect the relationship between gastric cancer and immunity. Our tool contained only four genes and had good stability. Comparing to the clinical prediction tools such as 21-gene score and 70-gene for breast cancer, the Ros-related signature had better economic viability and may use to improve the treatment of clinical GC in the future.
Data availability statement
The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.
Ethics statement
The studies involving human participants were reviewed and approved by Ethics Committee of Ningbo First Hospital (2022068A01). The patients/participants provided their written informed consent to participate in this study.
Author contributions
YD and YG. designed this study and drafted the manuscript. YM, ZW, and KH wrote the manuscript. KC and YG searched for publications and collected the data. YD and KH analyzed the data. The final paper was reviewed by all authors and approved for publication.
Funding
This study was supported by the Ningbo University Institute of Geriatrics (LNBYJS-2021).
Acknowledgments
We thank all members who contributed to this study.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2023.1074900/full#supplementary-material
SUPPLEMENTARY FILE S1 | The ROS-related genes extracted from the MSigDB.
SUPPLEMENTARY FILE S2 | The clinical information of TCGA and GEO datasets.
SUPPLEMENTARY FILE S3 | The primer sequences used for RT-qPCR.
SUPPLEMENTARY FILE S4 | Original R code of the analysis.
SUPPLEMENTARY FIGURE S1 | The flowchart of the analysis.
SUPPLEMENTARY FIGURE S2 | The PFS analysis between C1 and C2 clusters.
SUPPLEMENTARY FIGURE S3 | The K-M survival analysis between low- and high-risk GC.
SUPPLEMENTARY FIGURE S4 | The comparation analysis between ROS-related signature and other three models.
SUPPLEMENTARY FIGURE S5 | The mutation profile of low- and high-risk GC.
References
Alzahrani, A. S. (2019). PI3K/Akt/mTOR inhibitors in cancer: At the bench and bedside. Semin. Cancer Biol. 59, 125–132. doi:10.1016/j.semcancer.2019.07.009
Aran, D., Hu, Z., and Butte, A. J. (2017). xCell: digitally portraying the tissue cellular heterogeneity landscape. Genome Biol. 18, 220. doi:10.1186/s13059-017-1349-1
Becht, E., Giraldo, N. A., Lacroix, L., Buttard, B., Elarouci, N., Petitprez, F., et al. (2016). Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol. 17, 218. doi:10.1186/s13059-016-1070-5
Berretta, S., Berretta, M., Fiorica, F., Di Francia, R., Magistri, P., Bertola, G., et al. (2016). Multimodal approach of advanced gastric cancer: Based therapeutic algorithm. Eur. Rev. Med. Pharmacol. Sci. 20, 4018–4031.
Birben, E., Sahiner, U. M., Sackesen, C., Erzurum, S., and Kalayci, O. (2012). Oxidative stress and antioxidant defense. World Allergy Organ. J. 5, 9–19. doi:10.1097/WOX.0b013e3182439613
Bonneville, R., Krook, M. A., Kautto, E. A., Miya, J., Wing, M. R., Chen, H. Z., et al. (2017). Landscape of microsatellite instability across 39 cancer types. JCO Precis. Oncol. 2017, PO.17.00073. doi:10.1200/PO.17.00073
Cao, Y., Wang, J., Tian, H., and Fu, G. H. (2020). Mitochondrial ROS accumulation inhibiting JAK2/STAT3 pathway is a critical modulator of CYT997-induced autophagy and apoptosis in gastric cancer. J. Exp. Clin. Cancer Res. 39, 119. doi:10.1186/s13046-020-01621-y
Catalano, V., Labianca, R., Beretta, G. D., Gatta, G., de Braud, F., and Van Cutsem, E. (2005). Gastric cancer. Crit. Rev. Oncol. Hematol. 54, 209–241. doi:10.1016/j.critrevonc.2005.01.002
Chen, B., Khodadoust, M. S., Liu, C. L., Newman, A. M., and Alizadeh, A. A. (2018). Profiling tumor infiltrating immune cells with CIBERSORT. Methods Mol. Biol. 1711, 243–259. doi:10.1007/978-1-4939-7493-1_12
Chen, D., Zhang, X., Li, Z., and Zhu, B. (2021). Metabolic regulatory crosstalk between tumor microenvironment and tumor-associated macrophages. Theranostics 11, 1016–1030. doi:10.7150/thno.51777
Cieslak, M. C., Castelfranco, A. M., Roncalli, V., Lenz, P. H., and Hartline, D. K. (2020). t-Distributed Stochastic Neighbor Embedding (t-SNE): A tool for eco-physiological transcriptomic analysis. Mar. Genomics. 51, 100723. doi:10.1016/j.margen.2019.100723
Collin, F. (2019). Chemical basis of reactive oxygen species reactivity and involvement in neurodegenerative diseases. Int. J. Mol. Sci. 20, 2407. doi:10.3390/ijms20102407
Csardi, G., and Nepusz, T. (2006). The igraph software package for complex network research. InterJournal. Complex Syst. 1695, 1–9.
Ediriweera, M. K., Tennekoon, K. H., and Samarakoon, S. R. (2019). Role of the PI3K/AKT/mTOR signaling pathway in ovarian cancer: Biological and therapeutic significance. Semin. Cancer Biol. 59, 147–160. doi:10.1016/j.semcancer.2019.05.012
Geeleher, P., Cox, N., and Huang, R. S. (2014). pRRophetic: an R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLOS ONE 9, e107468. doi:10.1371/journal.pone.0107468
Hänzelmann, S., Castelo, R., and Guinney, J. (2013). GSVA: Gene set variation analysis for microarray and RNA-seq data. BMC Bioinforma. 14, 7. doi:10.1186/1471-2105-14-7
Iasonos, A., Schrag, D., Raj, G. V., and Panageas, K. S. (2008). How to build and interpret a nomogram for cancer prognosis. J. Clin. Oncol. 26, 1364–1370. doi:10.1200/JCO.2007.12.9791
Kanehisa, M., and Goto, S. (2000). KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 28, 27–30. doi:10.1093/nar/28.1.27
Kerr, K. F., Brown, M. D., Zhu, K., and Janes, H. (2016). Assessing the clinical impact of risk prediction models with decision curves: Guidance for correct interpretation and appropriate use. J. Clin. Oncol. 34, 2534–2540. doi:10.1200/JCO.2015.65.5654
Kim, J., Kim, J., and Bae, J. S. (2016). ROS homeostasis and metabolism: A critical liaison for cancer therapy. Exp. Mol. Med. 48, e269. doi:10.1038/emm.2016.119
Laha, D., Grant, R., Mishra, P., and Nilubol, N. (2021). The role of tumor necrosis factor in manipulating the immunological response of tumor microenvironment. Front. Immunol. 12, 656908. doi:10.3389/fimmu.2021.656908
Lee, C. H., and Yu, H. S. (2016). Role of mitochondria, ROS, and DNA damage in arsenic induced carcinogenesis. Front. Biosci. Sch. Ed. 8, 312–320. doi:10.2741/s465
Li, T., Fan, J., Wang, B., Traugh, N., Chen, Q., Liu, J. S., et al. (2017). TIMER: A web server for comprehensive analysis of tumor-infiltrating immune cells. Cancer Res. 77, e108–e110. doi:10.1158/0008-5472.CAN-17-0307
Liao, Y. C., Wu, S. Y., Huang, Y. F., Lo, P. C., Chan, T. Y., Chen, C. A., et al. (2021). NOX2-deficient neutrophils facilitate joint inflammation through higher pro-inflammatory and weakened immune checkpoint activities. Front. Immunol. 12, 743030. doi:10.3389/fimmu.2021.743030
Liu, J. Z., Hu, Y. L., Feng, Y., Jiang, Y., Guo, Y. B., Liu, Y. F., et al. (2020). BDH2 triggers ROS-induced cell death and autophagy by promoting Nrf2 ubiquitination in gastric cancer. J. Exp. Clin. Cancer Res. 39, 123. doi:10.1186/s13046-020-01620-z
Liu, Q., Yang, C., Wang, S., Shi, D., Wei, C., Song, J., et al. (2020). Wnt5a-induced M2 polarization of tumor-associated macrophages via IL-10 promotes colorectal cancer progression. Cell Commun. Signal. 18, 51. doi:10.1186/s12964-020-00557-2
Liu, Y., Wu, M., Cao, J., Zhu, Y., Ma, Y., Pu, Y., et al. (2022). Identification and verification of a glycolysis-related gene signature for gastric cancer. Ann. Transl. Med. 10, 1010. doi:10.21037/atm-22-3980
Mak, T. K., Li, X., Huang, H., Wu, K., Huang, Z., He, Y., et al. (2022). The cancer-associated fibroblast-related signature predicts prognosis and indicates immune microenvironment infiltration in gastric cancer. Front. Immunol. 13, 951214. doi:10.3389/fimmu.2022.951214
Mandrekar, J. N. (2010). Receiver operating characteristic curve in diagnostic test assessment. J. Thorac. Oncol. 5, 1315–1316. doi:10.1097/JTO.0b013e3181ec173d
Mayakonda, A., Lin, D. C., Assenov, Y., Plass, C., and Koeffler, H. P. (2018). Maftools: Efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 28, 1747–1756. doi:10.1101/gr.239244.118
Nissanka, N., and Moraes, C. T. (2018). Mitochondrial DNA damage and reactive oxygen species in neurodegenerative disease. FEBS Lett. 592, 728–742. doi:10.1002/1873-3468.12956
Noorolyai, S., Shajari, N., Baghbani, E., Sadreddini, S., and Baradaran, B. (2019). The relation between PI3K/AKT signalling pathway and cancer. Gene 698, 120–128. doi:10.1016/j.gene.2019.02.076
Pan, Y., Yu, Y., Wang, X., and Zhang, T. (2020). Tumor-associated macrophages in tumor immunity. Front. Immunol. 11, 583084. doi:10.3389/fimmu.2020.583084
Pizzino, G., Irrera, N., Cucinotta, M., Pallio, G., Mannino, F., Arcoraci, V., et al. (2017). Oxidative stress: Harms and benefits for human health. Oxid. Med. Cell. Longev. 2017, 8416763. doi:10.1155/2017/8416763
Plattner, C., Finotello, F., and Rieder, D. (2020). Deconvoluting tumor-infiltrating immune cells from RNA-seq data using quanTIseq. Methods Enzymol. 636, 261–285. doi:10.1016/bs.mie.2019.05.056
Raskov, H., Orhan, A., Gaggar, S., and Gögenur, I. (2021). Cancer-associated fibroblasts and tumor-associated macrophages in cancer and cancer immunotherapy. Front. Oncol. 11, 668731. doi:10.3389/fonc.2021.668731
Rich, J. T., Neely, J. G., Paniello, R. C., Voelker, C. C., Nussenbaum, B., and Wang, E. W. (2010). A practical guide to understanding Kaplan-Meier curves. Otolaryngol. Head. Neck Surg. 143, 331–336. doi:10.1016/j.otohns.2010.05.007
Ringnér, M. (2008). What is principal component analysis? Nat. Biotechnol. 26, 303–304. doi:10.1038/nbt0308-303
Ritchie, M. E., Phipson, B., Wu, D., Hu, Y., Law, C. W., Shi, W., et al. (2015). Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43, e47. doi:10.1093/nar/gkv007
Smyth, E. C., Nilsson, M., Grabsch, H. I., van Grieken, N. C., and Lordick, F. (2020). Gastric cancer. Lancet 396, 635–648. doi:10.1016/S0140-6736(20)31288-5
Subramanian, A., Tamayo, P., Mootha, V. K., Mukherjee, S., Ebert, B. L., Gillette, M. A., et al. (2005). Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. U. S. A. 102, 15545–15550. doi:10.1073/pnas.0506580102
Sun, Y., Zhao, Z., Zhang, H., Li, J., Chen, J., Luan, X., et al. (2020). The interaction of lead exposure and CCM3 defect plays an important role in regulating angiogenesis through eNOS/NO pathway. Environ. Toxicol. Pharmacol. 79, 103407. doi:10.1016/j.etap.2020.103407
Teng, F., Xu, Z., Chen, J., Zheng, G., Zheng, G., Lv, H., et al. (2018). DUSP1 induces apatinib resistance by activating the MAPK pathway in gastric cancer. Oncol. Rep. 40, 1203–1222. doi:10.3892/or.2018.6520
Tenopoulou, M., and Doulias, P. T. (2020). Endothelial nitric oxide synthase-derived nitric oxide in the regulation of metabolism. F1000Res. 9, F1000 Faculty Rev-1190. doi:10.12688/f1000research.19998.1
The Gene Ontology Consortium (2019). The gene Ontology resource: 20 years and still GOing strong. Nucleic Acids Res. 47, D330–D338. doi:10.1093/nar/gky1055
Tibshirani, R. (1997). The lasso method for variable selection in the Cox model. Stat. Med. 16, 385–395. doi:10.1002/(sici)1097-0258(19970228)16:4<385:aid-sim380>3.0.co;2-3
van Boxel, G. I., Ruurda, J. P., and van Hillegersberg, R. (2019). Robotic-assisted gastrectomy for gastric cancer: A European perspective. Gastric Cancer 22, 909–919. doi:10.1007/s10120-019-00979-z
Van Calster, B., Nieboer, D., Vergouwe, Y., De Cock, B., Pencina, M. J., and Steyerberg, E. W. (2016). A calibration hierarchy for risk models was defined: From utopia to empirical data. J. Clin. Epidemiol. 74, 167–176. doi:10.1016/j.jclinepi.2015.12.005
Wang, L., Wang, C., Tao, Z., Zhao, L., Zhu, Z., Wu, W., et al. (2019). Curcumin derivative WZ35 inhibits tumor cell growth via ROS-YAP-JNK signaling pathway in breast cancer. J. Exp. Clin. Cancer Res. 38, 460. doi:10.1186/s13046-019-1424-4
Wilkerson, M. D., and Hayes, D. N. (2010). ConsensusClusterPlus: A class discovery tool with confidence assessments and item tracking. Bioinformatics 26, 1572–1573. doi:10.1093/bioinformatics/btq170
Wu, D., and Cederbaum, A. I. (2003). Alcohol, oxidative stress, and free radical damage. Alcohol Res. Health. 27, 277–284.
Xie, J., Fu, L., and Zhang, J. (2022). Analysis of influencing factors on the occurrence and development of gastric cancer in high-incidence areas of digestive tract tumors based on high methylation of GPX3 gene. J. Oncol. 2022, 3094881. doi:10.1155/2022/3094881
Yoon, S. J., Park, J., Shin, Y., Choi, Y., Park, S. W., Kang, S. G., et al. (2020). Deconvolution of diffuse gastric cancer and the suppression of CD34 on the BALB/c nude mice model. BMC Cancer 20, 314. doi:10.1186/s12885-020-06814-4
Yu, G., Wang, L. G., Han, Y., and He, Q. Y. (2012). clusterProfiler: an R package for comparing biological themes among gene clusters. Omics 16, 284–287. doi:10.1089/omi.2011.0118
Zhao, S., Mi, Y., Guan, B., Zheng, B., Wei, P., Gu, Y., et al. (2020). Tumor-derived exosomal miR-934 induces macrophage M2 polarization to promote liver metastasis of colorectal cancer. J. Hematol. Oncol. 13, 156. doi:10.1186/s13045-020-00991-2
Zhou, C., Zhang, Y., Yan, R., Huang, L., Mellor, A. L., Yang, Y., et al. (2021). Exosome-derived miR-142-5p remodels lymphatic vessels and induces IDO to promote immune privilege in the tumour microenvironment. Cell Death Differ. 28, 715–729. doi:10.1038/s41418-020-00618-6
Zhou, X., Zhang, B., Zheng, G., Zhang, Z., Wu, J., Du, K., et al. (2022). Novel necroptosis-related gene signature for predicting early diagnosis and prognosis and immunotherapy of gastric cancer. Cancers (Basel) 14, 3891. doi:10.3390/cancers14163891
Zhu, L., Yang, F., Wang, L., Dong, L., Huang, Z., Wang, G., et al. (2021). Identification the ferroptosis-related gene signature in patients with esophageal adenocarcinoma. Cancer Cell Int. 21, 124. doi:10.1186/s12935-021-01821-2
Zhu, X., Guo, Y., Liu, Z., Yang, J., Tang, H., and Wang, Y. (2021). Itaconic acid exerts anti-inflammatory and antibacterial effects via promoting pentose phosphate pathway to produce ROS. Sci. Rep. 11, 18173. doi:10.1038/s41598-021-97352-x
Keywords: reactive oxygen species, tumor microenvironment, prognosis, gastric cancer, immunotherapy
Citation: Cen K, Wu Z, Mai Y, Dai Y, Hong K and Guo Y (2023) Identification of a novel reactive oxygen species (ROS)-related genes model combined with RT-qPCR experiments for prognosis and immunotherapy in gastric cancer. Front. Genet. 14:1074900. doi: 10.3389/fgene.2023.1074900
Received: 20 October 2022; Accepted: 03 April 2023;
Published: 14 April 2023.
Edited by:
Tao Chen, Tongji University, ChinaReviewed by:
Aimin Jiang, The First Affiliated Hospital of Xi’an Jiaotong University, ChinaJamshid Hadjati, Tehran University of Medical Sciences, Iran
Copyright © 2023 Cen, Wu, Mai, Dai, Hong and Guo. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Kai Hong, aG9uZ2thaTA2MjlAMTYzLmNvbQ==; Yangyang Guo, OTg0MDU0ODYzQHFxLmNvbQ==
†These authors have contributed equally to this work