Yuanyuan Feng1†
Wenkai Wang1†Shasha Jiang1†
Yongming Liu2
Yan Wang1Xiangyang Zhan3Huirong Zhu1*
Guoqing Du2*- 1Department of Oncology, Shuguang Hospital Affiliated to Shanghai University of Traditional Chinese Medicine, Shanghai, China
- 2Shi’s Center of Orthopedics and Traumatology, Shuguang Hospital Affiliated to Shanghai University of Traditional Chinese Medicine, Shanghai, China
- 3Urology Centre, Shuguang Hospital Afliated to Shanghai University of Traditional Chinese Medicine, Shanghai, China
Background: Kidney renal clear cell carcinoma (KIRC) is a common and clinically significant subtype of kidney cancer. A potential therapeutic target in KIRC is disulfidptosis, a novel mode of cell death induced by disulfide stress. The aim of this study was to develop a prognostic model to explore the clinical significance of different disulfidptosis gene typings from KIRC.
Methods: A comprehensive analysis of the chromosomal localization, expression patterns, mutational landscape, copy number variations, and prognostic significance of 10 disulfide death genes was conducted. Patients were categorized into distinct subtypes using the Non-negative Matrix Factorization (NMF) typing method based on disulfidptosis gene expression patterns. Weighted Gene Co-expression Network Analysis (WGCNA) was used on the KIRC dataset to identify differentially expressed genes between subtype clusters. A risk signature was created using LASSO-Cox regression and validated by survival analysis. An interaction between risk score and immune cell infiltration, tumor microenvironment characteristics and pathway enrichment analysis were investigated.
Results: Initial findings highlight the differential expression of specific DRGs in KIRC, with genomic instability and somatic mutation analysis revealing key insights into their role in cancer progression. NMF clustering differentiates KIRC patients into subgroups with distinct survival outcomes and immune profiles, and hierarchical clustering identifies gene modules associated with key biological and clinical parameters, leading to the development of a risk stratification model (LRP8, RNASE2, CLIP4, HAS2, SLC22A11, and KCTD12) validated by survival analysis and predictive of immune infiltration and drug sensitivity. Pathway enrichment analysis further delineates the differential molecular pathways between high-risk and low-risk patients, offering potential targets for personalized treatment. Lastly, differential expression analysis of model genes between normal and KIRC cells provides insights into the molecular mechanisms underlying KIRC, highlighting potential biomarkers and therapeutic targets.
Conclusion: This study contributes to the understanding of KIRC and provides a potential prognostic model using disulfidptosis gene for personalized management in KIRC patients. The risk signature shows clinical applicability and sheds light on the biological mechanisms associated with disulfide-induced cell death.
Introduction
Renal cell carcinoma (RCC) is the second-most common cause of death associated with malignant neoplasms of the urinary tract (Li et al., 2017). A series of possible features have been recognised, such as environmental factors (e.g., smoking, obesity, high blood pressure, long-term use of painkillers) and genetic factors, which are related to the renal cell carcinoma (Hsieh et al., 2017). KIRC is one type of RCC and constitutes 80% of all renal malignancies (Klümper et al., 2020; Chu et al., 2023). Although, the therapeutic efficacy of KIRC has improved dramatically due to rapid advances in medical treatment and research. However, the majority of patients with metastatic disease have a poor prognosis and limited therapeutic advances due to few or no early symptoms of KIRC (Zhuang et al., 2019). In this way, early detection of KIRC contributes to a better prognosis (Li et al., 2019). Understanding KIRC is of paramount importance for both clinical practice and research aimed at improving therapeutic outcomes for patients with this disease, and the identification of clear and reliable novel prognostic molecular biomarkers can accurately predict the prognosis of KIRC patients.
Disulfidptosis, a recently proposed mode of cell death, is a novel phenomenon with important implications for understanding and treating cancer. When cells are starved of glucose and solute carrier family 7 member 11 (SLC7A11) is expressed at high levels, it induces disulfidptosis and initiates the cell death process. The reduction-oxidation (REDOX) reaction and the formation of disulfide bonds are the mechanisms that trigger cell death (Wang et al., 2023). In a groundbreaking study by Liu et al., it was demonstrated that disulfidptosis is induced by disulfide stress and rapid cell death in SLC7A11 high-expression cells, specifically under conditions of glucose starvation (Liu et al., 2023). Abnormal accumulation of cystine in the cytoplasm leads to disulfide stress, forcing the reduction of cystine to the more soluble cysteine, which produces relatively strong cytotoxicity (Liu X. et al., 2021). The reducing of cystine to cysteine requires the consumption of reducing equivalents of nicotinamide adenine dinucleotide phosphate (NADPH), so SLC7A11-expressing cells with a high rate of cystine uptake demand a high level of NADPH to continuously reducing cystine to cysteine to support intracellular homeostasis. NADPH production is impaired by the intracellular accumulation of cystine and other disulfide molecules and rapid cell death in response to glucose starvation, known as disulfidemia (disulfidptosis) (Liu et al., 2020; Yan et al., 2023). Cancer development and progression are closely linked to the formation and cleavage of disulfide bonds. Sulfur (S)-based chemical bonds have been used to develop tumour-specific, redox-responsive DDSs, which include the thioether bond, the disulfide bond and the thioketal bond (Sun et al., 2019). Notably, this unique pattern of cell death was identified in a kidney cancer cell line, suggesting its relevance and possible potential as a treatment target in the context of KIRC (Zhang et al., 2022; Kang et al., 2023). In spite of these ground-breaking findings, the full impact of disulfide on the pathogenesis, progression, and treatment of KIRC remains to a large extent unexplored.
The study was designed to determine the clinical significance and predictive potential of the risk profile based on the expression pattern of the disulfide disulfidptosis genes in patients with KIRC. To gain insight into potential mechanisms of disulfide-induced KIRC cell death, the study also aimed to investigate the association of risk scores with immune cell infiltration, tumor microenvironment characteristics, drug sensitivity, and pathway enrichment analyses.
Methods
TCGA data acquisition and clinical details
RNA-seq transcriptome data and relevant clinical details, including sex, age, subtype, IDH status, and survival information for KIRC, were taken from the Cancer Genome Atlas (TCGA)-KIRC database, available at https://portal.gdc.cancer.gov/. In addition, somatic mutation counts and copy number variation (CNV) data were also obtained from the TCGA repository.
Comprehensive analysis of disulfidptosis-related genes (DRGs)
This study first comprehensively examines the chromosomal localization, differential expression patterns, mutational landscape, copy number variations, and prognostic significance of 10 disulfide death genes in KIRC cells. And the 10 disulfidptosis-associated genes (GYS1, NDUFS1, OXSM, LRPPRC, NDUFA11, NUBPL, NCKAP1, RPN1, SLC3A2, and SLC7A11) was obtained from previous literature (Liu et al., 2023). To accomplish this, the genomic positions of the 10 disulfide death genes on chromosomes were precisely mapped. The chromosome and gene positions were obtained from the NCBI database (https://www.ncbi.nlm.nih.gov/). And the changes in gene expression levels were compared between KIRC cells and healthy controls. Moreover, the mutational profiles and copy number variations of these genes specifically within KIRC cells were investigated. Lastly, the potential prognostic value of the ten genes in KIRC was explored. These data are all from the TCGA-KIRC project and were visualized using ggplot2 R package (version 3.4.4).
NMF typing of clear cell carcinoma patients based on ten disulfide death genes
KIRC patients were classified into distinct subtypes using the NMF typing method, which uses a panel of 10 disulfide death genes. The NMF clustering approach, implemented by the NMF R package (version 0.26), was used for this purpose. The study applied specific criteria to determine the optimal configuration, including the evaluation of the consistency map, the comprehensive coefficient and the silhouette coefficient. In particular, the analysis included ten runs with the number of clusters (k) set to 2. The resulting patient subtypes were defined according to the distribution of 10 disulfide death genes. Kaplan-Meier analysis was then performed to generate survival curves (Xu et al., 2022) by the survival R package (version 3.5–7). The log-rank test was applied to compare disease-specific survival differences between the two risk groups and to assess their significance by the stats R package (Nomura et al., 2020). The differences in immune cell infiltration, tumour microenvironment and stemness index between the two clusters were evaluated by the estimate R package (version 1.0.13).
Differential expression and WGCNA analysis of DEGs between clusters
Differential expression analysis of genes (DEGs) between cluster C1 and cluster C2 was performed using the limma R package (version 3.50.3). This study utilized WGCNA to buid a gene co-expression network and divide genes into five modules, and studied the gene significance and the relationship between the module and phenotype by using the WGCNA R package (version1.72–5) (Gao et al., 2022).
Development and validation of a risk signature for survival prediction in KIRC
Based on univariate Cox analysis, we identified survival-associated genes from the most important module and established a risk signature consisting of six genes (LRP8, RNASE2, CLIP4, HAS2, SLC22A11, and KCTD12) using LASSO-Cox regression with the glmnet R package (version 4.1–8) and survival R package (version 3.5–7). After constructing the risk score, we employed Cox regression models, including univariate and multivariate Cox proportional hazards regression models, to investigate the effect of various factors, including the risk score, on survival time (Xue et al., 2022). The TCGA samples were divided into a training set and a test set to validate our model, and with reference to the median score, we classified the KIRC patients into a low-risk group and a high-risk group. The Kaplan-Meier survival curve was then used to compare the overall survival (OS) and progression-free survival (PFS) of the two groups of patients (Ohashi et al., 2020). The concordance index, which represents the proportion of concordant pairs between predicted and observed outcomes across all patients, was used to evaluate the predictive ability of the model. The time-dependent ROC curve (timeROC) R package (version 0.4) was employed to calculate time-dependent receiver operating characteristic (ROC) curves for sensitivity and specificity based on the signature (Park et al., 2022). Using the rms R package (version 6.7-1), we constructed calibration plots, calibration curve plots, and nomograms for 1-year, 3-year, and 5-year survival based on the risk score and other significant clinical features. These evaluations were conducted to assess the clinical applicability of the risk model.
Estimation of immune and stromal cell abundance in the tumor microenvironment
The Microenvironment Cell Populations-counter (MCP counter) was used to estimate the absolute abundance of eight immune cell types in the tumour microenvironment, including T cells, CD8 T cells, cytotoxic lymphocytes, B lineage, NK cells, monocytic lineage, myeloid dendritic cells, and neutrophils, and two stromal cell types (endothelial cells and fibroblasts) (Horr and Buechler, 2021). CIBERSORT-ABS was employed to calculate the abundance of immune-related cells (Feng et al., 2022). The relationship between risk levels of immune cell infiltration and clusters was investigated using the limma R package. Furthermore, the linkET package was used to estimate the association between risk scores, immune checkpoint-related genes, immune cell infiltration, and immune function using Spearman’s correlation analysis (Wang et al., 2022). Additionally, this study compared the tumor mutation burden and mutation status of genes between high-risk and low-risk groups. The impact of tumor mutation burden on prognosis was explored through Kaplan-Meier (KM) analysis (Zhao et al., 2021).
Investigation of immunophenoscores and drug sensitivity
The Immunophenoscores (IPS) of KIRC patients were collected from the Cancer Immunome Atlas (TCIA, https://tcia.at/home) database, and the correlation of IPS with risk features was investigated using the Wilcoxon rank-sum test (Liu Y. et al., 2021). The oncoPredict R software package (version 0.2) was employed to predict drug susceptibility scores for KIRC. The correlation between risk scores and drug susceptibility scores was then assessed using the limma R package and plotted using the ggplot2 R package to explore the clinical value of this risk model in the management of KIRC.
Gene function and pathway analysis of disulfide-induced cell death
To elucidate the gene functions and biological pathways associated with disulfide-induced cell death, differential expression analysis of annotated genes was conducted using the R packages limma and clusterProfiler (version 4.2.2). Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were performed on the identified differentially expressed genes. Subsequently, gene set enrichment analysis (GSEA) was conducted using the GSEA software (http://www.broadinstitute.org/gsea/index.jsp) on the gene sets from MSigDB (c2. cp.kegg.v7.4. symbols.gmt) version 4.0 (Jia et al., 2023). Furthermore, to investigate the biological functional differences among the groups, Gene Set Variation Analysis (GSVA) enrichment analysis was performed using the GSVA R package (Tang et al., 2023).
Cell culture
786-O cells are derived from tissue of primary clear cell carcinoma, while CCC-HEK-1 cells are derived from normal human embryonic kidney tissue. 786-O cells were grown in complete medium (RPMI-1640 containing 10% fetal bovine serum (FBS) and 100 U mL−1 penicillin-streptomycin) in humidified air with CO2 at 37°C for 24 h. CCC-HEK-1 cells were grown in complete medium (DMEM containing 10% fetal bovine serum (FBS) and 100 U mL−1 penicillin-streptomycin) in humidified air with CO2 at 37°C for 24 h.
Patients’ samples collection
The normal kidney and KIRC samples from patients were obtained from the Shuguang Hospital Affiliated to Shanghai University of Traditional Chinese Medicine, China. The normal kidney tissue samples were obtained from healthy controls (sample size = 3), and KIRC tissue samples were obtained from patients diagnosed with KIRC (sample size = 5). The Ethics Committee of the Shuguang Hospital Affiliated to Shanghai University of Traditional Chinese Medicine approved the protocol for collecting human samples (approval number 2023-1311-78-01).
Quantitative real-time PCR
Total mRNA was extracted from cells and clinic samples using RNA Isolation Kit according to the manufacturer’s instructions. Then cDNA was synthesised using the First Strand cDNA Synthesis Kit according to the manufacturer’s instructions, which was treated as template DNA for the qRT-PCR assay. The qRT-PCR of six module genes (LRP8, RNASE2, CLIP4, HAS2, SLC22A11, and KCTD12) was performed using SYBR Green qPCR Mix on a real-time detector. Each experiment was repeated three times to calculate the mean and standard deviation (SD). The details of reagents and primers for the genes are in Supplementary Table S1.
Data statistical analysis
All data in this study are expressed as mean standard deviation (SD). All bioinformatics analysis is performed using R software (v.4.2.2). p values were calculated using SPSS Statistics 25 (*p < 0.05, **p < 0.01, and ***p < 0.001).
Results
Identification and expression analysis of DRGs
In this investigation, we meticulously analyzed a cohort of 10 DRGs to elucidate their role in KIRC. Chromosomal mapping, as depicted in Figure 1A, delineates the precise genomic locations of these genes. Notably, expression profiling reveals a distinct downregulation of NDUFA11 in tumor tissues, in contrast to the upregulation observed for RPN1; the expression of the other genes remained statistically unaltered between normal and tumor samples (Figure 1B). Our study further extends into the realm of genomic instability, where Figure 1C highlights the prevalence of CNV-related mutations across these DRGs. In this context, OXSM is notably predisposed to frequent CNV deletions, whereas the remaining genes exhibit a spectrum of CNV deletions and amplifications. Additionally, Figure 1D provides an insightful overview of the somatic mutation landscape in these genes within KIRC, illustrating a uniformly low mutation rate. Complementing these findings, Figure 1E uses network visualisation techniques to intricately map the connections between DRGs and their prognostic significance in KIRC, providing a comprehensive and multifaceted view of their biological relevance.
FIGURE 1. Identification and Expression Analysis of DRGs. (A) Chromosomal mapping. (B) Expression profiling. (C) CNV-related mutations. (D) Somatic mutation landscape. (E) Network Visual Analytics. (*p < 0.05, **p < 0.01, and ***p < 0.001).
NMF clustering and tumor microenvironment analysis
We employed the NMF algorithm to stratify KIRC patients into two distinct subgroups, designated as C1 and C2. Our analysis elucidated a robust consistency matrix for these patient clusters (Figure 2A). Additionally, we profiled the differential expression patterns of DRGs within these subgroups, as visualized in a comprehensive heat map (Figure 2B). Crucially, Kaplan-Meier survival analysis (Figure 2C) underscored a statistically significant disparity in patient survival rates between clusters C1 and C2 (p = 0.009), indicating potential prognostic implications. In terms of tumor microenvironment (TME) characterization, the C1 cluster demonstrated a marked enrichment of specific immune cell infiltrates, including CD8 T cells, monocytes, and neutrophils (Supplementary Figure S1A), suggesting a unique immunological landscape. Further, TME and stemness indices were analyzed in the C1 and C2 cohorts (Supplementary Figure S1B). To better understand the impact of clusters on immune cell infiltration in KIRC, we examined the differences in the tumor microenvironment between the C1 and C2 subtypes (Supplementary Figure S1C). These findings provide critical insights into the heterogeneity of the tumor microenvironment and its potential impact on the clinical trajectory of KIRC patients, as detailed in Supplementary Material S2.
FIGURE 2. NMF Clustering and survival analysis. (A) The NMF algorithm clusters patients into two subgroups. (B) The expression patterns associated with DRGs in the C1 and C2 clusters. (C) Difference in survival between patients in clusters C1 and C2. (*p < 0.05, **p < 0.01, and ***p < 0.001).
Hierarchical clustering of genes and module analysis
Our study adopts a hierarchical clustering approach to partition genes into distinct modules, and the resulting cluster dendrogram is illustrated in Figure 3A. To elucidate the relationship between these modules and key biological parameters, Figure 3B provides a graphical representation of the associations between the five identified modules and Stromal Score, Immune Score, as well as ESTIMATE Score. It can be found that MEblue and MEbrown modules are negatively correlated with Stromal Score, Immune Score, and ESTIMATE Score, while the other three modules are positively correlated with Stromal Score, Immune Score, as well as ESTIMATE Score. Among them, MEbrown has the closest relationship with Stromal Score. Furthermore, Figure 3C presents a comprehensive assessment of the linkages between individual modules and various phenotypic traits, shedding light on their potential functional relevance. It was found that the relationship between MEbrown module and phenotypic traits is the most significant, followed by the MEblue module as the second most significant. To gain further insights into the prognostic significance of the identified gene modules, univariate Cox analysis was conducted for genes within each module, and the results are visually represented in Figure 4A Subsequently, a risk scoring system was constructed through the application of LASSO-Cox regression, wherein six pivotal genes were selected from the most crucial module, as elucidated in Figures 3D, E. In order to provide a more in-depth understanding of the predictive capabilities of the model genes, Figure 4B showcases the outcomes of univariate Cox analysis, while Figure 4C elucidates the corresponding hazard ratio (HR) values and their associated p-values, thus offering valuable insights into the potential prognostic significance of these genes. Additionally, Figure 3F depicts the KM survival analysis results for individual genes, further illuminating their impact on patient outcomes. These analytical findings collectively contribute to a comprehensive assessment of the molecular landscape under investigation, offering critical insights into the potential clinical relevance of the identified gene modules and their constituent genes.
FIGURE 3. Hierarchical clustering methods for dividing genes into modules. (A) The cluster dendrogram. (B) Module-trait relationships. (C) Gene significance across modules. (D,E) Lasso model for screening model genes. (F) The KM survival analysis results for individual genes. (*p < 0.05, **p < 0.01, and ***p < 0.001).
FIGURE 4. Predictive power of modeled genes. (A) Univariate Cox analysis of genes within each model. (B) Univariate Cox analysis of model genes. (C) The corresponding hazard ratio values and their associated p-values. (*p < 0.05, **p < 0.01, and ***p < 0.001).
Risk stratification and validation of the risk model
Patients were categorized into low-risk and high-risk groups based on median risk scores. This stratification was pivotal in assessing the prognostic impact of the risk score. Figure 5A elucidates the outcomes of our univariate Cox regression analysis. It reveals a hazard ratio (HR) of 1.077 for the risk score [95% Confidence Interval (CI): 1.048–1.107, p < 0.001], signifying its significant prognostic value in the context of KIRC. Further, multivariate Cox regression analysis, detailed in Figure 5B, demonstrates an HR of 1.086 (95% CI: 1.050–1.122, p < 0.001) for the risk score. This analysis, accounting for various covariates, reinforces the risk score’s robustness as an independent prognostic indicator. To validate the prognostic model, TCGA samples were segregated into training and testing sets. Kaplan-Meier analysis, presented in Figure 5C, indicates a direct correlation between elevated risk scores and decreased overall survival (OS) and progression-free survival (PFS) in both sets, highlighting the predictive accuracy of our risk stratification approach. The calibration plot (Figure 5D) further underscores the clinical applicability of the risk score. It displays a remarkable concordance between predicted and observed survival outcomes, as evidenced by the calibration curves for 1-year, 3-year, and 5-year survival (Figure 5E). These curves exhibit a striking consistency between estimated and actual OS rates. ROC curves (Figure 5F) and decision curves (Figure 5G) further corroborate the risk score’s reliability. The ROC values (ROC train 0.845, test 0.738, all 0.788) demonstrate the model’s robust predictive capability. Lastly, the overall concordance index (c-index), depicted in Figure 5H, along with the c-index for the training (Figure 5I) and validation (Figure 5J) sets, confirms the model’s high predictive accuracy. These indices have reached satisfactory levels, reinforcing the model’s efficacy in forecasting clinical outcomes in KIRC patients. This comprehensive analysis solidifies the risk score as a vital tool in the prognostication of KIRC.
FIGURE 5. Risk Stratification and Validation of the Risk Model. (A) Univariate Cox regression analysis. (B) Multivariate Cox regression analysis. (C) Kaplan-Meier analysis. (D) The calibration plot. (E) Survival curves. (F) ROC curves. (G) Decision curves. (H) C-index. (I,J) Training and validation sets. (*p < 0.05, **p < 0.01, and ***p < 0.001).
Tumor microenvironment and immune infiltration analysis
We analyzed the TME and genetic profiles of patients with KIRC stratified into high-risk and low-risk groups based on their prognostic risk scores. Supplementary Figure S2A-C highlight that the high-risk group is characterized by a higher TME score, suggesting a more favorable prognosis for immunotherapy responsiveness. This observation underscores the potential for personalized therapeutic strategies based on risk stratification. In terms of the stemness index, the low-risk group exhibits an EREG-mRNAsi value closer to 1, indicative of a higher resemblance of tumor cells to stem cells. This finding may have significant implications for understanding the biology of lower-risk KIRC tumors. Supplementary Figure S2D further elucidates a strong association between the high-risk group and six types of immune infiltrating cells. This correlation may provide insights into the immunological underpinnings of high-risk KIRC. Tumor mutation burden (TMB) analysis shown no significant difference between high-risk and low-risk groups. However, Supplementary Figures S3B, C offer an in-depth view of the tumor gene mutation profiles for these groups, while Supplementary Figure S3D provides a comprehensive summary of the mutations, including variant classification and the top 10 mutated genes. Kaplan-Meier analysis (Supplementary Figure S3E) demonstrates that patients in the high H-TMB group exhibit a lower survival rate compared to those in the low H-TMB group (p < 0.001). Similarly, Supplementary Figure S3F shows that within the high-risk category, patients with high TMB have a lower survival rate than those in the low-risk group with low TMB (p < 0.01). Additionally, the differential expression of immune checkpoint-related genes is shown in Supplementary Figure S3J, with notable upregulation of TNF and HAVCR2 in the high-risk group. Supplementary Figure S3K presents a heatmap of the expression of these genes, while Supplementary Figure S3L highlights the differential expression of immune infiltration-related genes between the risk groups.
To further assess the value of the risk score for predicting immune checkpoint blockade (ICB) efficacy, The analysis includes PD1 and CTLA4 in the IPS analysis. However, no significant differences are observed in the average IPS values for PD1 and CTLA4 between low-risk and high-risk score groups, irrespective of their predicted response status. Finally, the high-risk group exhibits higher levels of immune dysregulation and immune exclusion reactions, providing critical insights into the immune landscape of high-risk KIRC (Supplementary Figure S4).
Drug sensitivity analysis based on risk scores
Spearman correlation techniques were employed to examine the association between risk scores and drug sensitivity within the context of KIRC. This analysis was conducted using data from the Genomics of Drug Sensitivity in Cancer (GDSC) database. Our findings, as illustrated in Figure 6, reveal that out of a panel of drugs evaluated, 40 exhibited significant correlations with the calculated risk scores. Intriguingly, the risk score demonstrated negative correlations with 20 drugs, including ULK1_ 4989, Topotecan, and Camptothecin. This negative correlation implies that higher risk scores in KIRC patients are associated with increased sensitivity to these drugs. Conversely, the risk score showed positive correlations with the sensitivities to 20 other drugs, such as OF-1, Sinularin, and Osimertinib. This positive correlation suggests that patients with higher risk scores are less sensitive to these therapeutic agents. Among these correlative drugs, ULK1_ 4989, and OF-1 emerged as having the most pronounced impact on drug sensitivity. This finding is particularly noteworthy as it suggests that these drugs could potentially serve as key therapeutic options, with their efficacy modulated by the patient’s risk score. This information could be pivotal in guiding personalized treatment strategies for patients with KIRC, tailoring drug choices based on individual risk profiles to optimize therapeutic outcomes.
FIGURE 6. The association between risk scores and drug sensitivity within the context of KIRC.
Pathway enrichment analysis and differential pathway expression
We employed GSEA to delineate the differential pathway activation between high-risk and low-risk patient groups. Figure 7A elucidates the GSEA results for the high-risk group, revealing significant enrichment in several critical biological pathways. Notably, this group demonstrated pronounced enrichment in pathways associated with immune response and autoimmune disorders, such as AUTOIMMUNE THYROID DISEASE, CYTOKINE-CYTOKINE RECEPTOR INTERACTION, PRIMARY IMMUNODEFICIENCY, SYSTEMIC LUPUS ERYTHEMATOSUS, and TYPE I DIABETES MELLITUS. These findings suggest a complex interplay between KIRC progression and immune dysregulation in high-risk patients. Conversely, the low-risk group, as depicted in Figure 7B, showed distinct pathway enrichment patterns. Key pathways such as NEUROACTIVE LIGAND-RECEPTOR INTERACTION, PEROXISOME, PPAR SIGNALING PATHWAY, RENIN-ANGIOTENSIN SYSTEM, and VALINE, LEUCINE, and ISOLEUCINE DEGRADATION were notably enriched. This enrichment indicates a divergent molecular profile in low-risk KIRC patients, potentially reflecting differences in metabolic processes and hormonal interactions. Further, our GSVA presented in Figures 7C, D, highlighted additional disparities between the groups. Noteworthy differences were observed in pathways such as PROXIMAL TUBULE BICARBONATE RECLAMATION and PHOTOPERIODISM. These results underscore the heterogeneity in biological processes and pathways between high-risk and low-risk KIRC patients, offering insights that could be pivotal in tailoring patient-specific therapeutic strategies and understanding the disease’s pathophysiology.
FIGURE 7. The differential pathway activation between high-risk and low-risk patient groups.(A) GSEA results for the high-risk group. (B) GSEA results for the low-risk group. (C,D) GSVA analysis.
Differential expression of model genes between normal and KIRC cells
A comprehensive analysis was conducted to examine the differential expression of specific model genes in tumor cells and clinical samples of KIRC. The results showed that CLIP4, HAS2, RNASE2, and KCTD12 genes were statistically significant in both tumor cells and clinical samples in Figure 8. Among them, RNASE2 was significantly upregulated (p < 0.01, Figure 8A), and the expression results of related mRNA showed a similar trend (p < 0.01, Figure 8B), which suggests a potential role for these genes in the tumorigenic processes or the adaptation of the tumor environment in KIRC. The enhanced expression of these genes might be reflective of alterations in cellular processes such as metabolism, signaling, or immune response within the tumor microenvironment. Conversely, a set of genes, namely CLIP4, HAS2, and KCTD12 (p < 0.01), were found to be notably downregulated in tumor cells. Besides, the relevant mRNA expression in clinical samples also revealed a significant downregulation (p < 0.05). The decreased expression of these genes could indicate their involvement in tumor suppressive pathways or mechanisms that are disrupted in the pathogenesis of KIRC. Crucially, the differential expression patterns of these model genes were found to be consistent with the results obtained from the gene differential expression analysis based on data from TCGA database. This consistency underscores the robustness of our model to a certain extent, reinforcing the reliability of our findings. The observed expression trends of these genes in tumor versus normal cells provide valuable insights into the molecular landscape of KIRC. These findings could be instrumental in identifying potential biomarkers for diagnosis or targets for therapeutic intervention, enhancing our understanding of the disease’s molecular underpinnings.
FIGURE 8. The expression level of LRP8, RNASE2, CLIP4, HAS2, SLC22A11, KCTD12 genes in cells and clinical samples. (A) Relative expression of 6 genes in normal human kidney cell lines (CCC-HEK-1) and KIRC kidney cancer cell lines (786-O cells) (B) Relative expression of 6 genes in clinical samples (*p < 0.05, **p < 0.01, and ***p < 0.001).
Discussion
KIRC is the most frequent histological subtype of renal cell carcinoma and is characterised by a relatively low prognosis and large metastatic potential. Understanding the underlying cellular and molecular mechanisms to guide biomarker discovery and intervention clinically is now a major challenge (Long et al., 2022). In recent years, as prognostic studies of KIRC have continued, the number of published articles has increased. Up to now, there have been only a few studies on disulfidptosis, which is a recently discovered new mode of cell death. Novel forms of cell death triggered by disulfide stress have important implications for cancer therapy. In contrast to common types of cell death, key genes involved in the regulation of disulfides include SLC7A11, SLC3A2, RPN1, and NCKAP1, as well as the Rac-WRC-Arp2/3 signaling pathway (Zheng et al., 2023). The investigation of the role of disulfidptosis in the prognosis of KIRC and the discovery of biomarkers is a therapeutic strategy worthy of investigation.
In this study, we performed a comprehensive analysis of disulfidptosis-related genes, tumor microenvironment, gene modules, risk stratification, immune infiltration, drug sensitivity, pathway enrichment, and differential gene expression in RCC based on the relevant literature research methods of Xu K et al. (Xu et al., 2023). Our findings reveal distinct chromosomal mappings and expression profiles of DRGs, notably the downregulation of NDUFA11 and upregulation of RPN1 in tumor tissues. The study further delves into the complexities of TME, highlighting differences in immune cell infiltrates and stemness indices between patient subgroups classified using the NMF algorithm based on the expression levels of disulfidptosis-related genes. Through hierarchical clustering and module analysis, we identified key genes and pathways associated with various biological scores and patient survival, underscoring their potential as prognostic markers. The development of a risk stratification model based on these findings provides a robust prognostic tool, as validated by Kaplan-Meier analysis, ROC curves, and Cox regression analysis using TCGA samples. Additionally, our exploration of drug sensitivities, utilizing Spearman correlation techniques, offers valuable insights for personalized treatment strategies by linking risk scores with drug responses. The study also employs GSEA to elucidate differential pathway activations between high-risk and low-risk patient groups, revealing a complex interplay between KIRC progression and immune dysregulation.
However, the study is not without limitations. The generalizability of the findings, especially regarding gene expression and drug sensitivity, may be constrained. The intricate nature of TME in KIRC warrants a more detailed analysis, and the conclusions drawn, particularly in subgroup analyses, could be limited by the cohort size. Furthermore, the correlation between risk scores and drug sensitivity does not necessarily imply causation and requires deeper investigation. Translating these findings into clinical practice faces challenges, including the need for validation in larger, more diverse cohorts and the integration of other clinical factors.
Despite these limitations, our study provides a detailed overview of the genomic and molecular landscape of KIRC, identifies potential biomarkers for prognosis and therapeutic targets, and underscores the significance of TME and immune infiltration in KIRC. What’s more, this study provides valuable insights for prognostic and therapeutic guidance by revealing the role of genes involved in the process of disulfuration in the prognosis and immunity of KIRC. However, our study is limited, and further experiments are needed to validate our findings and elucidate the underlying mechanisms. These findings lay the groundwork for developing personalized medicine approaches for KIRC treatment and enhance our understanding of its genetic underpinnings and complexity.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.
Ethics statement
The studies involving humans were approved by the Ethics Committee of the Shuguang Hospital Affiliated to Shanghai University of Traditional Chinese Medicine. The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.
Author contributions
YF: Conceptualization, Data curation, Investigation, Software, Writing–original draft. WW: Conceptualization, Data curation, Investigation, Methodology, Software, Writing–original draft. SJ: Formal Analysis, Methodology, Validation, Visualization, Writing–original draft. YL: Data curation, Formal Analysis, Methodology, Validation, Writing–original draft. YW: Project administration, Resources, Supervision, Writing–original draft, Writing–review and editing. XZ: Data curation, Methodology, Conceptualization, Formal analysis, Writing–review and editing. HZ: Funding acquisition, Project administration, Resources, Supervision, Writing–review and editing. GD: Funding acquisition, Project administration, Resources, Supervision, Writing–original draft, Writing–review and editing.
Funding
The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This work was funded by the National Natural Science Foundation of China (82174403, 82122075, 82074232, and 82030118), Shanghai Frontier Research Base of Disease and Syndrome Bi-ology of Inflammatory Cancer Trans-formation (2021KJ03-12), “Shu Guang” project supported by Shanghai Municipal Education Commission and Shanghai Education Development Foundation (21SG43), Three-year Plan Project of Shanghai Traditional Chinese Medicine (ZY (2021-2023)-0208).
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/fphar.2024.1343819/full#supplementary-material
References
Chu, J., Qian, X.-X., Zhang, X.-M., Jiang, T., Li, X.-J., and Sun, W. (2023). ETNK2 low-expression predicts poor prognosis in renal cell carcinoma with immunosuppressive tumor microenvironment. J. Oncol. 2023, 1743357. doi:10.1155/2023/1743357
Feng, L., Yang, K., Kuang, Q., Peng, M., Li, L., and Luo, P. (2022). A novel risk model for lncRNAs associated with oxidative stress predicts prognosis of bladder cancer. J. Oncol. 2022, 8408328. doi:10.1155/2022/8408328
Gao, Y., Cui, Y., Zhao, R., Chen, X., Zhang, J., Zhao, J., et al. (2022). Cryo-treatment enhances the embryogenicity of mature somatic embryos via the lncRNA–miRNA–mRNA network in white spruce. Int. J. Mol. Sci. 23, 1111. doi:10.3390/ijms23031111
Horr, C., and Buechler, S. A. (2021). Breast Cancer Consensus Subtypes: a system for subtyping breast cancer tumors based on gene expression. NPJ Breast Cancer 7, 136. doi:10.1038/s41523-021-00345-2
Hsieh, J. J., Purdue, M. P., Signoretti, S., Swanton, C., Albiges, L., Schmidinger, M., et al. (2017). Renal cell carcinoma. Nat. Rev. Dis. Prim. 3, 17009. doi:10.1038/nrdp.2017.9
Jia, Y., Li, H.-Y., Wang, J., Chen, X., Lou, L., Wei, Y.-Y., et al. (2023). Tripartite motif containing 69 elicits ERK2-dependent EYA4 turnover to impart pancreatic tumorigenesis. J. Cancer 14, 200–218. doi:10.7150/jca.79905
Kang, L., Wang, D., Shen, T., Liu, X., Dai, B., Zhou, D., et al. (2023). PDIA4 confers resistance to ferroptosis via induction of ATF4/SLC7A11 in renal cell carcinoma. Cell Death Dis. 14, 193. doi:10.1038/s41419-023-05719-x
Klümper, N., Ralser, D. J., Bawden, E. G., Landsberg, J., Zarbl, R., Kristiansen, G., et al. (2020). LAG3 (LAG-3, CD223) DNA methylation correlates with LAG3 expression by tumor and immune cells, immune cell infiltration, and overall survival in clear cell renal cell carcinoma. J. Immunother. Cancer 8, e000552. doi:10.1136/jitc-2020-000552
Li, J.-K., Chen, C., Liu, J.-Y., Shi, J.-Z., Liu, S.-P., Liu, B., et al. (2017). Long noncoding RNA MRCCAT1 promotes metastasis of clear cell renal cell carcinoma via inhibiting NPR3 and activating p38-MAPK signaling. Mol. Cancer 16, 111. doi:10.1186/s12943-017-0681-0
Li, Q. K., Pavlovich, C. P., Zhang, H., Kinsinger, C. R., and Chan, D. W. (2019). Challenges and opportunities in the proteomic characterization of clear cell renal cell carcinoma (ccRCC): a critical step towards the personalized care of renal cancers. Semin. Cancer Biol. 55, 8–15. doi:10.1016/j.semcancer.2018.06.004
Liu, X., Nie, L., Zhang, Y., Yan, Y., Wang, C., Colic, M., et al. (2023). Actin cytoskeleton vulnerability to disulfide stress mediates disulfidptosis. Nat. Cell Biol. 25, 404–414. doi:10.1038/s41556-023-01091-2
Liu, X., Olszewski, K., Zhang, Y., Lim, E. W., Shi, J., Zhang, X., et al. (2020). Cystine transporter regulation of pentose phosphate pathway dependency and disulfide stress exposes a targetable metabolic vulnerability in cancer. Nat. Cell Biol. 22, 476–486. doi:10.1038/s41556-020-0496-x
Liu, X., Zhang, Y., Zhuang, L., Olszewski, K., and Gan, B. (2021a). NADPH debt drives redox bankruptcy: SLC7A11/xCT-mediated cystine uptake as a double-edged sword in cellular redox regulation. Genes & Dis. 8, 731–745. doi:10.1016/j.gendis.2020.11.010
Liu, Y., Wu, Q., Fan, X., Li, W., Li, X., Zhu, H., et al. (2021b). A novel prognostic signature of immune-related lncRNA pairs in lung adenocarcinoma. Sci. Rep. 11, 16794. doi:10.1038/s41598-021-96236-4
Long, Z., Sun, C., Tang, M., Wang, Y., Ma, J., Yu, J., et al. (2022). Single-cell multiomics analysis reveals regulatory programs in clear cell renal cell carcinoma. Cell Discov. 8, 68. doi:10.1038/s41421-022-00415-0
Nomura, Y., Kakuta, E., Okada, A., Otsuka, R., Shimada, M., Tomizawa, Y., et al. (2020). Impact of the serum level of albumin and self-assessed chewing ability on mortality, qol, and ADLs for community-dwelling older adults at the age of 85: a 15 Year follow up study. Nutrients 12, 3315. doi:10.3390/nu12113315
Ohashi, R., Umezu, H., Sato, A., Abé, T., Kondo, S., Daigo, K., et al. (2020). Frequent germline and somatic single nucleotide variants in the promoter region of the ribosomal RNA gene in Japanese lung adenocarcinoma patients. Cells 9, 2409. doi:10.3390/cells9112409
Park, M.-K., Lim, J.-M., Jeong, J., Jang, Y., Lee, J.-W., Lee, J.-C., et al. (2022). Deep-learning algorithm and concomitant biomarker identification for NSCLC prediction using multi-omics data integration. Biomolecules 12, 1839. doi:10.3390/biom12121839
Sun, B., Luo, C., Zhang, X., Guo, M., Sun, M., Yu, H., et al. (2019). Probing the impact of sulfur/selenium/carbon linkages on prodrug nanoassemblies for cancer therapy. Nat. Commun. 10, 3211. doi:10.1038/s41467-019-11193-x
Tang, X., Ren, X., Huang, T., Miao, Y., Ha, W., Li, Z., et al. (2023). Prognostic and immunological significance of the molecular subtypes and risk signatures based on cuproptosis in hepatocellular carcinoma. Mediat. Inflamm. 2023, 3951940. doi:10.1155/2023/3951940
Wang, X., Lin, J., Li, Z., and Wang, M. (2023). In what area of biology has a “new” type of cell death been discovered? Biochim. Biophys. Acta Rev. Cancer 1878, 188955. doi:10.1016/j.bbcan.2023.188955
Wang, Z.-Z., Meng, T., Yang, M.-Y., Wang, W., Zhang, Y., Liu, Y., et al. (2022). ALYREF associated with immune infiltration is a prognostic biomarker in hepatocellular carcinoma. Transl. Oncol. 21, 101441. doi:10.1016/j.tranon.2022.101441
Xu, K., Zhang, Y., Yan, Z., Wang, Y., Li, Y., Qiu, Q., et al. (2023). Identification of disulfidptosis related subtypes, characterization of tumor microenvironment infiltration, and development of DRG prognostic prediction model in RCC, in which MSH3 is a key gene during disulfidptosis. Front. Immunol. 14, 1205250. doi:10.3389/fimmu.2023.1205250
Xu, Q.-S., Shen, Z.-Z., and Yuan, L.-Q. (2022). Identification and validation of a novel cuproptosis-related lncRNA signature for prognosis and immunotherapy of head and neck squamous cell carcinoma. Front. Cell Dev. Biol. 10, 968590. doi:10.3389/fcell.2022.968590
Xue, W., Sun, C., Yuan, H., Yang, X., Zhang, Q., Liao, Y., et al. (2022). Establishment and analysis of an individualized EMT-related gene signature for the prognosis of breast cancer in female patients. Dis. Markers 2022, 1289445. doi:10.1155/2022/1289445
Yan, Y., Teng, H., Hang, Q., Kondiparthi, L., Lei, G., Horbath, A., et al. (2023). SLC7A11 expression level dictates differential responses to oxidative stress in cancer cells. Nat. Commun. 14, 3673. doi:10.1038/s41467-023-39401-9
Zhang, L., Hobeika, C. S., Khabibullin, D., Yu, D., Filippakis, H., Alchoueiry, M., et al. (2022). Hypersensitivity to ferroptosis in chromophobe RCC is mediated by a glutathione metabolic dependency and cystine import via solute carrier family 7 member 11. Proc. Natl. Acad. Sci. U.S.A. 119, doi:10.1073/pnas.2122840119
Zhao, F., Wang, P., Guo, Y., Lu, Q., Kong, X., Su, D., et al. (2021). Identification of the potential roles of ring finger protein 8 in TP53-mutant breast cancer. Oncol. Lett. 21, 42. doi:10.3892/ol.2020.12303
Zheng, T., Liu, Q., Xing, F., Zeng, C., and Wang, W. (2023). Disulfidptosis: a new form of programmed cell death. J. Exp. Clin. Cancer Res. 42, 137. doi:10.1186/s13046-023-02712-2
Keywords: kidney renal clear cell carcinoma (KIRC), disulfidptosis, prognostic model, immune cell infiltration, tumor microenvironment, risk signature
Citation: Feng Y, Wang W, Jiang S, Liu Y, Wang Y, Zhan X, Zhu H and Du G (2024) Development of a prognostic model based on different disulfidptosis related genes typing for kidney renal clear cell carcinoma. Front. Pharmacol. 15:1343819. doi: 10.3389/fphar.2024.1343819
Received: 24 November 2023; Accepted: 21 February 2024;
Published: 13 March 2024.
Edited by:
Shuangiqan Yan, Fujian Normal University, ChinaReviewed by:
Wei Jie, Hainan Medical University, ChinaXiaoying Zhou, Guangxi Medical University, China
Copyright © 2024 Feng, Wang, Jiang, Liu, Wang, Zhan, Zhu and Du. 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: Guoqing Du, 530427540@qq.com; Huirong Zhu, stellargu0518@126.com
†These authors have contributed equally to this work