Identification and Validation of a Necroptosis-Related Prognostic Model in Tumor Recurrence and Tumor Immune Microenvironment in Breast Cancer Management

Background

Breast cancer (BRCA), as the most diagnosed cancer type and the leading cause of cancer-associated deaths in women, accounts for 701,000 deaths per year globally.¹ Due to the varied molecular subtypes, breast cancer is featured with high heterogeneity, hence the varied prognosis. Therefore, the importance of personalized medical therapies is emphasized in breast cancer management. In early breast cancer, surgery remains to be the first approach to treat the tumor, supplemented by chemotherapy or radiotherapy to prevent metastasis and recurrence. Though the prognosis has been markedly improved with the combination of chemotherapy, ER/PR-targeted or HER2-targeted therapy and immunotherapy, tumor relapse still poses a life-threatening event in breast cancer. Recurrent breast tumors are often unresponsive to treatments, hence the incurable outcome. As opposed to reactive medicine treatment after the onset of tumor relapse, the predictive and personalized paradigm focuses on the prediction and effective prevention before the clinical manifestation.² Yet, such predictive and personalized approaches towards breast cancer recurrence remain unsatisfactory. Signatures to guide effective approaches are under urgent demand.

Necrosis is proposed to occur when tumors outgrow the blood supply and is generally observed in aggressive solid tumors including breast cancer. It was regarded as an unregulated accidental cell death process before it was re-defined as a molecularly controlled regulated form of cell death with biochemical, genetic, and functional evidence.³ Necrosis is featured by rapid loss of plasma membrane integrity, organelle welling and mitochondrial dysfunction, while lacking typical apoptotic features such as internucleosomal DNA cleavage and nuclear condensation. Several cell-death modalities are included in regulated necrosis, such as necroptosis, parthanatos, ferroptosis, pyroptosis, pyronecrosis, mitochondrial permeability transition-dependent necrosis and NETosis.⁴ It has been illustrated that necroptosis of tumor cells under physio-pathological conditions is responsible for tumor necrosis and critical for metastasis in murine and human breast cancers.⁵

Necroptosis signaling is reported to be activated in breast cancer and correlated with its malignancy,⁶ and blockade of necroptosis signaling or pro-necroptotic proteins has been demonstrated to promote breast tumor recurrence.^6,7 Necroptosis, also known as programmed necroptosis, could be engaged by ligation of death receptors, which recruits RIPK1 to the intracellular region of the death receptor, followed by a series of ubiquitination and phosphorylation, resulting in necrotic cell death.⁸ It has been demonstrated that the interaction between RIPK1 and RIPK3 via the RHIM (RIP homotypic interaction motifs) domain is required for the initiation of necroptosis.⁹ Afterwards, necroptosis is executed by downstream components as reactive oxygen species (ROS), mitochondria, autophagy and so forth.¹⁰

Rising studies have been focusing on the role of necroptosis in immune response and immunotherapy.^11,12 It is implied that necroptosis might defend against tumor progression and enhance anti-tumor immunity via eliciting strong adaptive immune responses.¹³ In the process of necroptosis, cell-membrane ruptures as a result of necrotic cell death, and releases immunostimulatory intracellular components.¹⁴ In addition, primary macrophages are suggested to activate following necroptosis.¹⁴ Pro-necrotic signaling is also indicated to promote intratumoral immune response in breast cancer.¹⁵ Therefore, necroptosis might exert its positive impact on anti-tumor immunity through various mechanisms. In addition, studies concerning tumor immune microenvironment in breast cancer have unraveled that immune signature is associated with tumor recurrence in ER-negative breast tumors and enriched adaptive immunity portends satisfactory 5-year relapse-free survival (RFS).¹⁶

Despite the rising evidence emphasizing the correlation between necroptosis and tumor relapse, the detailed role of necroptosis in tumor immunity and relapse-free prognosis in breast cancer remains unclear. The potential value of necroptosis-targeted therapy as a predictive and personalized approach in breast cancer recurrence management is therefore undergoing exploration.

Therefore, we applied two independent prognosis models in the present study to investigate the potential value of necroptosis-related signatures in recurrence prognosis and tumor immune microenvironment in BRCA patient cohort. And we successfully identified a necroptosis-anchored tumor immune signature for patients with breast cancer, which helps to anticipate the relapse-free prognosis for patients and might enable a more tailored therapeutic regimen.

Methods

Data Processing

A total of 63 necroptosis-related genes were concluded from the literature.^17–20 Gene expression profiles and full clinical information of breast cancer patients were obtained from GEO database (GSE20685, GSE21653 and GSE25055). GSE21653 contained 266 early breast cancer patients who underwent initial surgery, gene expression data of 266 breast cancers were quantified by using whole-genome DNA microarrays (HG-U133 plus 2.0, Affymetrix), while 14 patients were excluded due to lack of important clinical information. GSE20685 contained 327 breast cancer samples, 268 patients underwent adjuvant chemotherapy and 91 patients had a relapse. GSE25055 contained 310 HER2-negative breast cancer cases treated with taxane-anthracycline chemotherapy pre-operatively and endocrine therapy if ER-positive. The patients with complete survival information were included in our analysis. Among them, GSE21653 was assigned as the training set, while GSE20685 and GSE25055 were enrolled for validation. Probes were re-assigned with gene symbols according to the corresponding gene platform. Robust Multi-array Average (RMA) method was applied to normalize raw data between datasets, including background correction, log2 transformation and normalization. RNA expression data was scaled with a standard deviation of 1 and a mean of 0. The RNA sequencing (RNA-seq) data (FPKM) and survival information of breast cancer were derived from The Cancer Genome Atlas (TCGA; https://tcga-data.nci.nih.gov/tcga/). Gene expression data from different samples is combined into genomicMatrix; all data is then log2 transformed.

Functional Enrichment Analysis

GO (gene ontology) analysis and KEGG (Kyoto Encyclopedia of Genes and Genomes) analysis were applied to perform the functional enrichment analysis, with the usage of R package clusterProfiler.²¹ Gene set variation analysis (GSVA) was also carried out with the R package “GSVA” to calculate the enrichment score of each pathway. The gene sets of “h.all.v7.2.symbols” downloaded in MSigDB and the known gene sets constructed by Mariathasan et al were used for GSVA enrichment analysis.

Construction of the Prognostic Gene Signature with Unsupervised Consensus Clustering

Sixty-three necroptosis-genes were selected from previously published literature. The patients were classified into two groups using the optimal k-means clustering. ConsensusClusterPlus R package was used to perform the cluster analysis with cycle computation 1000 times to ensure stability and reliability.²² The relapse-free survival (RFS) was calculated with the Kaplan–Meier method.

Construction of the Prognostic Gene Signature with LASSO-COX Regression

Lasso-Cox regression analysis was applied to identify the prognostic signature of genes of interest. Firstly, Lasso-penalized Cox analysis was performed with 10-fold cross-validation to narrow the genes for prognostic prediction. Six out of 63 necroptosis-related genes were screened out as a result. Thereafter, a stepwise multivariate Cox regression analysis was conducted to assess the contribution of a gene as an independent prognostic factor to relapse-free survival. The predictive model was established based on the expression levels of the 6 necroptotic genes and their relative coefficient (β) derived from the multivariate Cox regression analysis. The prognostic score = (−0.121*expression of ALDH2) + (−0.163*expression of CASP1) + (0.167*expression of EZH2) + (−0.242*expression of NLRP3) + (−0.365*expression of RIPK3) + (−0.431*expression of TLR3).

ESTIMATE Analysis

ESTIMATE (Estimation of Stromal and Immune Cells in Malignant Tumor tissues using Expression data) analysis²³ was performed to predict tumor purity, infiltrating stromal cells and immune cells in BRCA patients. And three scores, STROMAScore, IMMUNEScore and ESTIMATEScore, were yielded to evaluate the tumor microenvironment.

ROC Analysis

To investigate the capability of the established prognostic model to distinguish between subgroups, time-dependent receiver operating characteristic (ROC) analysis was performed with the “survivalROC” package in R software. Patients were divided into low-risk and high-risk groups based on the optimal threshold of the prognostic score calculated by „survminer‟ package in R. Further Kaplan–Meier survival analysis was used to assess the prognostic differences between groups. A two-sided Log rank test was performed afterwards, using „survival‟ package in R. Student’s t test was used to compare the difference between groups. All statistical analysis was performed using R software v4.0.1 and P < 0.05 was considered as the determining value of statistical significance.

Immunohistochemistry Analysis and Quantification

One hundred and seventy-one Breast cancer tissues embedded in paraffin were selected from Department of Pathology, Xiangya Hospital Central South University. The collection of tissue samples complied with the Declaration of Helsinki. The donors received the breast cancer resection during the span from 2010 to 2022 and were carefully followed-up. Deparaffination and rehydration of the tissue section were performed according to routine methods and followed by antigen retrieval with heat-induced epitope retrieval method. Endogenous peroxidase was blocked by incubation in methanol containing 3% hydrogen peroxide at room temperature. Blocking of nonspecific reactivity was performed with 10% goat serum in TBS. Primary antibodies were applied onto tissue sections overnight at 4°C. Second antibodies conjugated by horseradish peroxidase were applied afterwards. Color detection was performed by Diaminobenzidine and counterstaining was applied by Mayer’s hematoxylin. Dehydration and mounting were finally performed to complete the procedure. Permeabilization, if necessary, was applied by emerging sections in 0.1% TBS-Triton x100 buffer for 10 min.

For quantification, five pictures of random fields of view were taken for every slide with a Carl Zeiss microscope. Cells were counted by two independent researchers with the QuPath(v0.4.0). The primary antibody were ALDH2(Proteintech 15310-1-AP), CASP1(Proteintech 22915-1-AP), EZH2(Sigma, ZRB1095), NLRP3(Proteintech 19771-1-AP), RIPK3(Proteintech 17563-1-AP), TLR3(Invitrogen, PA5-81074).

siRNA and Cell Transfection

siRNA for genes (negative control, ALDH2, CASP1, EZH2, NLRP3, RIPK3, TLR3) were obtained from the commercial catalog in Ribobio, Guangzhou. siRNA transfections were performed using NEOFECT™ siRNA transfection reagent from Neofect (Beijing) biotech under instructions.

Cell Colony Formation Assay

Breast cancer cell line MD-MB-231 was purchased from Hunan Fenghui Biotechnology Co., Ltd., China. For the colony formation assay, equal amount of MDA-MB-231 cells were seeded into 6-well plate at a density of 500 cells per well. Cells were then cultured for 14 days before they were fixed and stained with crystal violet afterwards.

Combination Treatment with siRNA-Interference and Immunotherapy Therapeutic Study in vivo

All animal experiments were approved by the Ethics Committee at Xiangya Hospital Central South University. Mouse breast cancer cells E0771 (2x10^6 cells) were resuspended in 50 μL PBS (mixed with Matrigel at 1:1 ratio) and injected subcutaneously into the fourth pair of mammary fat pads of eight-week female c57BL/6 mouse. Tumor size was monitored every 4 days. siRNA of genes of interest and control siRNA (Ribobio, 5 nmol/kg) were dissolved in distilled water and were injected intratumorally every 3 days until the end time point of the experiment. Drugs of immunotherapy, anti-mouse TIM-3 (clone BT3-23, BioXCell, BE0115) plus anti-mouse PD-1 (clone RMP1-14, BioXCell, BE0146), were administrated with intraperitoneal injection at a dose of 10mg/kg once a week. The administration of siRNA and immunotherapy started 12 days after xenograft implantation.

Results

The Functional Pattern of Necroptosis-Related Genes in BRCA Cohort

In this study, gene expression and corresponding clinical data were achieved from GEO datasets, and a workflow of this study is shown in Figure 1A.

Figure 1 (A)Study flow chart. (B) GO enrichment analysis of the functional pattern of necroptosis-related genes in BRCA cohort. (C) KEGG enrichment analysis of differential necroptosis-related genes in BRCA cohort. (D) Forest plot analysis of individual gene of necroptosis-related gene set. (E) The interactions between the selected genes. This circle size represented the effect of each gene on the prognosis. Green dots stand for the favorable factors, while the purple dots stand for risk factors. Values were calculated by Log rank test.

Three GEO BRCA datasets, including 878 patients in total, with available survival data and clinicopathological information (GSE20685, GSE21653 and GSE25055) were enrolled in our study. GSE21635 with 252 patients was introduced as the training set (Supplementary Table S1). And GSE20685 and GSE25055 were assigned into the validation sets.

A total of 63-related genes were concluded from related literature. GO analysis was performed to validate the molecular functions of the genes of interest in BRCA patients. As expected, the selected genes were enriched in necroptotic processes and cell membrane structure. Intriguing, the synthesis and regulation of IL-1 (interleukin-1) were also abundantly enriched in this gene set, indicating the potential function of the selected necroptosis-related genes in immunity and inflammation (Figure 1B and C). Single impact of each gene in our selected gene set was evaluated with univariate Cox regression analysis and 12 of them were shown with significant impact on the patients’ prognosis in terms of tumor recurrence. Among them, ALDH2, CASP1, CASP9, IRF1, IRF2, NDRG2, NLRP3, RIPK3 TLR3 and TLR4 were shown to favor the prognosis, while EZH2 and PGAM5 presented negative role (Figure 1D). Further interaction analysis has also confirmed positive correlations between favorable factors as well as the negative correlation between risk and favorable factors (Figure 1E).

Consensus Clustering of BRCA Patients Based on Necroptotic Signature

Unsupervised consensus clustering was performed to examine the necroptotic pattern in breast cancer patients. Unsupervised consensus clustering serves as a functional method to proceed unsupervised class discovery, which works as a data mining technique to detect possible groups based on intrinsic biological characteristics. It provides quantitative and visual stability evidence to estimate the unsupervised subgroups in a dataset.²²

The consensus clustering matrix yielded 2 clusters with distinguishable differences. Further, principal component analysis (PCA analysis) confirmed a remarkable difference between the 2 clusters (Figure 2A and B and Table 1). It is worth mentioning that a worse relapse-free survival was indicated for the patients in cluster A compared to those in cluster B (Figure 2C). Necroptosis-related signature is therefore indicated to correlate with the tumor recurrence in breast cancer.

Table 1 Clinicopathologic Characteristics of Breast Cancer Patients According to the Necroptosis Pattern

Figure 2 Consensus clustering of BRCA patients according to necroptotic gene set. (A)Consensus clustering matrix for k=2. (B) PCA analysis of the 2 clusters was shown with low overlapped samples. (C) Kaplan-Meier analysis for the 2 clusters in BRCA. Patients in cluster A have worse relapse-free survival. (D) Unsupervised clustering of necroptosis-related genes in the BRCA cohort. BRCA molecular subtype, survival analysis, expression of KI67/P53/ERBB2/ER/PR, tumor Grade, primary node metastasis, primary tumor stage and histology were used as patient annotations. High expression was presented in red, while low in blue.

Unsupervised clustering of necroptosis-related genes has suggested that the molecular subtypes and PR/ER signatures have a clear differentiated distribution among the two clusters. Basal breast cancer samples with PR negative and ER negative were mainly distributed in cluster B, which was consistently featured by higher tumor grade. (Figure 2D, Supplementary Table S2)

Necroptosis Signatures Assign Breast Cancer with Differentiated Profiles of Tumor Immune Microenvironment

Gene set variation analysis (GSVA) analysis was performed to investigate the difference between necroptotic sub-clusters. We found that the immune-related signaling pathways were indicated to be the most altered pathways between the 2 clusters, suggesting the impact of necroptosis on immune system. Chemokine signaling was listed as one important pathway that has significant difference between clusters (Figure 3A). Consistently, patients included in cluster B were shown with significantly upregulated activity of the major immune cells (Figure 3B and C). ESITMATE analysis were performed between the two clusters and patients B cluster were shown to harbor significantly higher infiltration of immune cells (Figure 3D). Further analysis of immune checkpoint signature had additionally provided convincing data that patients in cluster B were rendered with higher sensitivity to immunotherapy with immune checkpoint blockade (Figure 3E).

Figure 3 Features of tumor immune microenvironment in two necroptotic clusters. (A) KEGG analysis of the differential genes between two clusters. (B and C) Immune cell infiltration analysis between two clusters. (D) ESTIMATE analysis for the comparisons between cluster A and B were presented in histogram. (E) Immune checkpoint analysis between two clusters.

Construction of Risk Prognostic Model Based on Necroptosis Signature

Lasso Cox regression analysis was performed on 63 necroptosis-related genes in GSE21653 cohort, to construct a necroptosis-related risk model. LASSO Cox regression serves as a method for variable selection and shrinkage in Cox’s proportional hazards model. A more refined model could be generated after a penalty function constructed by LASSO Cox regression analysis. It has been widely applied in screening genes in cancer research.²⁴

6 genes were filtered out by the Lasso analysis and were therefore selected to build up the prognostic model, ie TLR3, RIPK3, NLRP3, CASP1, ALDH2 and EZH2 (Figure 4A and B). In line with former Results, the 6 genes were also included in the 12 genes that were shown to influence tumor relapse in Figure 1. A subset of necroptotic gene set was therefore identified.

Figure 4 Risk model derived from necrotic-cell-death related genes. (A and B) Lasso Cox regression analysis of 63 necroptosis-related genes. (C and D) Kaplan-Meier Survival analysis for patients in high/low risk and ROC curve of 1-, 3-, 5-year survival in GSE21653. (E and F) Kaplan-Meier Survival analysis for patients in high/low risk and ROC curve of 1-, 3-, 5-year survival in GSE20685. (G and H) Kaplan-Meier Survival analysis for patients in high/low risk and ROC curve of 1-, 3-, 5-year survival in GSE25055.

Prognostic score was calculated based on the expression data of the 6 genes for each patient included. The score was further weighted by the multivariate Cox regression coefficients as following formula: (−0.121*expression of ALDH2) + (−0.163*expression of CASP1) + (0.167*expression of EZH2) + (−0.242*expression of NLRP3) + (−0.365*expression of RIPK3) + (−0.431*expression of TLR3). Afterwards, we obtained the optimum cut-off score using “survminer” package in R, patients had a score of 0.05 or higher were separated into high-risk group, and those with a risk score lower than 0.05 were divided into low-risk group (Table 2). A dramatically shortened relapse-free survival (RFS) was observed in the high-risk group (Figure 4C). Further subsequent time-dependent ROC analysis joined to confirm the high capability of such risk model to distinguish between the two subgroups. The separability of the necroptosis-related risk model was verified with the high AUCs of 0.759, 0.799, 0.791 at 1, 3, 5 year-RFS, respectively (Figure 4D).

Table 2 Clinicopathologic Characteristics of Breast Cancer Patients According to the Necroptosis-Related Signature

Validation of the Necroptosis-Based Prognostic Risk Model

GSE20685 and GSE25055 was applied into our study as the validation sets to assess the prognostic value of aforementioned necroptosis-based risk model. Patients in each validation set were stratified into low-risk and high-risk subgroups according to the necroptotic signature identified in the training set. In line, patients included in the high-risk group suffered a significantly shortened RFS compared to those in the low-risk subgroup (Figure 4E and G). Additional time-independent ROS analysis provided convincing evidence for the satisfactorily differentiable model of the necroptotic signature (Figure 4F and H).

Univariate and multivariate Cox analysis were performed in the training and validation sets. Our results indicated that the necroptotic signature acted as an independent prognostic factor in the training group (Supplementary Figure S1A, B). In line with the result in the training group, the stratification of patients into low- and high-risk subpopulations was validated as an independent prognostic factor in both validation sets. Additionally, the presence of lymph node metastasis and organ metastasis also play the roles as independent risk factor in the risk models in the validation database (Supplementary Figure S1C-F).

To further validate the predictive power of the Necroptosis-related signature for breast cancer patients, we tested the signature in TCGA dataset. According to the signature identified above, patients in the lower-risk group had significantly better OS and RFS In TCGA dataset, AUCs at 1, 3 and 5 years were 0.606, 0.594 and 0.645 for OS. AUCs of RFS at 1, 3 and 5 years were 0.580, 0.593 and 0.536 (Supplementary Figure S2).

To testify the value of our prediction model based on the necroptosis-signature for patients with breast cancer, in addition to univariate and multivariate Cox regressions to determine its accuracy in predicting clinical outcomes, we further compared our signature with other cell death-related prediction models (ferroptosis and pyroptosis).^25–30 As shown in Figure 5 and Supplementary Figure S3, our signature had better diagnostic efficiency than previously reported signatures.

Figure 5 Features of tumor immune microenvironment in low-risk and high-risk groups. (A) KEGG analysis of the differential genes between two subgroups. (B and C) Immune cell infiltration analysis between two subgroups. (D) ESTIMATE analysis for the comparisons between cluster A and B were presented in histogram. (E) Immune checkpoint analysis between two subgroups.

Furthermore, we investigated whether the risk value of our prognosis model was informative in treatment regiments with the dataset GSE25055. GSE25055 is a neoadjuvant study of 310 HER2-negative breast cancer cases treated with taxane-anthracycline chemotherapy pre-operatively and endocrine therapy if ER-positive. We observed that the risk scores were positively associated with the residual cancer burden (RCB) after neoadjuvant chemotherapy. Higher risk scores indicated poor pathological response to neoadjuvant chemotherapy (Supplementary Figure S4).

Clinicopathological Patterns and Tumor Immune Microenvironment Patterns in Necroptotic Risk Models

We further stratified patients according to clinicopathological risk factors and examine the prognostic value of necrotic cell death signature in the subpopulation of BRCA cohort with Kaplan-Meier analysis. Regardless of various clinicopathological factors, patients in the high-risk group suffered an unfavorable prognosis compared to those in the low-risk group (Supplementary Figure S5).

Based on the analysis of GSVA, chemokine signaling and some infection-related signaling pathways were among the important differential pathways between two risk groups (Figure 5A). Patients in the low-risk subgroup were shown with significantly upregulated infiltration of B cells, CD8+ T cells, eosinophils, macrophages, monocytes, NK cells, plasmacytoid dendritic cells and so forth (Figure 5B and C). An active immune microenvironment was therefore indicated. ESITMATE analysis also joined to confirm that patients in the low-risk group were shown to harbor significant higher infiltration of immune cells, however significantly elevated stroma infiltration was also suggested in the low-risk group (Figure 5D). Further analysis of immune checkpoint signature had revealed elevated expression of CD200, CD27, CD28, CD40, CD86, HAVCR2, NRP1, LAIR1, LGALS9, TNFRSF14, TNFSF14 in low-risk group, depicting the activated immune checkpoint patterns in the low-risk group (Figure 5E).

In addition, consistent immune profiles based on necroptosis-related signature were confirmed independent of prior therapies in the analysis of GSE25055 (Supplementary Figure S6). Though with different prior therapies, activated B cells, activated CD8+ T cells, eosinophils, immature B cells, MDSCS, macrophages, mast cell, monocytes, natural killer T cells, T follicular helper cells, type 1 helper cells stayed commonly more infiltrated in the low-risk group with better prognosis. ImmuneScore and ESTIMATEScore stayed higher in low-risk compared to high-risk group. Therefore, no significant changes of necroptosis-related immune signature were observed with prior therapies.

Individual Role of Each Gene from Necroptosis-Related Subset in Predicting Clinicopathological Features and Tumor Microenvironment

We identified a subset of necroptotic gene set from the Lasso-Cox regression model, composed by TLR3, RIPK3, NLRP3, CASP1, ALDH2 and EZH2. The expression of each gene in this subset was examined in patient sub-cohorts stratified by different clinicopathological risk factors. Tumor size, lymph node status, tumor grade, TRP53 mutation and molecular subtypes are included in the stratification. ALDH2 and TLR3 express differently in tumors of different sizes. ALDH2, EZH2, RIPK3 and TLR3 express differentially in diverse tumor grades. And only EZH2 and RIPK3 present different expression between patients harboring wild type TRP53 and those harboring mutated TRP53. Overall, the more advanced the tumor grade is and the bigger the tumor size is, the lower expression of ALDH2, CASP1, RIPK3 and TLR3 were detected as well as the higher expression of EZH2. While the expression of genes in subset stays unaltered in samples with different lymph node status, all of the 6 genes exhibit significantly different expression in molecular subtypes of breast cancer (Supplementary Figure S7).

Gene expression of the indicated six necroptotic subset genes was validated in both prognostic models. And it was indicated that the subset genes were correspondently elevated in the group with high immune cell infiltration in both models (Figure 6A and B). The Alluvial plot revealed the variation of necroptotic gene cluster, necroptotic risk group and tumor recurrence in individual BRCA patient and suggested the consistency of established prognostic models (Figure 6C).

Figure 6 (A) Gene expression of necroptotic subset genes in two clusters after consensus clustering. (B) Gene expression of necroptotic subset genes in high-/low- risk groups derived from risk models with Lasso-Cox analysis. (C) Alluvial diagram of BRCA patient necroptotic gene cluster, necroptotic risk group and tumor recurrence. (D-I) Prognostic validation of the six necroptosis-related genes, i.g. ALDH2, CASP1, EZH2, NLRP3, RIPK3 and TLR3. (J) ESTIMATE analysis of the selected 6 necroptotic genes and ESTIMATEScore, ImmuneScore and StromaScore was presented respectively. (K) The correlation between the selected 6 necroptotic genes and infiltrated immune cells of TME. (L) Immune checkpoint analysis based on the selected 6 necroptotic genes.

Prognostic validation was performed by Kaplan–Meier analysis for the indicated six necroptosis-related genes. Except for EZH2, the expression of which was negatively correlated with survival time, the rest 5 genes of necroptotic signatures were suggested to act as favorable prognostic factors (Figure 6D-I).

Consistently, TLR3, RIPK3, NLRP3, CASP1 and ALDH2 were exemplified to present positive association with ESTIMATEScore, ImmuneScore, StromaScore as well as immune cell infiltrations, whilst EZH2 was negatively correlated with the aforementioned scores (Figure 6J, K). Especially, EZH2 was shown to have negative correlation with infiltrated Eosinophils, MDSCs, macrophages, NK cells, Neutrophils and plasmacytoid dendritic cells. Correlation analysis between immune checkpoint signature and each gene in the selected 6 gene set had also joined to support the activated pattern of immune checkpoint upon the onset of the necroptotic subset (Figure 6L).

In addition, a small cohort of clinical BRCA specimens from our institute was included to join the verification of the selected necroptosis-related genes. A total of 171 BRCA samples with clinical information were collected from the department of pathology and breast oncology, Xiangya Hospital. Immunochemistry analysis was applied on specimen samples. With the help of the quantification of IHC staining, samples were classified into low- and high-expression groups for each protein of interest. Kaplan–Meier analysis was performed afterwards for each gene. In consistence, patients with higher expression of ALDH2, CASP1, NALRP3, RIPK3 and TLR3 have longer disease-free survival, while patients with high EZH2 expression suffered earlier tumor recurrence (Figure 7A).

Figure 7 (A) Prognosis of the patient in low- and high-necroptotic genes expression groups. (B) Efficiency of siRNA for necroptotic genes of interest. (C) Colony formation of MDA-MB-231 upon siRNA targeted necroptotic genes of interest. (D) Alteration of immune checkpoints after silencing of necroptotic genes. (E) Tumor growth of xenograft with the administration of co-blockade of anti-PD-1 and anti-TIM3 along with different siRNA injection. (ns, not significant; *, p<0.05; **, p<0.01; ***, p<0.001).

The Potential of Targeting Necroptosis in Immunotherapy in Breast Cancer Recurrence Management

As mentioned above, we have identified C10orf54, CD200, CD200R1, CD27, CD28, CD40, CD40LG, CD48, CD86, HAVCR2, LAIR1, LGALS9, TIGIT, TNFRSF14, TNFSF14 as the immune checkpoints uniformly altered in both stratification models. To validate the genetic changes, we silenced the key genes that we have selected out from the two stratification models in breast cancer cell line MD-MB-231 (ie, ALDH2, CASP1, EZH2, NLRP3, RIPK3, TLR3) and testified the impact on biological function and immune checkpoints (Figure 7B and D). We observed that the expressions of immune checkpoints of interest were significantly suppressed upon the silencing of CASP1 and NLRP3, while silencing of EZH2 led to a reverse result. Moreover, Cal-9 (LGALS9) were indicated to be consistently regulated along with the alteration of the necroptotic genes (Figure 7C).

Interestingly, HAVCR2(TIM-3) and its ligand LGALS9 were also observed to be significantly elevated in the low-risk group in both stratification models. As anti-TIM-3 therapy has been proved to facilitate the response to anti-PD-1 therapy due to its important role in T-cell exhaustion,^31,32 we further investigated if necroptotic signature could help to predict the outcome of breast cancer immunotherapy with co-blockade of TIM-3 and PD-1.

Therefore, the response to the co-blockade immunotherapy with anti-TIM-3 and anti-PD-1 under different manipulation of necroptotic genes was examined in vivo. As expected, treatment with siRNA of CASP1 and NLRP3 suppressed the response to co-blockade treatment, whereas siRNA-EZH2 significantly sensitized tumor response to the co-blockade immunotherapy (Figure 7E).

Discussion

Necroptosis, serving as an alternative form of programmed cell death, is the best-characterized form of regulated necrotic cell death.⁴ The role of necroptosis in tumor immunity remains controversial. On the one hand, it could evoke strong adaptive immune responses to defend against tumor promotion. On the other hand, the following recruited inflammatory response might also promote tumorigenesis and metastasis. An immunosuppressive tumor environment and promoted oncogenesis were also implied upon the necroptotic RIP1/RIP3 signaling via CXCL1.³³

We demonstrated that our selected gene set of necroptosis enriched in not only the programmed necrotic cell death but also in IL-1 synthesis and regulation pathway, suggesting the potential function of necroptosis in immunity and inflammation, which is consistent with previous reports.³⁴ Further, our current study applied two models to identify the subpopulations of BRCA cohorts. Unsupervised consensus clustering had resulted in two clusters with well-separable tumor immune microenvironment, in which patients in cluster B present much higher immune cell infiltration and dramatically enriched immunity and inflammation pathways, thus the better prognosis for patients in this cluster. Meanwhile, Lasso-Cox regression has assigned patients into low-risk and high-risk subgroup of breast cancer patients with explicitly-separable survival difference between groups. More immune cell infiltration and high immune score were concluded in patients in the low-risk group as expected, though the difference between the two risk subgroups was not as significant as that in two consensus clustered subgroups.

RIPK1 has been putatively recognized as one of the key necroptosis regulators. But it was also argued that ZBP1, instead of RIPK1, mediated tumor necroptosis in breast cancer.³⁵ Histologic tumor necrosis has been presumed to occur as tumor progresses aggressively. And it was suggested that a progressive switch from predominantly apoptotic to necrotic tumor cell death was involved in the development of invasive cancer. Thus the presence of necroptosis has been indicated to serve as an unfavorable prognostic factor for some patients.³⁶ However, the majority of key regulators of the necroptotic pathway are generally downregulated in tumor cells.¹³

In our current analysis, a subset of necroptotic gene set from the Lasso-Cox regression model was identified, composed of TLR3, RIPK3, NLRP3, CASP1, ALDH2 and EZH2. Among them, 5 necroptotic genes, TLR3, RIPK3, NLRP3, CASP1 and ALDH2, were suggested to be favorable prognostic factor, while EZH2 might predict poor outcomes for BRCA patients. To examine the detailed value of the selected necroptotic markers might help to make better therapeutic decisions.

Recent studies concerning genes involved in necroptosis have discovered its distinct role of RIPK3 in necroptosis.^37,38 Researchers have shown that necroptotic gene MLKL presents a positive correlation with B, NK and T lymphocytes in multiple breast cancer subtypes, and RIPK3 exhibits positive correlations with lymphoid cells only in HER2+ and TNBC breast cancer, which is consistent with our data. While DFNA5, a molecule mediating post-apoptotic secondary necrosis, is associated with the monocytic lineage and macrophages in ER+ breast cancer.¹⁵ However, RIPK3 was also reported to robustly increase in recurrent breast tumor cells and promote productive cell cycle.³⁹ It is also reported that necroptosis could be induced by the combination of activation of innate immune signaling through Toll-like receptor 3 (TLR3) and IFNγ treatment.¹⁴

Consensus clustering has provided us with separated clusters with drastically differential tumor immune microenvironments based on different necroptotic traits. BRCA patients in low-risk group, identified by Lasso-Cox regression model, exhibits consistent more infiltration of immune cell in tumor tissue. In addition, the selected subset of necroptosis-related genes was intriguingly shown to be significantly associated with immune cell infiltration in tumor, as well as the immune checkpoint profiles. Taken together, the necroptotic signature of breast cancer patient could be used as a potential biomarker to predict tumor immunity and patient’s response to immune checkpoint blockade immunotherapy.

Recent studies have revealed the important role of necroptosis in oncogenesis and tumor metastasis. Necroptosis-related lncRNAs were also indicated to correlate to the infiltration of activated CD4+ memory T cells, M1 macrophages, and resting mast cells, as well as some common immune checkpoint molecules in triple negative breast cancer.⁴⁰ Tumor necroptosis-mediated shedding of cell surface proteins was suggested to suppress anti-tumor immunity induced by T cells.⁴¹ Our study joined to provide the tumor immune microenvironment profile based on necroptosis signature.

In the current study, we have identified a few immune checkpoints that were consistently altered in both models, i.g. C10orf54, CD200, CD200R1, CD27, CD28, CD40, CD40LG, CD48, CD86, HAVCR2, LAIR1, LGALS9, TIGIT, TNFRSF14, TNFSF14. Among them, C10orf54 (also known as VISTA), Siglecs family (CD200, CD200R1), HAVCR2 (also known as TIM-3) as well as its ligand LGALS9, and NECL subset (TIGIT) were regarded as novel targets of immunotherapy, besides PD-1/PD-L1 and CTLA-4. In addition, besides the immune response-related pathways, JAK-STAT signaling pathway was indicated to be the commonly altered signaling pathway between two subgroups after stratification, inferring the potential regulation mechanism behind necroptosis and anti-tumor immunity microenvironment.

Further, necroptosis has been agreed as a novel target for cancer treatment, and growing discoveries of therapeutic compounds have been reported to defend against tumor by inducing necroptosis.⁴² Some compounds targeting RIP1 and RIP3 were verified to exert their effect of anti-tumor through the induction of necroptosis.⁴³ And the induction of necroptosis was also reported to circumvent cancer drug resistance in breast cancer.⁴⁴ Recent study focusing on necroptosis-related lncRNAs also implied that necroptosis-related lncRNAs signature could provide patients with a reference to make chemotherapy decision.⁴⁵ Our data also illustrated that low risk score in our necroptotic prognosis model could indicate a low RCB and thus a better pathological response to neoadjuvant chemotherapy. A hint of personalized neoadjuvant chemotherapy was therefore indicated for patients with low risk.

In addition, we validated that regulation of some key necroptotic genes could affect the response to immunotherapy with combined inhibitors of TIM-3 and PD-1. Patients in the necroptotic low-risk group were thereafter indicated with better response to co-blockade of TIM-3/PD-1 immunotherapy. The potential of targeting necroptosis as a novel therapy for malignancies is therefore fundamentally implied.

By analyzing necroptosis signatures in breast cancer, it was implied that the necroptotic signatures could have predictive value in terms of cancer relapse-related prognosis as well as immunotherapy response. Based on the molecular profile of necroptosis signature, we could identify patients who were more likely to present early cancer relapse or respond poorly to immunotherapy. These could further help in stratifying patients into different risk groups and guiding treatment decisions. In our study, highly activated necroptotic signatures indicated late recurrence occurrence. In addition, more immunogenic malignancy with activated immune microenvironment was also implied in the activated necroptotic signature, suggesting the more efficiency of immunotherapy in this scenario. Specifically, distinctly altered immune checkpoints upon varied necroptotic signatures were also inferred in our study, which could aid in the prediction of response to immunotherapy with certain immune checkpoint blockades.

Further, the necroptotic signature is suggested to contribute to the personalization of medical services in recurrent breast cancer management. By integrating information about the necroptotic profile with other clinical and molecular data, healthcare providers could develop individualized treatment plans for patients. Taken the results from our study, personalized strategies could include selecting therapies that target necroptosis gene subset based on the patient’s molecular profile. Patients with highly activated necroptotic signatures and corresponding activated immune microenvironments could be assigned with immunotherapy like immune checkpoint blockades. Whereas patients, with suppressed necroptotic signatures should be tailored with specific necroptosis-targeting interventions to aid in immunotherapy. Personalization of medical services can lead to improved treatment outcomes, reduced adverse effects, and optimized healthcare resource utilization.

Conclusion

We demonstrated two strategies to stratify breast patients based on their necroptotic profiles and showed that necroptosis-related signature could assign patients with distinct tumor immune microenvironment patterns and corresponding relapse-free prognosis. A subset of necroptotic gene set, composed of TLR3, RIPK3, NLRP3, CASP1, ALDH2 and EZH2, was identified and could be used as a biomarker set for predicting tumor recurrence and immunotherapy response in breast cancer, and might help to tailor personalized approach in breast cancer recurrence management.