Single-cell and spatially resolved analysis uncovers cell heterogeneity of breast cancer

Abstract

The heterogeneity and the complex cellular architecture have a crucial effect on breast cancer progression and response to treatment. However, deciphering the neoplastic subtypes and their spatial organization is still challenging. Here, we combine single-nucleus RNA sequencing (snRNA-seq) with a microarray-based spatial transcriptomics (ST) to identify cell populations and their spatial distribution in breast cancer tissues. Malignant cells are clustered into distinct subpopulations. These cell clusters not only have diverse features, origins and functions, but also emerge to the crosstalk within subtypes. Furthermore, we find that these subclusters are mapped in distinct tissue regions, where discrepant enrichment of stromal cell types are observed. We also inferred the abundance of these tumorous subpopulations by deconvolution of large breast cancer RNA-seq cohorts, revealing differential association with patient survival and therapeutic response. Our study provides a novel insight for the cellular architecture of breast cancer and potential therapeutic strategies.

Keywords: Breast cancer; Heterogeneity; Single-nucleus RNA sequencing; Spatial transcriptomics; Tissue architecture.

Fig. 1

snRNA-seq analysis of tumor from sample BC-A. a Schematic of the single-nucleus RNA-seq and ST experiment and analysis. b UMAP visualization of 4,093 nuclei from BC-A tumor analyzed by snRNA-seq showing nine major cell types. c UMAP visualization of inferred epithelial cells from BC-A tumor analyzed by snRNA-seq. Clusters are colored and labeled according to their inferred cell subtypes. d Feature plot of subcluster-specific marker genes in epithelial cells. e Heatmap of differentially expressed genes in each epithelial subcluster. The color bars above the heatmap reflects the subcluster PAM50 subtype estimated by ‘pseudobulk.’ f Heatmap showing the pathway enrichment of each epithelial subcluster using MSigDB HALLMARK gene sets. Mean score of GSVA was z-score transformed. g Heatmap of the area under the curve (AUC) scores of TF motifs estimated per cell by SCENIC. h Differentiation trajectories of epithelial cells by Monocle 2. i Circos plot showing the interactions between ligands and receptors across cell types

Fig. 2

ST analysis of BC-A and cell type deconvolution. a Scaled deconvolution values for six epithelial subclusters overlaid onto tissue spots. b Different spatial distribution of six groups as defined in a. c Heatmap of estimated scores of immune cells in each group by MCPcounter. d Heatmap of the ssGSEA score of six gene signatures estimated in each METABRIC sample. e Kaplan–Meier survival curve for METABRIC cohort in three groups. P value was calculated with log-rank test. Log-rank p value < 0.05 was considered as statistically significant. f Box plot of the estimated proportion of six gene signatures in two NAC cohorts using CIBERSORTx. Statistical significance was determined using a two-sided t-test in a pairwise comparison of means between groups, with P values adjusted using the Benjamini–Hochberg procedure. *P < 0.05, **P < 0.005, ***P < 0.005