Radiosensitivity is associated with antitumor immunity in estrogen receptor-negative breast cancer

This study evaluated radiosensitivity and the tumor microenvironment (TME) to identify characteristics of breast cancer patients who would benefit most from radiation therapy. We analyzed 1903 records from the Molecular Taxonomy of Breast Cancer International Consortium cohort using the radiosensitivity index and gene expression deconvolution algorithms, CIBERSORT and xCell, that estimates the TME composition of tumor samples. In this study, patients were stratified according to TME and radiosensitivity. We performed integrative analyses of clinical and immuno-genomic data to characterize molecular features associated with radiosensitivity. Radiosensitivity was significantly associated with activation of antitumor immunity. In contrast, radioresistance was associated with a reactive stromal microenvironment. The immuno-genomic analysis revealed that estrogen receptor (ER) pathway activity was correlated with suppression of antitumor immunity. In ER-negative disease, the best prognosis was shown in the immune-high and radiosensitive group patients, and the lowest was in the immune-low and radioresistant group patients. In ER-positive disease, immune signature and radiosensitivity had no prognostic significance. Taken together, these results suggest that tumor radiosensitivity is associated with activation of antitumor immunity and a better prognosis, particularly in patients with ER-negative breast cancer.


Introduction
Immunotherapy has revolutionized cancer treatment paradigms since approval of immune checkpoint inhibitors (ICIs) for treatment of metastatic melanoma [1]. Traditionally, radiation therapy has been used to treat localized lesions via the direct cytotoxic effects of ionizing radiation. However, recent advances in ICI use reveal the importance of immunologic cell death and effects of radiotherapy on the tumor microenvironment (TME) [2,3]. ICIs have led to a paradigm shift in cancer therapy, but there are challenges associated with biomarker development. In particular, numerous researchers are working to predict treatment response when ICIs are combined with other treatments, such as chemotherapy or radiotherapy.
Radiotherapy is an integral component of loco-regional management of breast cancer. It may be used during almost every stage of breast cancer. However, it is not easy to evaluate radiation response in patients with tumors treated using a multimodal approach, because the aim of breast radiotherapy is usually adjuvant rather than radical. Different groups have developed genomic assays that predict radiation response in patients with breast cancer [4][5][6][7]. Although prospective validation is necessary, the results of retrospective studies suggest that radiotherapy-specific genomic assays can consistently identify patients who benefit from radiotherapy [8]. The radiosensitivity index (RSI) is a 10-gene rank-based signature developed to predict the radiation sensitivity of 48 cancer cell lines [9,10]. The prognostic value of the RSI has been independently validated for in multiple disease sites, including 1 3 breast cancer, head and neck cancer, non-small-cell lung cancer, pancreatic cancer, endometrial cancer, melanoma, and glioma [11,12]. However, results relating to the prognostic value of the RSI have been inconsistent between studies [11,13] that examined different molecular subtypes of breast cancer.
Breast cancer is a heterogeneous disease. Treatment response is highly variable across molecular subtypes. The use of ICIs has resulted in survival improvements in patients with triple-negative breast cancer [14,15], but has had minimal therapeutic effects in patients with estrogen receptor (ER)-positive breast cancer. Tumor-infiltrating lymphocytes (TILs) have critical roles in mediating responses to chemotherapy and improving clinical outcomes in breast cancer [16][17][18][19]. While TILs have prognostic value for ER-negative breast cancer, their roles in ER-positive breast cancer remain unclear [20,21]. Study findings suggest that radiosensitivity is associated with enhanced antitumor immunity in solid tumors [22][23][24]. However, the relationships between molecular subtype, TME, and radiosensitivity have not been defined for breast cancer.
This study was designed to evaluate radiosensitivity and TME according to breast cancer molecular subtype. It also examined the role of ER for antitumor immunity in relation to radiosensitivity.

Data source
We used data from the Molecular Taxonomy of Breast Cancer International Consortium Molecular Taxonomy of Breast Cancer International Consortium (METABRIC) study [25]. The clinical information and multi-omics data for 1,903 patients were downloaded from the cBioPortal for cancer genomics [26]. This information included gene expression microarray (Illumina Human v3 microarray; log-transformed mRNA z-scores) and clinical data including age, menopausal status, treatment information, tumor size, lymph node status, histology, PAM50 intrinsic subtype [27,28], and immunohistochemistry receptor status. The median follow-up time was 12.3 years for patients alive at censoring.

Transcriptomic analysis by TME and RSI
Radiation response was measured using the RSI, a rank-based linear algorithm [9][10][11]. Each of 10 RSI genes were ranked using the gene expression value (genes with the lowest and highest expression values were ranked as 1 and 10, respectively). The RSI value was determined using: Patients were categorized into radiosensitive (RSIlow) or radioresistant (RSI-high) tumor groups, using a previously defined cut-off point 0.3745 [22]. To classify patients according to the TME, compositions of 64 immune and stromal cells were estimated using the xCell package [29]. For xCell signatures, we performed k-means (k = 2) clustering using Cluster 3.0 and plotted heatmaps using Java Treeview [30]. We nominated two clusters as immune-high and immune-low groups. We also used the CIBERSORT deconvolution method to estimate immune cell proportions [31]. Differentially expressed genes (DEGs) were identified using eBayes in the R package, limma (version 3.46.0) [32]. KEGG pathway and GO enrichment analyses were performed using DAVID [33]. Gene Set Variance Analysis (GSVA) was applied to measure cancer hallmark pathway activity using the R package, 'GSVA' [34].

Statistical analysis
T-tests, Wilcoxon rank-sum tests, or chi-square tests were used to compare clinical characteristics according to radiosensitivity and immune subtype. Kaplan-Meier curves with log-rank tests and multivariate stepwise Cox regression analyses were used to estimate the prognostic roles of radiosensitivity and immune cell signatures. Statistical analyses were performed using R software (version 3.5.3; R Foundation for Statistical Computing, Vienna, Austria).

Radiosensitivity was associated with activation of antitumor immunity and a better prognosis
To examine RSI-associated pathways and the TME component, we performed Spearman's correlation analyses of RSI and enrichment scores for cancer hallmark pathway and 64 human cell types (Fig. 1a). RSI was negatively correlated with the antitumor immunity-related cell component and DNA damage repair-associated pathways ( Fig. 1a left), while stromal cells and metastatic properties were positively correlated with RSI ( Fig. 1a right). Using a pre-defined RSI cut-off value of 0.3745, each patient was classified into one of two tumor groups; radiosensitive or radioresistant. Results for baseline characteristics for each group are presented in Supplementary Table 1.
The radiosensitive subtype was more prevalent in tumors from older patients, and that had HER2 (human epidermal growth factor receptor 2)-negative, lymph node positive, and high cellularity results. We also screened molecular subtypes according to radiation response. There were more patients with claudin-low and luminal-B subtype tumors in the radiosensitive group, while there were more patients with basal, normal-like, and HER2-type tumors in the radioresistant group (Fig. 1b). Next, we used a pathway analysis to compare DEGs between the radiosensitive and radioresistant groups (Fig. 1c). The results indicated that patients with radiosensitive tumors had upregulated genes for oxidative phosphorylation, Myc target, and immune activation-related pathways. Patients with radioresistant tumors had up-regulation of genes associated with epithelial-mesenchymal transition (EMT), Wnt/ß-Catenin, PI3K-AKT-mTOR, and angiogenesis pathways. In the ER-negative group, the Kaplan-Meier curve and logrank test results revealed that patients with radiosensitive tumors had significantly better overall survival (log-rank hazard ratio (HR) 0.67, p = 0.001) and recurrence-free survival (log-rank HR 0.63, p < 0.001) (Fig. 1d). However, there was no significant survival difference between groups stratified by radiation response in patients with ER-positive tumors. A survival analysis of patients who received radiotherapy versus and those who did not found that the patients with radiosensitive tumors had significantly better outcomes if they had ER-negative tumors and received adjuvant radiotherapy ( Supplementary Fig. 1). There was no significant difference in survival for patients who did not receive radiotherapy. This result indicated that RSI scores predicted radiation response in patients with ER-negative tumors. In ER-negative group patients, radiosensitivity was associated with activation of antitumor immunity and a better prognosis. a Spearman's correlation between RSI and enrichment scores for cancer hallmark pathway, immune cells, and stromal cells. The top 10 positively correlated (right) and negatively correlated (left) cancer pathways or cell types, with RSI value. b Proportion of radiosensitiv-ity across PAM50 molecular subtypes. c Differentially expressed gene analysis between RS (radiosensitive; RSI-low) and RR (radioresistant; RSI-high) groups. Biologic processes enriched in the RS (upper) and RR (lower) groups. d Overall survival (upper) and recurrencefree survival (lower) stratified by RSI in the ER-positive (left) and ER-negative (right) group patients

Immune-high group had better prognosis than immune-low group in patients with ER-negative tumors
To reveal clinical characteristics associated with the TME, we defined two groups (immune-high and immune-low) using k-means clustering and an xCell deconvolution algorithm that predicted 64 immune and stromal cell types (Fig. 2a). When patients were clustered by immune and stromal cell profiles, there were 650 patients in the immunehigh group and 1,253 patients in the immune-low group. To verify whether immune subgroups reflected the TME, we performed a pathway analysis of DEGs between the two subgroups (Fig. 2b). We found that immune-activating pathways were upregulated in the immune-high subgroup. We also found that the cytoskeleton process and ER response were upregulated in the immune-low subgroup. Survival analysis revealed that immune infiltration status was a significant prognostic factor only in the ER-negative patients. For the ER-negative patients, the immune-high subgroup had better overall survival (log-rank HR 0.69, p = 0.002) and recurrence-free survival (log-rank HR 0.68, p = 0.003) than the immune-low subgroup (Fig. 2c). The immune-high patients showed significantly higher age, tumor size, number of lymph node, and histologic grade compared with the immune-low patients (Supplementary Table 2). We also found that radiosensitive patients were more prevalent in the immune-high subgroup. There were more luminal A/B patients in the immune-low subgroup and more claudin-low patients in the immune-high subgroup (Fig. 2d). The xCell Immune-high group had better prognosis than immune-low group in ER-negative group patients. a Immune and stromal cell enrichment scores are clustered into two groups using K-means clustering. The group in which patients had higher immunescore were classified as the immune-high group; otherwise, they were assigned to the immune-low group. b Differentially expressed gene analysis between immune-high and immune-low groups. Biological pro-cesses enriched in the immune-high (upper) and immune-low (lower) groups. c Overall survival (upper) and recurrence-free survival (lower) stratified by immune cell infiltration signatures in the ERpositive (left) and ER-negative (right) group patients. d Proportion of immune subgroups across PAM50 molecular subtypes. e Violin plot showing the distribution of immunescore according to PAM50 molecular subtypes immunescore was compared with the PAM50 molecular subtypes (Fig. 2e), and the claudin-low subtype patients had higher immunescore and luminal A/B patients had lower immunescore. Although the distributions of immunescore were significantly different across molecular subtypes, we found that the immuenscore was remarkably heterogeneous across patients, even within the same molecular subtypes.

ER pathway activity was associated with suppression of antitumor immunity
We found that ER-positive patients had no significant survival differences when patients were stratified by radiation response or immune profile (Figs. 1d, 2c). In order to define the underlying mechanism, we compared the clinical characteristics of ER-positive and ER-negative patients (Supplementary Table 3). The ER-negative patients were more prevalent among young, high grade, node positive, HER2-positive and immune-high patients. ER-negative patients were more frequently treated with chemotherapy and radiotherapy. We also measured the degree of immune cell infiltration status using the CIBERSORT algorithm; antitumor immune markers were significantly enriched in the ER-negative and radiosensitive groups compared to the ER-positive or radioresistant groups (Fig. 3a). On the other hand, immune suppressive markers, M2 or M0 macrophages, were more depleted in the ER-negative and radiosensitive groups, compared with the ER-positive or radioresistant groups. We speculated that ER pathway activity was associated with the immune suppressive environment. Antitumor immunity might not be observed in spite of immune cell infiltration in patients with ER-positive breast cancer. To examine our hypothesis, we identified DEGs between patients with ER-positive versus ER-negative tumors within the immune-high subgroup (Fig. 3b). A pathway analysis revealed that patients with ER-positive tumors had downregulated antitumor immunity-associated pathways (e.g., inflammatory response, interferon gamma response). The immune suppressive biomarker M2 macrophage was higher in patients with ER-positive tumors; PD-1 and PD-L1 expression, which represented immune activation was significantly lower in patients with ER-positive tumors (Fig. 3c upper; adjusted p < 0.001). Moreover, within the ER-positive group, the ESR1 expression level also represented immunologic status. Patients with higher ESR1 levels were significantly more likely to have the higher M2 macrophage and the lower PD-1 and PD-L1 expression indicated immune suppression (Fig. 3c lower; adjusted p < 0.001).

Radioresistance was associated with reactive stromal microenvironment
Unsupervised clustering of 650 immune-high group patients with cancer hallmark enrichment scores subdivided the patients into two groups. One group had an enrichment of EMT and angiogenesis-related pathways; this group was also significantly enriched with radioresistant group patients (Fig. 4a). Moreover, a gene set enrichment analysis revealed a similar trend toward the enrichment of these stromal and angiogenic gene sets in radioresistant tumors (Fig. 4b). We evaluated the degree of stromal cell infiltration status estimated from the xCell algorithm and found that stromal cells were significantly enriched in the radioresistant and immune-high group patients (Fig. 4c). The results for patients with ER-negative tumors stratified by radiosensitivity and immune status indicated that recurrence-free survival was highest in the immune-high and radiosensitive patients and lowest in the immune-low and radioresistant patients (Fig. 4d). Multivariable Cox regression analysis revealed that immune infiltration status and radiosensitivity were independent predictors of overall survival and recurrencefree survival ( Fig. 4e; Supplementary Table 3). The radiosensitive patients had a better prognosis within the immunehigh group patients.

Discussion
Immunotherapy has had promising results as a treatment for various types of solid tumors. However, its therapeutic role in breast cancer is limited to the treatment of triple-negative breast cancer [14,15]. Although breast cancer immunophenotypes have been actively studied, there is limited understanding of the radiotherapy-specific immune contexture in breast cancer. Using the METABRIC dataset, we found that radiosensitive tumors were characterized by enhanced antitumor immunity, and radioresistant tumors had reactive stromal microenvironments. Because ER pathway activity was associated with immune suppressive properties, radiosensitivity was significantly associated with a better prognosis, particularly in patients with ER-negative and radiotherapytreated breast cancer patients.
Our study identified that the immune status and radiosensitivity were prognostic only for ER-negative breast cancer in METABRIC dataset. Similar results were noted in the study by Torres-Roca et al., reported that RSI and RSI-based genomically adjusted radiation dose (GARD) are associated with local recurrence risk following radiotherapy, particularly in triple-negative patients [13,35]. In contrast, however, there was a report that RSI was an independent predictor of patient outcome in radiotherapy-treated ERpositive patients. Subgroup analysis was performed to define the mechanism via which an immune-based biomarker was prognostic only for patients with ER-negative tumors. In the METABRIC cohort, most patients with ER-positive tumors were not treated with chemotherapy, and about one-half of the patients with ER-negative tumors received chemotherapy. Therefore, we hypothesized that immune-based biomarkers were prognostic only in patients treated using chemotherapy. When patients were stratified by chemotherapy status, RSI and immune profile results were prognostic for patients with ER-negative tumors, regardless of chemotherapy status (Supplementary Fig. 2). Next, immunogenomic profiles were analyzed to examine the hypothesis that ER pathway activity suppressed antitumor immunity and that ER-positive tumors were unable to induce immunogenic cell death, in spite of TILs. We found that not only ER positivity, but also higher ER expressed, group patients had associated elevated M2 macrophage levels and immune suppression. Among TILs, macrophages are a major component Biologic processes enriched in the ER-(upper) and ER + (lower) group patients. c Distribution of M2 macrophage enrichment, PD-1 and PD-L1 expression, stratified by ER positivity (upper) or ESR1 expression(lower) level. *FDR < 0.05; ***FDR < 0.001; ER + estrogen receptor-positive, ER − estrogen receptor-negative, RS radiosensitive, RR radioresistant of the TME and orchestrate the immune response to tumor cell death [36]. Macrophages develop into one of two main subtypes, the classically activated macrophage (M1) or the alternatively activated macrophage (M2) [37]. M1 and M2 macrophages have different immune surveillance functions.
M2 macrophages secrete anti-inflammatory cytokines to induce immune tolerance, which have been studied as a prognostic biomarker [38,39]. Osteoporosis studies found that estrogen decreases M1 subtype activity and facilitates M2 macrophage polarization [40,41]. How estrogen induces Fig. 4 Radioresistance was associated stromal cells and metastatic properties. a Heatmap showing gene set variation analysis enrichment scores for immune-high patients. Unsupervised clustering identified that RR group patients were associated with angiogenesis and epithelial-mesenchymal transition-related pathways. b In the immune-high group, the top enriched gene sets in radioresistant (RSI-high) group patients are shown, compared with radiosensitive (RSI-low) group patients. c Distribution of endothelial cells, mesangial cells, fibroblasts, and overall enrichment of stromal cells (stro-mal score) estimated using xCell. Boxplots show stromal cell subset proportions according to ER status and radiation response. d Recurrence-free survival stratified by immune cell infiltration signatures and radiation response in ER-positive (left) and ER-negative (right) group patients. e Hazard ratio plot for recurrence-free survival (left) and overall survival (right) using multivariable Cox regression analysis. ***FDR < 0.001; ER + estrogen receptor-positive, ER − estrogen receptor-negative, RS radiosensitive, RR radioresistant alteration of cancer immunity remains unclear. Further studies are needed to define molecular mechanisms of estrogenmediated immune evasion.
Immunohistochemistry and transcriptomic analyses revealed three major tumor immunophenotypes, "immunedesert," "immune-excluded," and "immune-inflamed" [42,43]. Immune-inflamed tumors are characterized by high infiltration of cytotoxic lymphocytes that express PD-1 and PD-L1 expressing tumor cells. TILs are generally present in immune-excluded tumors, but cytotoxic lymphocytes are excluded from tumor cores and are caught in fibrotic tumor stroma. Immune-excluded tumors are characterized by a reactive stromal microenvironment with high TGF-ß signaling, increased presence of myeloid-derived suppressor cells, and tumor angiogenesis [44][45][46]. In this study, we analyzed the clinical characteristics by classifying the immune status into two groups (immune-high, immune-low) according to immune and stromal cell enrichment scores. Using cancer hallmark clustering analysis, we confirmed that stromal environment was important for radiation response (Fig. 4). The radioresistant and immune-high patients were identified as having an immune-excluded phenotype, while the radiosensitive and immune-high patients were characterized as having an immune-inflamed tumor phenotype. Dai et al. found that RSI value was associated with IFN-γ-mediated signaling and predicted therapeutic efficacy of PD-1 blockade in 12,832 primary tumors across 11 major cancer types in a merged microarray-acquired dataset and a The Cancer Genome Atlas (TCGA) database [24]. Using data from 10,240 genomically profiled distinct solid primary tumors from the Total Cancer Care database, Strom et al. found that low RSI values are significantly associated with immune activation [22]. The RSI was initially designed from cell line radiosensitivity data, which did not include TME information. Although further studies are needed, evidence that cellular biomarkers for radiosensitivity represent the TME suggests that intrinsic cancer cell phenotypes are critical for TME formation.
Breast cancer is a heterogeneous disease, and stratification by molecular subtype is critical to achieve better clinical outcomes. However, the results from single-cell RNA sequencing studies of the TME and intratumoral heterogeneity suggest it is not sufficient to assess treatment response based on molecular subtyping results [47]. The cellular analysis by Wu et al. revealed that epithelial, immune, and mesenchymal phenotypes existing within every tumor have remarkable heterogeneity [48]. In our study, we showed that distribution of immunescore varied across individual patients, within the same tumor molecular subtype (Fig. 2e). Therefore, the TME and treatment response will likely vary between individual patients.
This study had some limitations. First, because use of the RSI has not been validated in prospective clinical trials, it may not be justified to conclude that RSI-low group patients will be good responders to radiotherapy. However, we had additional findings that the RSI was well-correlated with the TME, which might represent a good response to radiotherapy. Second, recent evidence suggests that GARD (RSI adjusted radiation dose) is significantly associated with local recurrence and overall survival for patients with cancer treated with radiotherapy [12,35]. However, in the META-BRIC study, there are no data providing clinical variables including radiation dose, the type of nodal surgery, and lymphovascular invasion status, which may affect loco-regional recurrence. It would be an interesting study to determine whether GARD can predict loco-regional recurrence in large cohorts and how it differs depending on the molecular subtypes. Third, we found that ER pathway activity was associated with suppressed antitumor immunity, which was represented by M2 macrophage infiltration. Further studies are needed to elucidate the detailed mechanism for the role of this estrogen in immune regulation. Fourth, we did not differentiate the immune-desert subgroup in this study. Considering that the TME is an important factor affecting treatment response and survival in patients, future studies with more detailed immunologic subtyping may provide novel insight in the precision treatment of breast cancer. Finally, the METABRIC cohort was collected before 2012, treatment has evolved significantly last 20 years. In the METABRIC dataset, none of the HER2-positive patients received trastuzumab, and most of the ER-positive and/or node-negative patients did not receive chemotherapy. Since this study included the patients who received only radiotherapy without hormone therapy or chemotherapy, it can be a reasonable source of valuable information for responses to radiotherapy. However, RSI may not have the same prognostic role at this point as the systemic treatment has significantly changed. We performed external validation using TCGA breast cancer dataset and confirmed that the progression free interval was significantly different according to radiosensitivity in the patients who received radiotherapy ( Supplementary Fig. 4). However, radiosensitivity and immune infiltration status had the prognostic value for progression free interval in both ER-positive and ER-negative patients. Given that there is no long-term follow-up data for meaningful survival analyses for ER-positive TCGA patients, further studies would be necessary to predict the response of radiotherapy in the current era.

Conclusion
The results of this study suggest that tumor radiosensitivity is associated with activation of antitumor immunity and a better prognosis, particularly in patients with ER-negative breast cancer. By classifying tumor subsets with radiosensitivity and tumor-infiltrating immune and stromal cell signatures, we demonstrate a strong association among radiosensitivity, TME, ER pathway activity, and patient survival.