Identification of Disulfidptosis-Related Subtypes, Prognosis Model and Cerna Network in Breast Cancer

In females, breast cancer is the most common cancer with the highest incidence. In addition, breast cancer is the second leading cause of cancer death among women [1]. At present, the treatment of breast cancer includes surgery, chemotherapy, biologically targeted therapy, and other treatments. However, the resistance of some patients continuously increases after treatment, and the recurrence rate is also high, making maintaining the long-term survival of patients a significant challenge. According to a statistical study from 2022, the mortality pattern reflects the incidence trend and the mortality rate of breast cancer has declined more slowly in recent years, which indicates that the progress of breast cancer research has been stagnant [2]. The latest statistics from 2023 suggest that although the mortality rate of breast cancer patients is decreasing year by year, the increasing incidence rate may attenuate future progress [3]. Therefore, both the therapeutic efficacy and prognosis of BRCA are not ideal. New biomarkers for immunotherapy and the prognosis of BRCA are urgently needed to improve the treatment and prognosis of BRCA patients.

A recent study [4] has defined a new way of cell death modality, disulfidptosis, which is distinct from necrosis, besides cuproptosis, and may provide a new strategy for cancer therapy. This study demonstrates that tumor cells usually require the high expression of the SLC7A11 at post-transcription level, taking up more cystine to synthesize glutathione to balance the oxidative effects due to their highly active metabolic actions. But ingested cystine must be promptly reduced to cysteine to eliminate cysteine toxicity. Otherwise, the over-accumulated cystine may cause tumor cell death. They reported [5] that clinical studies found that metabolic therapy via glucose transporter inhibitors could trigger a disulfide fatty liver and inhibit cancer growth. In addition, many previous studies [6-8] have shown that disulfides play a significant regulatory role in immune cell recognition, tumor cell activity, and drug resistance in BRCA. Although there have been many previous studies exploring the mechanism of disulfide in breast cancer, it is worth noting that there have not been any studies on the mechanism of disulfidptosis in breast cancer.

Since breast cancer is a heterogeneous disease, tumor typing has crucial implications for achieving better clinical outcomes [9]. A study [10] identified nine breast cancer prognosis-related subtypes using the consensus clustering method based on DNA methylation patterns. Consensus clustering analysis identified three immune infiltration ways with different prognostic and biological features based on ICI patterns [11]. However, no studies have proposed a breast cancer-associated subtype under the disulfidptosis mechanism. LncRNA is a functional regulator that mediates multiple mechanisms and plays critical regulatory roles in various types of cancer growth [12]. For instance, it was found that lncRNA RP11-551L14.4 could inhibit the malignant development of BRCA by inhibiting the expression of miRNA miR-4472 [13]. A study also constructed a prognosis model through six genomic instabilityassociated lncRNAs, which could provide survival prediction for breast cancer patients [14]. Under the mechanism of disulfidptosis, the role of lncRNA in processes such as breast carcinogenesis and metastasis is still a gap. LncRNA TINCR weakens the efficacy of immunotherapy against BRCA via sponging DNMT1 and downregulating MiR-199a-5p through the STAT1-TINCR-USP20- PD-L1 axis [15]. The lncRNA BC069792 functions as a tumor suppressor gene by targeting KCNQ4 in BRCA, which can be done as a new diagnostic indicator and therapeutic target for BRCA [16]. However, the regulatory mechanism of the ceRNA network associated with disulfidptosis in BRCA remains unclear.

In summary, to explore the role of disulfidptosis in BRCA occurrence and development, three modular analyses were mainly performed in this study. First, we used the consensus clustering method to identify disulfidptosis-related subtypes in breast cancer based on three distinct patterns. Next, machine learning methods constructed a breast cancer prognosis model based on 10 DRLs. Finally, the DT algorithm was used to screen the DRCLs, the downstream miRNAs and mRNAs were screened through online databases, and a ceRNA regulatory network was constructed. It is incredibly novel that the contents of this study have never been reported before concerning each other. Our findings will provide academic help in understanding the mechanism of disulfidptosis in BRCA and offer new theoretical support for the prognosis and immunotherapy of BRCA.

Source and Preprocessing of Data

Two datasets, TCGA-BRCA and GSE42568, were used in this study. Among them, TCGA-BRCA datasets and clinical data were retrieved from the Cancer Genome Atlas (TCGA) database (https:// portal.gdc.cancer.gov/), totaling 1231 samples, containing 113 normal breast samples and 1118 tumor breast samples. After data in TCGA-BRCA were normalized to fragments per kilobase million (FPKM), the mRNA and lncRNA expression matrices were extracted. GSE42568 gene expression profiles were retrieved from Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm. nih.gov/geo/). This dataset has 121 breast samples, including 17 normal and 104 tumor samples. GSE42568 is from the GPL570 platform, and probe IDs for all genes were converted to gene names with reference to the platform file. After excluding normal samples of TCGA-BRCA and GSE42568, all tumor samples were merged via the “combat” method in the sva (v3.42.0) R package. The merged gene expression matrix normalized the data by the limma (v3.50.3) R package.

Acquisition of Disulfidptosis-Related Genes and Identification of Subtypes

According to the study [4], there are ten disulfidptosisrelated genes (DRGs), namely SLC7A11, SLC3A2, RPN1, NCKAP1, NUBPL, NDUFA11, LRPPRC, OXSM, NDUFS1 and GYS1. Using the Consensusclusterplus (v1.58.0) R package, consensus unsupervised clustering analysis was carried out, and classify patients into two subtypes according to the expression levels of DRGs. Correlations within groups increased after clustering, whereas correlations between groups decreased. Differences in the enrichment of KEGG pathways between the two subtypes were analyzed by GSVA and presented as heatmaps. In addition, ssGSEA analysis was performed to demonstrate which immune cells were differentially expressed in the two subtypes using boxplots. Then, principal component analysis (PCA) was performed according to the expression levels of DRGs. Finally, the difference in DRGs between the two subtypes was obtained by differential analysis, and the survival curves of the two subtypes were obtained by survival analysis.

Identification of Differential Genes in Different DRGs Subtypes, Consensus Clustering Analysis

First, we identified genes that were differentially expressed in the two subtypes of DRGs, and the genes with P < 0.05 and |logFC| > 1.5 were regarded as differentially expressed genes (DEGs) of the two subtypes of DRGs. Through consensus clustering analysis, patients were divided into two subtypes based on the expression levels of DEGs. Then, the differential expression of 10 DRGs on these two subtypes was exhibited by boxplots. Finally, GSVA, ssGSEA, and PCA were performed for these two subtypes of DEGs using the same analysis steps as in the previous section.

Analysis of Disulfidptosis-Related mRNA-lncRNA Co- Expression and Construction of A Prognosis Model

Through co-expression analysis, using P < 0.05 and correlation coefficient |R| > 0.3 as the selecting criterion, we then obtained the disulfidptosis-related lncRNAs (DRLs) co-expressed with DRGs. We subsequently merged the expression profiles of these DRLs and patients’ survival information for subsequent analyses. The DRLs associated with the prognosis of BRCA patients were obtained using univariate Cox regression analysis. After several trials, we randomly split these DRLs in a 4:1 ratio into training and test sets and adopted the LASSO regression algorithm to compress the features to avoid overfitting. Finally, multivariate Cox regression analysis was performed on the training set to construct a disulfidptosisrelated prognosis model based on the ten prognostically valuable DRLs. The risk score for DRLs was reckoned via this formula:

In the formula, DRL_Ci, DRL_Ei and N represent the regression coefficient of each DRL in the multivariate Cox regression analysis, the expression value of each DRL, and the total number of prognosisrelated DRLs, respectively. The median risk score was used to classify patients into high- and low-risk groups. We then evaluated the predictive value of the model by Kaplan-Meier survival curves and scatter plots. Univariate Cox and multivariate Cox regression analyses were used to identify whether these models and other Clinical traits were independent prognostic factors.

Drug sensitivity Analysis

The half maximal inhibitory concentration (IC50) of some anti- BRCA drugs was reckoned using the prophetic (V0.5) software package, set at P < 0.001 as the screening threshold, which excavates potential sensitive medicines for breast cancer patients. Differential expressions of drug sensitivity between high and low-risk groups was presented as boxplots by the Wilcoxon rank-sum test.

Immune Infiltration Analysis and Immunophenotyping

The relative proportions of 22 immune cells in each sample were calculated using the CIBERSORT algorithm. According to the results of immune infiltration analysis, the correlation between immune cells and risk score was explored, exhibited via scatter plots. The correlation between immune cells and 10 DRLs is presented as a heatmap. Each sample’s tumor microenvironment (TME) score was evaluated using the ESTIMATE algorithm. Then the expression differences in the stromal score, immune score, and estimate score between patients in high- and low-risk groups were explored. Obtain two immune subtypes via consensus clustering analysis. Differential analysis was then performed according to the expressed amount of DRGs and immune cells in the two subtypes, both visualized via boxplots. Finally, GSVA and ssGSEA were analyzed on these two immune subtypes.

Identification of Disulfidptosis-Related Characteristic lncRNAs

Decision tree (DT), a commonly used machine learning algorithm in data mining, fits high-dimensional data, can alleviate overfitting globally, and generate local optimal solutions for the model. Feature selection was performed on 178 previously obtained DRLs using the DT algorithm to select the disulfidptosis-related characteristic lncRNAs (DRCLs). We set overall > 3.5 as the screening threshold, where overall is the importance score. By taking the intersection of lncRNAs screened out by the DT algorithm and the above ten lncRNAs used to construct the prognosis model, we finally screened out 2 DRCLs.

Construction of ceRNA Regulatory Network

After identifying 2 DRCLs in the previous step, we performed the differential analysis of DRCLs based on TCGA-BRCA datasets. After confirming that DRCLs were differentially expressed between tumor and normal samples, miRNAs of targeting relationship with these 2 DRCLs were screened out via miRCode (http:// www.mircode.org/) and StarBase (https://starbase.sysu.edu. cn/) databases, respectively. Next, the predicted results from the two databases were intersected to obtain the miRNAs of targeting relationship with DRCLs. Before proceeding to the subsequent analysis, we verified the expression of these miRNAs via the dbDEMC database (https://www.biosino.org/dbDEMC/). According to the scientific hypothesis of ceRNA, key miRNAs that satisfy this hypothesis were identified. The target mRNAs of these key miRNAs were screened out via TargetScan (https://www. targetscan.org/) and miRDB (http://www.mirdb.org/) databases, respectively, and the intersection genes of the results of the two databases were made mRNA targets. The differential expression of these mRNAs was verified using the UALCAN (https://ualcan.path. uab.edu/) database, and mRNAs that met the scientific hypothesis of ceRNA were screened out. Finally, the ceRNA regulatory network was drawn by Cytoscape software.

Statistical Analysis

In this study, Pearson’s method was employed for correlation analysis. The log-rank test was used to assess the significance between Kaplan-Meier survival curves. All statistical analyses were performed with R (V 4.1.2) software. *** denoting p < 0.001, ** denoting p < 0.01 and * denoting p < 0.05. P < 0.05 was considered statistical significance.

Notably, datasets in two databases, TCGA and GEO, in which the expression data of TCGA-BRCA were separated into mRNA and lncRNA expression profiles, were used in this study. Then, the expression matrix of DRGs was extracted in the merged matrix of TCGA-BRCA mRNA and GSE42568 for subsequent disulfidptosisrelated molecular typing studies. Through TCGA mRNA-lncRNA co-expression analysis, disulfidptosis-related lncRNAs (DRLs) were identified in the TCGA-BRCA lncRNA expression matrix for subsequent construction of disulfidptosis-related prognosis model and ceRNA network. The flow chart of the analysis of this study is shown in (Figure 1).

Identification of Two Subtypes of DRGs

The network diagram (Figure 2A) demonstrates a strong correlation among the 10 DRGs and, in patients with breast cancer, all eight DRGs are risk factors, except NUBPL and RPN1 are favorable factors. All patients were clustered according to the expression profiles of the 10 DRGs by the consensus clustering method. The ideal number of subtypes was two, confirmed by the relative change in consensus matrix plots, cumulative distribution function (CDF) plots, and the area under the CDF curve (Figures 2B-2D). These two subtypes have the names A and B. Through GSVA analysis (Figure 2E), subtype A was markedly enriched in pathways like glutathione metabolism and arachidonic acid metabolism in the KEGG pathway. In contrast, subtype B was mainly enriched in pathways like ubiquitin-mediated proteolysis and inositol phosphate metabolism. Subtype A was enriched primarily on glutathione metabolism. This result suggested that subtype A might be a key subtype involved in the disulfidptosis process in the development of breast cancer.

The ssGSEA results (Figure 2F) revealed significant differences in the expression of most immune cells between subtypes A and B. Interestingly, among the immune cells with differential expression, with the exception of Type 2 T helper cell, which was downregulated in subtype B, other immune cells were significantly upregulated in subtype B. Results of PCA (Figure 2G) showed that different subtypes could be clearly distinguished via the expression levels of DRGs. We then assessed the differential expression of the 10 DRGs between subtypes A and B, and the results (Figure 2H) showed that five DRGs, including SLC7A11 and NCKAP1, were significantly downregulated in subtype B, whereas three DRGs, including SLC3A2 and RPN1, were significantly upregulated in subtype B. Finally, subtype B was found to have a higher probability of survival than subtype A by survival analysis, indicating that subtype B has a better prognosis (Figure 2I).

Consensus Clustering Analysis of DEGs in Different DRGs Subtypes

Of the two subtypes of DRGs derived above, we Determined 191 genes that were differentially expressed between the two subtypes via differential analysis. So, we next performed the consensus clustering analysis of these 191 DEGs according to the expression levels of DEGs, which finally resulted in two subtypes, A and B (Figures 3A-3C). We then explored the differential expression (Figure 3D) of the 10 DRGs between these two subtypes and found that SLC3A2, RPN1, and NDUFA11 were significantly upregulated in subtype A, whereas NCKAP1, NUBPL, and NDUFS1 were downregulated considerably in subtype A. Surprisingly, the results (Figure 3E) of GSVA showed that subtype A was highly enriched in the top 20 differentially enriched pathways relative to subtype B, which had multiple immune cell signaling pathways like T and B cell receptor signaling pathways. Looking at boxplots of the ssGSEA results (Figure 3F), we could see that all immune cells were significantly different in relative content between subtypes A and B, and more interestingly, only CD56dim natural killer cell and Neutrophil were less abundant in subtype A than in subtype B. All other immune cells were more abundant in subtype A than in subtype B. Finally, we evaluated the expression patterns of subtypes A and B using the PCA method, which showed good discrimination between subtypes A and B according to the expression levels of DEGs (Figure 3G).

Construction of a Prognosis Model Based on Drls

Through TCGA mRNA-lncRNA co-expression analysis, we obtained 178 DRLs (Figure 4A). We next merged the expression profiles of DRLs and patients’ survival information. After several trials that randomly divided these 178 DRLs into training and test sets, we found the optimal ratio of 4:1 for the training and test sets. Fourteen DRLs (Figure 4B) associated with the prognosis of BRCA were screened by univariate Cox regression analysis in the training set samples. We then used LASSO regression to compress the feature variables (Figures 4C&4D), and the best fit was achieved when lambda = 0.01634, and we obtained 13 DRLs with significant prognostic significance. Finally, the DRLs used to build the prognosis model were screened from the above 13 DRLs using multivariate Cox regression analysis, and as a result, 10 DRLs were selected, which were EXOC3-AS1, LINC00115, SNHG15, RASGRP3- Through TCGA mRNA-lncRNA co-expression analysis, we obtained 178 DRLs (Figure 4A). We next merged the expression profiles of DRLs and patients’ survival information. After several trials that randomly divided these 178 DRLs into training and test sets, we found the optimal ratio of 4:1 for the training and test sets. Fourteen DRLs (Figure 4B) associated with the prognosis of BRCA were screened by univariate Cox regression analysis in the training set samples. We then used LASSO regression to compress the feature variables (Figures 4C&4D), and the best fit was achieved when lambda = 0.01634, and we obtained 13 DRLs with significant prognostic significance. Finally, the DRLs used to build the prognosis model were screened from the above 13 DRLs using multivariate Cox regression analysis, and as a result, 10 DRLs were selected, which were EXOC3-AS1, LINC00115, SNHG15, RASGRP3-

Survival Analysis

For the prognosis model constructed in the previous section, the expression of 10 DRLs of all breast cancer patients was next substituted into the risk score formula to obtain a risk score for each patient. A median risk score was calculated based on each breast cancer patient’s risk score, and using the median risk score as a cut-off, all BRCA patients were split into high- and low-risk groups. The result (Figure 4F) illustrates the expression differences of 10 DRGs between high- and low-risk patients, among which SLC7A11, NCKAP1, LRPPRC, and NDUFS1 were significantly upregulated in the high-risk group. At the same time, SLC3A2, NDUFA11, and GYS1 were significantly downregulated in the high-risk group. Through survival analysis, the overall survival (OS) curve was plotted to evaluate the survival difference between patients in high- and lowrisk groups. In both the whole sample and the training set samples (Figures 5A&5B), the results revealed that patients in the highrisk group had a considerably poorer prognosis than those in the low-risk group. As the risk score increases gradually, breast cancer patients’ survival time and survival status change accordingly. And they showed that the higher the risk score, the lower the average patient survival time and the greater the danger of death (Figures 5C&5D). We then used the test set samples to test the robustness of this prognosis model, and the evaluation method was consistent with the above. The survival curves (Figure 5E) showed that the high-risk patients in the test set had a significantly worse prognosis than the low-risk patients. The risk curve and scatter plot (Figure 5F) demonstrate that patients have shorter survival times and a higher risk of death as the risk score increases. From this, the test set’s outcomes matched those of the training set, demonstrating that the prognosis model constructed with 10 DRLs could well predict the prognosis of patients with BRCA (Figure 5).

Assessing the Performance of Prognosis Model

To evaluate the performance of our constructed prognosis model, whether the risk score could be made as an independent prognostic factor for BRCA was first examined by univariate and multivariate Cox regression analysis. Fortunately, the univariate and multivariate Cox regression analysis results (Figures 6A&6B) were consistent. Both indicated that the risk score could be a prognostic factor independent of other clinical traits, which illustrated that the risk score could be an independent prognostic factor for BRCA. Next, we evaluated whether this prognosis model could balance the specificity and sensitivity of breast cancer prognosis by ROC curves and calculated the area under the ROC curve (AUC). The results indicated that this prognosis model had high accuracy, and in all samples (Figure 6C), the AUC for predicting 1-, 3-, and 5-year OS were 0.669, 0.707, and 0.686, respectively. The AUC results for the training and test sets are shown in (Figure 6D) and (Figure 6E), respectively. In addition, compared with the clinical traits, the CI value (Figure 6F) of the prognosis model was the highest, and its AUC value (Figure 6G) was also relatively high, thereby demonstrating the superior prognosis performance of this prognosis model. Finally, we constructed a nomogram by combining partial clinical traits and risk scores, and the overall score was obtained by summing all the scores for each patient’s clinical traits and risk and then substituting them into the nomogram could predict patients’ survival. As shown in (Figure 6H), the total score of this patient was 149, and the probabilities of survival time greater than 1-, 3-, and 5-year were 0.988, 0.935, and 0.874, respectively. The calibration curve (Figure 6I) was utilized to assess the nomogram’s capacity for prediction. The results indicated a high fit of the predicted and actual observed values, thus demonstrating that the nomogram has some predictive power (Figure 6).

Identification of Breast Cancer Associated Sensitivity Drugs

Through drug sensitivity analysis, we calculated the IC50 of some anti-BRCA drugs and explored the expression difference of drug sensitivity between patients in high- and low-risk groups of BRCA. The results (Figure 7A) indicated that BMS.708163, BX.795, GDC0941, OSI.906, and PLX4720 were significantly downregulated in patients in the high-risk group. Whereas 13 drugs (Figure 7B), including Vinblastine, Tipifarnib, and Salubrinal, were significantly upregulated in the high-risk group. The above results show five drugs, BMS.708163, BX.795, GDC0941, OSI.906, and PLX4720, with lower IC50 values in the high-risk group, might be potential therapeutic medications for high-risk patients. In contrast, the other 13 drugs could be appropriate for treating low-risk patients (Figure 7).

Correlation Analysis of Immune Cells and The Risk Score

After completing the immune infiltration analysis, we explored the correlation between immune cells and the risk score. A total of four immune cells (Figures 8A-8D), i.e., Macrophages M2, Macrophages M0, Mast cells activated, and Neutrophils, were significantly positively correlated with the risk score. In comparison, seven immune cells (Figures 8E-8K), including T cells follicular helper and Macrophages M1, were significantly negatively correlated with the risk score. In addition, we also explored the correlation between immune cells and 10 DRLs used to construct the prognosis model, and the correlation heatmaps (Figure 8L) showed that all immune cells had some correlation with 10 DRLs. T cells CD8 and 10 DRLs were significantly correlated, except for this immune cell’s significant negative correlation with AC108693.2 and this immune cell’s significant positive correlation with the other nine DRLs. Finally, we evaluated the TME score of every sample via the ESTIMATE method. We explored the expression differences of the stromal, immune, and estimate scores between patients in highand low-risk groups via difference analysis. Unfortunately, only the stromal score showed a significant expression difference between high- and low-risk groups, but there was no discernible difference in the immune and estimate score (Figure 8M).

Immunophenotyping

According to the outcomes of immune infiltration, we obtained two immune subtypes, A and B, using consensus clustering analysis (Figures 9A-9C). Through differential expression analysis (Figure 9D) of DRGs in these two subtypes, we found that SLC7A11, SLC3A2, and RPN1 were significantly upregulated in subtype A, whereas NUBPL, OXSM, and NDUFS1 were downregulated. We subsequently explored the differences in expression between the two subtypes and showed that four immune cells (Figure 9E), including NK cells resting and Macrophages M0, were significantly upregulated in subtype A, while 15 immune cells, including Plasma cells and T cells CD8, were significantly downregulated in subtype A. Furthermore, immune scores differed considerably between subtypes A and B, but stromal scores did not. Excitingly, the results of GSVA (Figure 9F) showed striking similarities with those of GSVA of DEGs subtypes, such as subtype A was also highly enriched in the top 20 pathways with enrichment differences compared with subtype B, which also had multiple immune cell signaling pathways like T and B cell receptor signaling pathway, indicating that this immune subtype A and the subtype A of DEGs were closely related.

The results of ssGSEA (Figure 9G) showed that all immune cells were significantly different in relative content between subtypes A and B and that only immune cell Neutrophil was more abundant in subtype B than in subtype A, and all other immune cells were less abundant in subtype B than in subtype A. As it happens, Neutrophil was also more abundant in subtype B than in subtype A in the ssGSEA results of the DEGs subtype. So subtype B of immune typing may closely correlate with the subtype B of DEGs typing. Finally, we assessed the risk score difference among the three subtyping modalities (Figures 9H-9J). And they have indicated that the risk score was significantly different only between the two subtypes of DRGs, but not among the DEGs subtypes and Immune subtypes. The patient typing, patient risk grouping, and patient survival status are presented in the Sanger diagram (Figure 9K).

Identification of DRCLs and Construction of a ceRNA Regulatory Network

We screened 38 DRCLs (Figure 10A) using the feature selection of 178 DRLs derived from the result of co-expression analysis via the DT algorithm. Then the intersection of these 38 DRCLs and 10 DRLs in which the prognosis model was constructed, we got 2 DRCLs, i.e., LINC00115 and SNHG15 (Figure 10B). Next, we validated these two DRCLs’ expressions in TCGA-BRCA datasets and found that LINC00115 and SNHG15 were more highly expressed in tumor samples than in normal samples in TCGA-BRCA (Figure 10C). Twenty miRNAs with a target relationship with LINC00115 and ten miRNAs with a target relationship with SNHG15 were screened out using the miRCode database. Six miRNAs with a target relationship with LINC00115 were screened out using the Starbase database, and 75 miRNAs with a target relationship with SNHG15 were screened out. Then integrating the screening results from the two databases, we obtained one miRNA target (hsa-miR-146b-5p) of LINC00115 (Figure 10D) and three miRNA targets (hsa-miR-4735- 3p, hsa-miR-338-3p, and hsa-miR-490-3p) of SNHG15 (Figure 10E).

To meet the regulatory relationship between lncRNAs and miRNAs in the ceRNA scientific hypothesis, we validated the expression of these miRNAs using the dbDEMC database and showed that hsa-miR-146b-5p [17] and hsa-miR-338-3p [18] were significantly under expressed in cancer samples (Figure 10F), whereas hsa-miR-4735-3p and hsa-miR-490-3p were overexpressed considerably in cancer samples. According to the ceRNA scientific hypothesis, we should screen out miRNAs that are lowly expressed in tumor samples as miRNA targets for DRCLs, and thus identify hsa-miR-146b-5p and hsa-miR-338-3p as key miRNAs. We screened 17 target genes of hsa-miR-146b-5p and 40 target genes of hsa-miR-338-3p by comprehension of the result of TargetScan and miRDB databases. Then, the differential expression of these target genes was validated using the UALCAN database. According to the hypothesis of ceRNA science, we should screen out mRNAs overexpressed in tumor samples as the key target genes of key miRNAs. Finally, we got seven key target genes of hsa-miR- 146b-5p (Figure 11A) and nine key target genes of hsa-miR-338-3p (Figure 11B). Finally, we got two lncRNAs-miRNA-mRNA regulatory networks that conformed to the hypothesis of ceRNA science. The target relationship of the ceRNA regulatory network was visualized by Cystoscope software (Figure 11C).

There are several challenges in BRCA prognosis and treatment. This challenge is caused by numerous factors, including the limited treatment options for some BRCA subtypes, the rise in tumor resistance, and the subpar performance of prognosis models, all resulting in unsatisfactory clinical outcomes [19]. Therefore, it is necessary to study new subtypes of BRCA, discover new sensitive drugs, and construct new prognosis models for improving the prognosis and immunotherapy of BRCA. Disulfidptosis is a new way of cell death. So far, to the best of the authors’ knowledge, no research has revealed the role of disulfidptosis in BRCA. Therefore, we discussed the biological characteristics of three patterns of disulfidptosis-related subtypes in breast cancer to improve the prognosis of BRCA. In addition, we also constructed a prognosis model based on 10 DRLs and a ceRNA regulatory network based on 2 DRCLs. First of all, we obtained the first disulfidptosis-related typing based on the expression of DRGs. Through GSVA, we found that glutathione metabolism was one of the activated pathways of subtype A, and glutathione metabolism was an essential substance in disulfidptosis, indicating that subtype A may have a close link with disulfidptosis.

A study [20] showed that buthionine sulfoximine inhibits glutathione production to specifically cause tumor regression in PI(3)K pathway mutant breast cancer cells. And study [21] into redox pathways in breast cancer and the startling discovery that glutathione peroxidase 2 is important in patient survival and oncogenic signaling. These two studies demonstrate a link between glutathione and breast cancer development. In the ssGSEA results, only Type 2 T helper cell has a higher abundance in subtype A than subtype B. Directly inhibiting spontaneous breast carcinogenesis, CD4+ T helper (Th2) cells cause the cancer cells to differentiate to their final state [22]. Therefore, under the mechanism of disulfidptosis, subtype A may play a significant role in the development of BRCA. We obtained the second disulfidptosisrelated typing according to the expression of DEGs. The result of GSVA showed that subtype A was highly enriched in multiple immune cell signaling pathways like T and B cell receptor signaling pathways.

A dominating gammadelta1 T cell population seen in lymphocytes infiltrating BRCA had a strong inhibitory influence on effector and dendritic cell development and function [23]. B lymphocytes that have infiltrated the tumor in cases of breast cancer group together and go through isotype switching, proliferation, and activation that are regulated by B-cell receptors [24]. It can be seen that subtype A may participate in activities such as breast cancer proliferation. After ssGSEA, we found a high abundance of CD56dim natural killer cell and Neutrophil in subtype B. Low levels of natural killer cell activity were seen in patients due to the targeted apoptosis of circulating CD56(dim) natural killer cells [25]. The primary factor contributing to breast cancer-related death is lung metastasis. The inflammatory microenvironment linked to neutrophils facilitates lung metastasis of tumor cells [26]. This suggests that subtype B may participate in the death process of BRCA cells. The third disulfidptosis-related typing was based on the immunization pattern.

Fortuitously, the immune subtypes shared strong similarities with the GSVA and ssGSEA results for the DEGs subtypes. Subtype A in both tapings was all significantly enriched in multiple immune cell signalings pathways like T and B cell receptor signaling pathways. In addition, Neutrophil in both typing’s was more abundant in subtype B. So the immune typing may have a close correlation with the DEGs typing. Then, we used univariate Cox, LASSO, and multivariate Cox regression analyses, which yielded 10 DRLs with predictive value, i.e., EXOC3-AS1, LINC00115, SNHG15, RASGRP3-AS1, AC108693.2, AC016065.1, AL157756.1, YTHDF3- AS1, AL139124.1, and AL121832.2. Study [27] identified seven m7G-related lncRNAs with predictive value, including LINC00115, and then built a BRCA prognosis model according to these seven lncRNAs. Through sponging miR-381, SNHG15 knockdown eliminated the cisplatin resistance of BRCA, offering a brand-new therapeutic target for the disease [28]. Apoptosis was increased in MCF7 cells, cell proliferation was reduced, and T-47D cells were made more vulnerable to the effects of tamoxifen and trastuzumab.

These effects were all caused by down-regulating RasGRP3 expression in BRCA [29]. Clinically, brain metastases in patients with BRCA are correlated with YTHDF3 overexpression [30]. Study [31] identified three m5C-related lncRNAs with predictive value, including AL121832.2, and a BRCA prognosis model was constructed based on these three lncRNAs. In previous reports, five of the 10 DRLs we screened have shown an association with breast cancer, while the other five have not been reported in any study. Of note, LINC00115 and AL121832.2 of the 10 DRLs were previously reported to be used to construct a prognosis model, which showed that our prognosis model built via the 10 DRLs had sure theoretical support. A risk score was reckoned for every patient based on the model formula, and the median risk score was used to classify patients into high- and low-risk groups. We next explored the expression differences of the 10 DRGs between highand low-risk groups. BRCA cells overexpress SLC7A11 to avoid cell death [32]. When NCKAP1 is silenced, the WASF3 complex becomes unstable, which reduces the ability of BRCA cells to invade [33]. Cell proliferation of SUM159 and MCF-7 cells was significantly enhanced by overexpression of NDUFS1 protein, promoting breast cancer progression [34].

These are consistent with our findings that SLC7A11, NCKAP1, and NDUFS1 are high-risk genes. The ROC curve could show the excellent predictive ability of the model. In all samples, the AUC for predicting 1-, 3-, and 5-year OS were 0.669, 0.707, and 0.686, respectively. In the training set sample, the AUC for predicting 1-, 3-, and 5-year OS were 0.654, 0.705, and 0.706, respectively. In the test set sample, the AUC for predicting 1-, 3-, and 5-year OS were 0.706, 0.719, and 0.619, respectively. It is worth mentioning that this is the first time a BRCA-related prognosis model has been constructed under a disulfidptosis mechanism. In addition, we screened for BRCA-related sensitivity drugs, in which BMS.708163, BX.795, GDC0941, OSI.906, and PLX4720 were determined as therapeutic drugs for high-risk breast cancer patients. Four inhibitors, such as OTSSP167 and BX.795, showed anti-neuroblastoma activity and anti-cancer properties [35]. Potently sensitized by GDC0941 to ABT-737 in BRCA, which increased ABT-737’s anti-cancer effectiveness [36]. Hormone-independent tumor development was more successfully controlled by the ER inhibitor, fulvestrant and OSI.906 when used together than when used separately [37].

However, BMS.708163 and PLX4720 have yet to be reported for the treatment of BRCA, indicating that they may be new drugs for treating BRCA. Finally, we identified 38 DRCLs using the DT algorithm and obtained 2 DRCLs by intersecting with the 10 DRLs used to construct the prognosis model. Two studies [38,39] reported that LINC00115 and SNHG15 were significantly up-regulated in BRCA. These are in line with the outcomes we validated in TCGABRCA. It was innovative to slow BRCA development by targeting the lncRNA SNHG15/miR-451/c-Myc signaling cascade [40]. Then we identified hsa-miR-146b-5p as the downstream miRNA target of LINC00115 and hsa-miR-338-3p as the downstream miRNA target of SNHG15. Both miRNAs were validated to be underexpressed in breast cancer through the dbDEMC database. We continued to screen the mRNA targets of these two key miRNAs, which ultimately resulted in seven targets of hsa-miR-146b-5p and nine targets of hsa-miR-338-3p. Therefore, we obtained two ceRNA regulatory networks that satisfied the scientific hypothesis of ceRNA. These two ceRNA networks have never been reported in previous studies.

This study is highly innovative in examining the role and biological characteristics of disulfidptosis in breast cancer. The present study has three significant contributions: First, three disulfidptosis-related typing based on three different patterns were obtained. Second, a disulfidptosis-related prognosis model was constructed based on 10 DRLs. Third, the DT algorithm was used to screen DRCLs, and two ceRNA regulatory networks were built, which conformed to the hypothesis of ceRNA science. Although the current study is based mainly on public databases, we supply new theoretical evidence for the mechanism of the treatment and prognosis of BRCA.

In summary, three different patterns of disulfidptosis-related typings were identified in BRCA, and the relationships between these typings and DRGs, KEGG pathways, as well as immune cells were explored. In addition, a disulfidptosis-related prognosis model was constructed, which exhibited excellent predictive performance. We then mined two new therapeutic agents for breast cancer. Finally, we built two disulfidptosis-related ceRNA regulatory networks. These have yet to be previously studied. Therefore, the findings of this study can provide a theoretical basis for improving the prognosis of breast cancer and provide new hope for reducing the mortality rate of breast cancer.

None.

No Conflict of Interest.