Tumor Microenvironment Characterization in Breast Cancer Identifies Prognostic and Neoadjuvant Chemotherapy Relevant Signatures

Immune response which involves distinct immune cells is associated with prognosis of breast cancer. Nonetheless, less study have determined the associations of different types of immune cells with patient survival and treatment response. In this study, A total of 1,502 estrogen receptor(ER)-negative breast cancers from public databases were used to infer the proportions of 22 subsets of immune cells. Another 320 ER-negative breast cancer patients from Guangdong Provincial People’s Hospital were also included and divided into the testing and validation cohorts. CD8+ T cells, CD4+ T cells, B cells, and M1 macrophages were associated with favourable outcome (all p <0.01), whereas Treg cells were strongly associated with poor outcome (p = 0.005). Using the LASSO model, we classified patients into the stromal immunotype A and B subgroups according to immunoscores. The 10 years OS and DFS rates were significantly higher in the immunotype A subgroup than immunotype B subgroup. Stromal immunotype was identified as an independent prognostic indicator in multivariate analysis in all cohorts and was also related to pathological complete response(pCR) after neoadjuvant chemotherapy. The nomogram that integrated the immunotype and clinicopathologic features showed good predictive accuracy for pCR and discriminatory power. The stromal immunotype A subgroup had higher expression levels of immune checkpoint molecules (PD-L1, PD-1, and CTLA-4) and cytokines (IL-2, INF-γ, and TGF-β). In addition, patients with immunotype A and B diseases had distinct mutation signatures. Therefore, The stromal immunotypes could predict survival and responses of ER-negative breast cancer patients to neoadjuvant chemotherapy.


INTRODUCTION
Breast cancer is the most common cancer in women and a leading cancer worldwide. Genomic changes in cancer cells may be used to predict prognosis and treatment responses as well as to develop new targeted therapy (Cancer Genome Altas Network, 2012;Curtis et al., 2012;Ali et al., 2014;Ji et al., 2019). Recently, it has been reported that the tumor microenvironment also played an important role in tumor progression and chemotherapy efficacy. breast cancer is composed of a mixtures of cancer cells and non-cancer cells such as stromal cells, vascular cells, and tumor-infiltrating lymphocytes (TILs), with the roles of noncancer cells remain unclear. Among the non-cancer cells, TIL values have been reported to be associated with pathologically complete response (pCR) and prolonged overall survival (OS) in patients with breast cancer (Ladoire et al., 2008;Kitano et al., 2017;Luen et al., 2017), it also could be used as a drug target, but the roles of specific immune cells have not been well clarified. To better understand the diverse immune cells of breast cancer to the response of treatment and construct the classification of stromal immunotype for survival prediction, we enumerate immune cells in a way that accounts for the breadth of their specialized functions, and to reliably investigate the interaction of the immune response with breast cancer, and finally find the effective parameters to support the personalized therapy.
Therefore, the aim of this study was to quantify the composition of immune cells and investigate their relationships with responses to neoadjuvant chemotherapy and survival of breast cancer patients.

Study Population
The present study was approved by the Ethics Committee of the Guangdong Provincial People's Hospital. Three independent cohorts of patients with breast cancer were included. The training cohort was comprised of 1,502 estrogen receptor(ER)negative breast cancer patients with gene expression data from the public studies, details of patients from public studies can be found in corresponding publications (Supplementary Table S1), some data were compiled and created as previously described (Haibe-Kains et al., 2012;Ali et al., 2016). some data were downloaded from Gene Expression Omnibus and ArrayExpress(TCGA and METABRIC). After excluding patients with comorbidities (e.g., other malignant tumors), incomplete follow-up data, and metastatic disease, 320 patients with pathologically diagnosed ER-negative breast cancer using needle core biopsy at the Guangdong Provincial People's Hospital between June 2009 and December 2015 were selected and randomly divided into the testing cohort (n 218) and the validation cohort (n 102) by using computer-generated random numbers. pCR was defined as the absence of any residual invasive carcinoma or DCIS on pathologic review of a surgical specimen following neoadjuvant chemotherapy. The study design is shown in Supplementary Figure S1. Informed consent was signed by each patient to allow the use of their data in clinical researches.

Immunoscore Establishment
Normalized gene expression data were used to infer the relative proportions of 22 types of infiltrating immune cells using the CIBERSORT algorithm. For the TCGA dataset, RNA sequencing data were transformed using voom, converting count data to values more similar to those resulting from microarrays (Cancer Genome Altas Network, 2012;Law et al., 2014). The CIBERSORT is a welldesigned method for the analysis of microarray gene expression profiles (Ali et al., 2016;Li et al., 2019). The LM22 file covers 547 genes that can be used to accurately discriminate 22 distinct functional subsets of immune cells and activation states including seven T cell types, naive and memory B cells, plasma cells, NK cells, and myeloid subsets. So in our study, we used the CIBERSORT and the LM22 gene signature file to estimate the fractions of diverse immune cells in the patients. The associations of immune cells with the OS and diseasefree survival (DFS) of patients in the training cohort were analyzed using the estimated fractions by univariate analysis. The COX regression analysis with the least absolute shrinkage and selection operator (LASSO) model was used to develop an immunoscore for the classification of immunotypes A and B breast cancer (Figures 1A-C).

Statistical Analysis
The coefficients and partial likelihood deviance for each prognostic feature were calculated with the "glmnet" package in the R program. The LASSO Cox regression model was used to estimate the ideal coefficient and likelihood deviance. Restricted cubic spline regression was used to characterize the relationship between immune score and patient survival in all three cohorts. Besides, The nomogram integrating immune type and clinical parameters was constructed according to results of the

Stromal Immunoscore and Stromal Immunotype
The fractions of the 22 types of immune cells are showed in Figure 1. Among them, five types (CD8 + T cells, CD4 + T cells, B cells, M1 macrophages, and Treg cells) were significantly associated with OS and DFS (Supplementary Figure S2), so we chose these five types as research cells. In the testing and validation cohorts, the mean densities of B cells, CD8 + T cells, CD4 + T cells, Treg cells, and M1 macrophages in the stroma were higher than the densities in tumor core (38.320, 77.234, 45.291, 26.291, and 71.248 cells/mm 2 Vs 9.580, 25.745, 10.065,6.5728, and 24.568 cells/mm 2 , respectively). (Figure 2). The immune cell infiltration in the stroma was significantly associated with the OS of ER-negative breast cancer (Supplementary Table S2). The representative images of these immune cells are shown in Figure 2. We built prognostic classifiers by using the LASSO COX model in the training cohort (Supplementary Table S1; Figure 1) and calculated the coefficients using the following formulas: IS 2.202 × Treg cell count-3.455 × B cell count-2.559 × CD8 + T cell count-2.077 × CD4 + T cell count-2.074 × M1 macrophage count. Therefore, We classified patients into immunotype A and immunotype B groups based on the median immunoscore. Clinical characteristics of patients with breast cancer of stromal immunotypes A and B in the training cohort and the testing and validation cohorts did not vary significantly ( Table 1).

Association of Stromal Immunotype With Patient Survival
As shown in Figure 3, tumors with high immunoscore generally contained increased Treg cells and reduced B cells, CD8 + T cells,     (Figure 4).

Association of Stromal Immunotype With Response to Neoadjuvant Chemotherapy
To quantitatively predict responses to neoadjuvant chemotherapy, we constructed a nomogram which integrated both stromal immunotype and clinicopathological factors using the data from the testing cohort and validated it using the data from the validation cohort. First, we analyzed the relationships between pCR and clinicopathologic characteristics and found that the rate of pCR was associated with tumor size, histologic grade, HER2 status, TNM stage, immunotype in the testing cohort and the validation cohort (all p < 0.05) ( Table 2). Multivariate logistic analyses in both cohorts showed that histologic grade, TNM stage, HER2 status, and immunotype were independent predictors for pCR (all p < 0.05) ( Table 3). Therefore, we constructed a nomogram using these factors ( Figure 5). The calibration plots and ROC curve (AUC: 0.8128) showed that the derived nomogram performed well, with a high pCR-predictive ability in the validation cohort  (AUC: 0.7199). Interestingly, the predictive accuracy of our nomogram was higher than that of the TNM staging system in both cohorts (AUC: 0.6129 and 0.6043).

Identification of Stromal Immunotype-Associated Biological Pathways and Immune Checkpoint Molecules
By IHC staining, we found three important cytokines (IL-2, INFγ, and TGF-β) that were differently expressed between immunotype A and B subgroups. IL-2, INF-γ expression was significantly higher in the immunotype A subgroup, and TGF-β expression was significantly higher in the immunotype B subgroup ( Figure 6). Furthermore, the expression levels of several immune checkpoint molecules (PD-L1, PD-1, and CTLA-4) were significantly higher in the immunotype A subgroup than in the immunotype B subgroup ( Figure 6). The Gene Set Enrichment Analysis showed that the T-cell receptor signaling pathway, B-cell receptor signaling pathway, and NK cell-mediated immunity were highly enriched in the immunotype A subgroup ( Figure 6).

DISCUSSION
Tumor microenvironment is involved in not only the progression but also the responses to anticancer therapies and outcomes of breast cancer (Yu et al., 2017). TILs in the microenvironment can be used to monitor the immune response and are important in predicting treatment responses in many cancers (Asano et al., 2016). Denkert et al. (Denkert et al., 2018) have demonstrated that increased TIL concentration was associated with good response to neoadjuvant chemotherapy in all molecular subtypes of breast cancer and with a survival benefit in HER2-positive and triple-negative breast cancer. TILs include various types of immune cells. Some types of TILs have been shown to be associated with survival or proliferation of cancer cells, for example, high Treg cell count was related to poor clinical outcomes (McCoy et al., 2015). However, not all TILs have definitive prognostic values. Therefore, we applied the CIBERSORT algorithm (Ali et al., 2016) to estimate the relative proportions of 22 distinct functional subsets of immune cells in ER-negative breast cancer patients. CD8 + T cells, CD4 + T cells, B cells, M1 macrophages, and Treg cells with significant prognostic values were selected using univariate cox regression. Then five selected types of immune cells were used to build a stromal immunotype. The immunotype not only   Frontiers in Molecular Biosciences | www.frontiersin.org October 2021 | Volume 8 | Article 759495 8 immunotype B subgroup had significantly shorter OS and DFS than those in the immunotype A subgroup.
Interleukin-2 (IL-2), a cell growth factor in the immune system, can promote the proliferation of Th0 and CTL. IFN-γ, a multifunctional cytokine, plays an important role in promoting antigen presentation and enhancing the activity of macrophages and natural killer cells. TGF-β regulates tumor-stroma interactions, angiogenesis, and metastasis and shows an immunosuppressive function, inhibiting the proliferation and activity of CTL and B cells, while promoting the proliferation of M2 macrophages (Sia et al., 2017;Xue et al., 2020). Interestingly, in the present study, we found that IL-2 and IFN-γ were increased in the immunotype A patients, whereas the TGF-β-mediated immunosuppression pathway was elevated in immunotype B patients. What's more, the Gene Set Enrichment Analysis showed that the B-cell receptor signaling pathway, T cell-mediated pathway, and NK cellmediated immunity were highly enriched in immunotype A patients. All these findings suggest that immunotype A indicates a hyperimmune state, with a strong tumor-immune cell interaction. CD4 + T-cells in immunotype A may play an important role in recruiting and modulating cytotoxic T-cells in antitumor immunity and CD4 + T-cell activity contributes to the full function of CD8 + cytotoxic T-cells, which is effective for tumor control. While immunotype B indicates a status of immunosuppression for the level of Tregs. These may also explain the prognostic heterogeneity in patients with different immunotypes. In addition, we found that immunotype A was associated with high expression levels of immune checkpoint molecules and cytokines, suggesting that immune checkpoint inhibitor might work well in this subgroup. Therefore, the stromal immunotype might be used to predict responses to immunotherapy.
Recently, many studies focused on biomarkers that were highly associated with treatment responses such as chemotherapy (Lee et al., 2019;Al Amri et al., 2020;Jiao et al., 2020;Valdés-Ferrada et al., 2020;Zhao et al., 2020), and most of them focused on signatures of genes, microRNAs, lncRNAs, and epigenetic biomarkers for the prediction of long-term survival in patients with tumor. However, these signatures cannot be widely used in clinic because of the variability in gene sequencing methods, the inconvenience in the use of assay platforms, and the requirement for specialized analyses. In the present study, we quantified immune cells using IHC staining, which might be easily applied in clinical practice. Moreover, breast cancer has heterogeneous prognosis. To better predict the responses after neoadjuvant chemotherapy, we constructed a nomogram combining clinical characteristics and immunotype. The calibration plots and ROC curve showed that the derived nomogram outperformed the TNM staging system (the 7 th edition).
The present study had several limitations. First, as a retrospective study with relatively small sample sizes, our results need to be validated in a prospective study with larger cohorts. Second, not all of the 22 subsets of immune cells were analyzed. We selected only five subsets using the CIBERSORT method. Other important immunotypes such as neutrophil infiltration might also be valuable, but were not analyzed in the present study. The immunotype classification needs to be optimized with more important immune markers being identified. Third, the underlying mechanisms of the relationship between immunotype and patient prognosis were not well investigated.
The role of the immune profile in the development and invasion of breast cancer needs to be further investigated.
In conclusion, immune cell infiltration in breast cancer is associated with prognosis. We defined two immunotypes by integrating the indicators of the immune cell infiltration, which could be used to predict survival and responses of breast cancer patients to neoadjuvant treatment. The immunotypes might also have significant implications in immunotherapy for the patients who are insensitive to chemotherapy.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

ETHICS STATEMENT
The present study was approved by the Ethics Committee of the Guangdong Provincial People's Hospital. The patients/ participants provided their written informed consent to participate in this study.