The expression patterns of different cell types and their interactions in the tumor microenvironment are predictive of breast cancer patient response to neoadjuvant chemotherapy

Abstract

The tumor microenvironment (TME) is a complex ecosystem of diverse cell types whose interactions govern tumor growth and clinical outcome. While the TME's impact on immunotherapy has been extensively studied, its role in chemotherapy response remains less explored. To address this, we developed DECODEM (DEcoupling Cell-type-specific Outcomes using DEconvolution and Machine learning), a generic computational framework leveraging cellular deconvolution of bulk transcriptomics to associate the gene expression of individual cell types in the TME with clinical response. Employing DECODEM to analyze the gene expression of breast cancer (BC) patients treated with neoadjuvant chemotherapy, we find that the gene expression of specific immune cells (myeloid, plasmablasts, B-cells) and stromal cells (endothelial, normal epithelial, CAFs) are highly predictive of chemotherapy response, going beyond that of the malignant cells. These findings are further tested and validated in a single-cell cohort of triple negative breast cancer. To investigate the possible role of immune cell-cell interactions (CCIs) in mediating chemotherapy response, we extended DECODEM to DECODEMi to identify such CCIs, validated in single-cell data. Our findings highlight the importance of active pre-treatment immune infiltration for chemotherapy success. The tools developed here are made publicly available and are applicable for studying the role of the TME in mediating response from readily available bulk tumor expression in a wide range of cancer treatments and indications.

Figure 1:. The analysis pipeline, DECODEM and DECODEMi.

First, we apply CODEFACS on bulk expression to generate nine cell-type-specific expression profiles and further apply LIRICS to infer the cell-cell interactions (CCIs) present in the tumor microenvironment. Second, we train a multi-stage machine learning pipeline with either the nine cell-type-specific expression profiles to build nine cell-type-specific clinical response predictors (DECODEM) or the CCI profile to build a CCI-based predictor (DECODEMi). In validation, we deploy the abovementioned first stage to procure the cell-type-specific expression or CCIs and directly feed those as inputs to the trained predictors to evaluate model performance.

Figure 2:. The prominent cell types mediating chemotherapy response in breast cancer tumor microenvironment.

A-B. Comparison of prediction scores (A) and model performance (B) in cross-validation with TransNEO (n = 91) for nine cell-type-specific, bulk, and Sammut et al. predictors. R and NR stand for responders and non-responders. The differences between the prediction scores were computed by using a one-tailed Wilcoxon rank-sum test. AUC and AP stand for the area under the receiver operating characteristics curve and average precision (equivalent to the area under precision-recall curve), respectively. Cell types are ranked by their AUC values in a descending order. C. Comparison of model performance in external validation with ARTemis + PBCP (n = 55) for the seven prominent cell types, bulk, and Sammut et al. predictors. Cell types are ranked by their AUC values in a descending order. D. Comparison of model performance in external validation with BrighTNess (n = 106) for the seven prominent cell types and bulk predictors. Cell types are ranked by their AUC values in a descending order.

Figure 3:. The prominent multi-cell-type ensembles mediating chemotherapy response in breast cancer tumor microenvironment.

A-B. Validation performance with ARTemis + PBCP (n = 55) for the five most prominent two-cell-ensembles (A) and three-cell-ensembles (B) along with the bulk. AUC and AP stand for the area under the receiver operating characteristics curve and average precision (equivalent to the area under the precision-recall curve), respectively. Ensembles are ranked by their AUC values in a descending order. C-D. Validation performance with BrighTNess (n = 106) for the five most prominent two-cell-ensembles (C) and three-cell-ensembles (D) along with the bulk. Ensembles are ranked by their AUC values in a descending order. E. The enriched Reactome pathways (FDR-adjusted p-value ≤ 0.2) across the prominent cell types. For each cell type, only the top 10 relevant pathways are displayed. NES and logP stand for the normalized enrichment score from GSEA and log-scaled p-value, respectively.

Figure 4:. The prominent cell-cell interactions (CCIs) mediating chemotherapy response in breast cancer tumor microenvironment.

A. Prediction scores for DECODEMi in cross-validation with TransNEO (n = 94) and in external validation with ARTemis + PBCP (n = 55) and BrighTNess (n = 106). R and NR stand for responders and non-responders, respectively. The differences between the prediction scores were computed by using a one-tailed Wilcoxon rank-sum test. B. Model performance in cross-validation with TransNEO for CCI-based, the top cell-type-specific, and bulk predictors. AUC and AP stand for the area under the receiver operating characteristics curve and average precision (equivalent to the area under precision-recall curve), respectively. C-D. Model performance in external validation with ARTemis + PBCP (C) and BrighTNess (D) for CCI-based, the top cell-type-specific, and bulk predictors. E-F. The 10 most predictive CCIs in external validation with ARTemis + PBCP (E) and BrighTNess (F). Each CCI is displayed as a quadruplet where the first pair containing the ligand and receptor cell types (separated by ‘−’) is separated by ‘::’ from the second pair containing the corresponding ligand and receptor genes (separated by ‘−’). Mean decrease in Gini impurity is the built-in feature importance measure in random forest where a higher value indicates a higher importance and vice versa [98]. Feature directionalities were computed using Fisher’s enrichment analysis. G-H. Prediction scores (G) and AUC value (H) for chemotherapy response prediction in external validation with the SC-TNBC cohort (n = 200, sourced from Zhang et al. [48]) using the top 170 response-relevant CCIs (among B-cells, myeloid, and T-cells) identified in bulk by DECODEMi.

Figure 5:. Generalizability of DECODEM to single-cell (SC) transcriptomics and for TCGA patient survival stratification.

A-B. Prediction scores (A) and the area under the receiver operating characteristics curve (AUC) values (B) in external validation with SC expression of TNBC patients (SC-TNBC) treated with chemotherapy alone (n = 6) or combined with immunotherapy (n = 9) from Zhang et al. [48]. R and NR stand for responders and non-responders, respectively. The differences between the prediction scores were computed by using a one-tailed Wilcoxon rank-sum test. C-D. Kaplan-Meier plots depicting the overall survival (OS; C) and progression-free interval (PFI; D) of 705 BC patients from TCGA (TCGA-BRCA). Patients were stratified into low-risk, likely responder (n = 353) and high-risk, unlikely responder (n = 352) groups by their median DECODEM score computed using the expression of cancer epithelial cells, whereby the low-risk and high-risk groups comprised individuals with scores above and below the median, respectively. The differences between the two curves were computed by using a Log-rank test.