scRNA-seq analysis of the tumor ecosystem in primary breast cancer and paired PBMC
The overall research scheme is shown in Fig. 1. The Seurat package was utilized to construct a comprehensive tumor niche map for breast cancer, enabling us to perform cell classification and identify marker genes. The resulting data was visualized using the Uniform Manifold Approximation and Projection (UMAP) method, revealing 20 distinct clusters (Fig. 2A, D, and E). Despite the initial data being derived from flow sorting of CD45+ cells, we observed a small number of non-immune cells within these clusters, primarily fibroblasts (DCN and COL1A1) and epithelial cells (KRT19 and KRT18) (Fig. 2A, B, and F).
Among the immunocytes, we identified various cell subtypes, including myeloid-derived cells (LYZ and CSF1R), T cells (CD3D and CD3E), B cells (CD19 and CD79A), plasma cells (IGLL1 and MZB1), mast cells (CPA3), plasmacytoid dendritic cells (pDC) (LILRA4 and IRF7), and natural killer (NK) cells (KLRF1 and KLRD1) (Fig. 2B, C, and F). These cell subtypes were present across tumor-infiltrating immune cells (TILs), peripheral blood mononuclear cells (PBMCs), and exhibited variations between individual patients (Fig. 2G).
Importantly, we observed significant heterogeneity in the composition of myeloid cells between breast cancer tumor-infiltrating immune cells and PBMCs, with a higher proportion of myeloid cells found in TILs (Fig. 2G). Previous studies have highlighted the crucial role of myeloid cells as the major contributors to the tumor microenvironment. They are essential for the normal functioning of both the innate and adaptive immune systems and exhibit heterogeneity and dynamic properties, suggesting diverse roles within the breast cancer microenvironment 25.
PBMC and TIL have different myeloid cell immune ecosystems.
We performed unsupervised clustering of myeloid cells, resulting in a total of 11 cell clusters. One cluster with high expression of CD1C, a marker gene, was identified as dendritic cells and referred to as cDC2 based on literature reports. Through further analysis, we defined the clusters with high expression of APOE, C1QC, and APOC1 as macrophages (Fig. 3D). With the development of single cell analysis methods, the functions of various subsets of cells in tumor microenvironment have gradually attracted people's interest. For example, Rodrigo Nalio Ramos found in the functional analysis of tumor-associated macrophage subsets that FOLR2-labeled macrophages were located in the perivascular region of tumor matrix, where they interacted with CD8 T cells. Promote the anti-tumor effect of CD8 T cells. Finally, it was found that the density of FOLR2 macrophages in tumor was positively correlated with better patient survival rate 26. We also try to use this as a research route to find the immunocytological basis that may lead to immune differences by subdividing subgroups. So the remaining clusters were classified as monocytes, resulting in a total of 5 monocyte clusters, 5 macrophage clusters, and one dendritic cell cluster (Fig. 3A).
However, we observed that several cell clusters were missing in the PBMC samples compared to the TIL and PBMC samples (Fig. 3E). This suggests the presence of myeloid cell subpopulations specific to the tumor microenvironment. We speculate that different infiltrated immune cells can interact with other immune components in various ways and play multiple roles in tumor occurrence and development 27.
Regarding the myeloid cells, Mono1 cells exhibited high expression of S100A8/A9, VACN, and FCN1, while Mono2 cells showed high expression of LYPD2 and LST1, and Mono3 cells displayed high expression of GNLY and IL32. Mono4 cells exhibited high expression of CXCL10 and MX1, whereas Mono5 cells demonstrated high expression of THBS1 and EREG (Fig. 3C) (Supplementary material 1). The multifaceted functions of these genes suggest that more research is necessary to understand the role of monocytes in the tumor microenvironment. Take for example. CXCL10 regulates the persistence and heterogeneity of CD8+ T cells and positively promotes the killing effect of T cells 28.
Pseudo-time reveals myeloid cell transition trajectory.
To explore the developmental trajectory of myeloid cells in breast cancer, we employed two methods to infer cell developmental trajectories. Firstly, using Fitdevo, we predicted the developmental potential of subpopulations within the myeloid cells. Notably, Mono4 was positioned at the arrowhead end, suggesting that it represents a terminally differentiated monocyte subpopulation with limited developmental potential (Fig. 3F).
Next, utilizing Monocle, we inferred state trajectories (Fig. 3H) to investigate the dynamic immune and cell developmental states of myeloid cells in tumor-infiltrating immune cells and peripheral blood mononuclear cells (PBMCs) of breast cancer patients. Our analysis revealed that monocytes predominantly occupied the early stage of pseudo-time, while macrophages mainly populated the late stage of pseudo-time (Fig. 3I). This pattern suggests that monocytes undergo differentiation into macrophages, with the majority of monocytes in the early stage of differentiation and macrophages in the late stage. This transformation often involves functional changes in the cells and variations in immune infiltration. For example, previous studies have reported that cytokines such as IL-6 and GM-CSF can drive monocytes to differentiate into M2-type tumor-associated macrophages, ultimately promoting cancer metastasis 29. Furthermore, another study demonstrated that monocytes can differentiate into bone marrow-derived macrophages in the central nervous system, highlighting their potential as gene delivery vectors and therapeutic targets for refractory CNS diseases 30.
Additionally, we analyzed the distribution of cells in the pseudo-time density map and found that Macro1, Macro3, Macro5, and Mono4 were primarily located towards the end of the pseudo-time (Fig. 3G). Using Monocle, we performed BEAM (Branched Expression Analysis Modeling) analysis on branch node 1 and identified highly expressed MHC-II-related genes, such as HLA-DPA, HLA-DPB1, and HLA-DPA, at the right end of the branch point. This suggests their involvement in antigen assembly and presentation by macrophage MHC-II. Moreover, macrophages exhibiting high expression of MHC-II have been associated with nervous system inflammation and behavioral changes 31,32. The high expression of macrophage marker genes C1QC, APOE, and SPP1 on the right side of the branch point further supports the notion that macrophages predominantly reside in the late stage of pseudo-time 33. Conversely, monocyte marker genes, such as S100A8, FCN1, and VCAN, showed high expression on the left side of branch point 1 (Fig. 3J). Notably, the late-stage distribution of Mono4 in the pseudo-time period suggests that it represents a highly differentiated monocyte. The absence of this cell cluster in the PBMC samples suggests its potential importance within the tumor microenvironment (Fig. 3E).
These findings provide valuable insights into the developmental trajectory of myeloid cells in breast cancer, emphasizing the differentiation of monocytes into macrophages and their potential roles within the tumor microenvironment.
Mono4 cells are target cells that respond to interferons.
The highly expressed genes of Mono4 been analyzed to investigate its biological function and the findings revealed that Mono4 cells exhibited high expression of MX1, IFI44L, CXCL10, IFIT3, and IFI T1 (Fig. 4B). Importantly, these genes were almost not expressed in other myeloid cell subpopulations, indicating the strong heterogeneity between Mono4 and other myeloid cell subgroups.
Next, we performed gene enrichment analysis on the highly expressed genes of Mono4 (Supplementary Material 2) and found significant enrichment in biological functions related to response to interferon-γ (IFN-γ), cellular response to type I interferon, positive regulation of type I interferon production, interferon-gamma-mediated signaling pathway, and response to tumor cells (Fig. 4D). This enrichment analysis further supports the notion that Mono4 cells are target cells that respond to interferon action, aligning with previous reports highlighting the essential role of interferon in the antitumor effect 34. Moreover, our ROGUE analysis revealed that Mono4 cells had the highest homogeneity within the monocyte population (Fig. 4C), indicating their widespread presence in the tumor tissues of breast cancer patients.
Furthermore, we performed transcription factor activity analysis to investigate the transcriptional regulation in myeloid cell subpopulations. Our results demonstrated significant enrichment of IRF7, a transcription factor in Mono4 cells. IRF7 plays a critical role in various aspects of tumor immune response, including the differentiation and development of endogenous immune cells and regulation of the response to exogenous pathogens 35. Previous studies have also identified IRF7 as a key transcriptional regulator in interferon (IFN)-dependent immune responses (Fig. 4J) 36.
Pathway enrichment analysis using PROGENy was conducted to explore the activated pathways in Mono4 cells. We found significant enrichment of the JAK-STAT pathway in Mono4 cells (Fig. 4K). Moreover, STAT1, a key component of the JAK-STAT pathway, was found to be significantly enriched in Mono4 cells through transcription factor analysis (Fig. 4J). STAT1 plays a crucial role in tumor immunity and body homeostasis, protecting CD4+ T cells from NK cell clearance in the immune microenvironment to ensure the complete immune function of the body. The STAT1-related pathway is well-known for its potent pro-inflammatory, antiviral, and immune response functions in the peripheral circulation. The absence of STAT1 in hosts leads to rapid death in severe infections 37. These results indicate that Mono4 cells activate their pro-inflammatory response through interferon binding to its receptor and subsequent activation of the JAK-STAT pathway.
Mono4 cells are target cells that respond to interferon-γ (IFN-γ), rather than type I interferons.
To further explore whether Mono4 cells are target cells of type I interferons or interferon-γ (IFN-γ), we examined the expression of IFNA1 and IFNB1 in various cell types. IFNA1 and IFNB1 encode type I interferon alpha and interferon beta, respectively. We did not detect the expression of IFNA1, and IFNB1 showed minimal expression in myeloid cells, almost negligible (Fig. 4E).
Next, we investigated the expression of IFNG, which is the gene encoding interferon-γ (IFN-γ). We found that IFNG was mainly expressed in T cells and NK cells, with the highest expression observed in T cells (Fig. 4F). We then examined the expression of IFNG receptor coding genes, IFNGR1 and IFNGR2, in myeloid cell subpopulations. We observed relatively high expression levels of both genes in Mono4 cells (Fig. 4G and H).
Furthermore, we performed Gene Set Enrichment Analysis (GSEA) on the differentially expressed genes between Mono4 and other myeloid cell subpopulations. We found that the gene set "REACTOME REGULATION OF IFNA IFNB SIGNALING" had a p.adjust value greater than 0.05 in GSEA (Fig. 4I), indicating no statistically significant enrichment. In contrast, the gene set "SANA RESPONSE TO IFNG UP" had a p.adjust value much smaller than 0.05 (Fig. 4I). These results indicate that Mono4 cells are target cells of interferon-γ (IFN-γ) rather than type I interferons.
Cellular interaction of Mono4 with T cells.
T cells, being the principal source of interferon-γ (IFN-γ) secretion, were subjected to unsupervised clustering analysis (Fig. 4L). Subsequently, T cell subtypes were identified based on established cell markers. Notably, we characterized several distinct T cell subtypes, including naive T cells (Naive T,CCR7+), CD8+ effector T cells (CD8+ Teff, GZMB+), T helper cells (Th, CD4+), NKT cells (NK T, NCR3+), and a subset of T cells exhibiting elevated expression of IFIT3 (IFIT3+ T) (Fig. 4M).
To further investigate the cell-cell interactions between Mono4 cells and T cell subpopulations, we utilized CellChat for the analysis of the interaction between Mono4 cells and T cell subpopulations. We observed variations in the number and intensity of cell communication among these cell types (Fig. 5A and B). Subsequently, we analyzed the interaction strength between Mono4 cells as the ligand and other cells, and identified that Mono4 cells exhibited the strongest communication intensity with CD8+ Teff cells and IFIT3+ T cells (Fig. 5C). Moreover, we analyzed the interaction strength between CD8+ Teff cells as the ligand and other cells, and found that CD8+ Teff cells demonstrated the strongest communication intensity with Mono4 cells (Fig. 5D). These findings suggest the presence of cell-cell interactions between Mono4 cells and CD8+ Teff cells as well as IFIT3+ T cells.
In our previous study, we confirmed Mono4 cells as target cells responding to interferon-γ (IFN-γ). Therefore, we analyzed the expression of IFNG in T cell subpopulations, and observed its expression in CD8+ Teff cells and IFIT3+ T cells (Fig. 5E). These results indicate that the interferon-γ (IFN-γ) released by CD8+ Teff cells and IFIT3+ T cells further activates Mono4 cells through the interferon-γ (IFN-γ) receptor on Mono4 cells.
Additionally, we analyzed the CXCL signaling pathway network and identified Mono4 cells as the sender, while CD8+ Teff cells and IFIT3+ T cells act as the receivers in the pathway (Fig. 5F). Through further analysis of the receptor-ligand pairs, we found that Mono4 cells can secrete CXCL10, which binds to CXCR3 on CD8+ Teff cells and IFIT3+ T cells (Fig. 5I). Furthermore, we assessed the expression of CXCR3 in various T cell subpopulations and observed higher expression levels in CD8+ Teff cells and IFIT3+ T cells (Fig. 5J). Previous studies have shown that the chemokine receptor CXCR3 reacts with ligand CXCL10 to promote the persistence and heterogeneity of CD8 T cells. It plays an indispensable role in anti-tumor immunity 38,39.
What’s more, we analyzed the SPP1 signaling pathway network and found that Mono4 cells act as the source of the SPP1 signaling pathway network, targeting cDC2 cells (Fig. 5G). Through receptor-ligand analysis, we discovered that the SPP1 secreted by Mono4 cells can interact with the receptor CD44 on cDC2 cells. Previous studies have shown that OPN (SPP1) can stimulate DC cells to secrete IL12, which further induces T cells to produce interferon-γ 40,41.We detected the expression of IL12 in myeloid cells, and as expected, cDC2 cells had the highest expression (Fig. 5H).
To explore the functionality and biological roles of IFIT3+ T cells, we visualized the top 5 highly expressed genes in the identified T cell subpopulations (Fig. 5K). We observed that IFIT3+ T cells, similar to Mono4 cells, predominantly exhibited high expression of interferon response genes. Furthermore, we conducted gene enrichment analysis and pathway enrichment analysis for the highly expressed genes in IFIT3+ T cells. Similar to Mono4 cells, we found enrichment of the JAK-STAT pathway (Fig. 5M) and interferon-related biological functions (Fig. 5L). Based on these findings, we infer that IFIT3+ T cells are a distinct subset of T cells that respond to interferon-γ (IFN-γ) and possess the capability of autocrine interferon-γ (IFN-γ) secretion.
In summary, based on the aforementioned studies, we have discovered that Mono4 cells, through the interaction with interferon-γ (IFN-γ) secreted by IFIT3+ T cells and CD8+ Teff cells, activate the JAK-STAT pathway by binding to their interferon-γ (IFN-γ) receptors. Furthermore, Mono4 cells regulate the expression of SPP1 and CXCL10 through the transcription factors STAT1 and IRF7. On one hand, Mono4 cells recruit IFIT3+ T cells and CD8+ Teff cells to the tumor microenvironment through CXCL10, thereby exerting their anti-tumor effects. On the other hand, Mono4 cells activate the production of interferon-γ (IFN-γ) in T cells through SPP1. The interferon-γ (IFN-γ) to further activate Mono4 cells. Additionally, IFIT3+ T cells can produce interferon-γ (IFN-γ) through autocrine signaling. Thus, a positive feedback regulation is established between the pro-inflammatory Mono4 cells, IFIT3+ T cells, and CD8+ Teff cells. Based on the high expression of CXCL10 in MONO4, we defined this cell population as CXCL10+ monocytes (Fig. 4A).
Significant increase of CXCL10+monocyte after immunotherapy
To further investigate the effects of immunotherapy on CXCL10+monocyte, we performed single-cell sequencing analysis on myeloid cells from mice treated with aPD1+aCTLA4 immunotherapy. Using the Uniform Manifold Approximation and Projection (UMAP) method, we identified and visualized five clusters, including one cluster of highly expressed cell proliferation-related genes, such as MKI67, which we defined as Cycling (Fig. 6A) (Supplementary material 3). Our analysis also demonstrated the presence of CXCL10+ monocytes in breast cancer tumor tissues of mice, which were characterized by high expression levels of MX1, IFI44L, IFIT3, and IFIT1 (Fig. 6B). This finding is consistent with our previous human single-cell sequencing analysis, indicating a high degree of homology between CXCL10+ monocytes in humans and mice. Differential gene analysis showed that CXCL10+monocyte had significantly higher expression levels of CXCL10, CXCL9, CCL5, CCL8, and CCL6 after immunotherapy (Fig. 6C). Chemokines, which are proteins that promote induction, differentiation, and chemotaxis of immune cells, have gained increasing attention in tumor immune microenvironment studies. CXCL9 and CXCL10, which belong to the chemokine family, have been shown to have various functions, including promoting immunity, regulating immune cell activation and migration, and inhibiting tumors 42. Recent studies have also demonstrated that CXCL9 and CXCL10 can transform cold tumors into hot tumors, improving the efficacy of tumor immunotherapy 43. In addition, CC chemokines play an important role in antitumor activity, including chemotaxis of tumor-infiltrating lymphocytes (TIL).
Our results demonstrated a significant increase in CXCL10+monocyte after immunotherapy, suggesting that CXCL10+monocyte respond to immunotherapy and further activate antitumor immunity (Fig. 7C). These findings indicate that CXCL10+monocyte exhibit stronger antitumor effects following immunotherapy.
CXCL10 + monocyte is associated with patient prognosis.
To investigate the relationship between CXCL10+monocyte and patient prognosis, we conducted survival analysis on the Top 10 highly expressed genes (Supplementary material 1) of CXCL10+monocyte and two monocyte marker genes, CSF1R and LYZ, using the TCGA database 44. We discovered that patients with high expression of these 12 genes had a better prognosis (Fig. 7A) 45. Our single-cell data revealed that CXCL10+monocyte is often accompanied by high MX1, IFI44L, IFIT3 and IFIT1 expression. We analyzed the correlation of these genes in breast cancer and found that these genes were significantly and positively correlated with CXCL10(Fig. 7B). This suggests that CXCL10+monocyte is prevalent in breast cancer.
To verify these findings, we performed multi-color immunofluorescence analysis on tumor tissue sections from patients in early, late, and metastatic stages. We observed a difference in cell density of CXCL10+monocyte among the three groups and selected two groups of patients in the same pathological stage to determine the impact of different CXCL10+monocyte densities on patient prognosis. The results were independently judged by three pathologists at Zhejiang Provincial People's Hospital (Supplementary material 4), and we found that the density of CXCL10+monocyte was significantly higher in the early group than in the other two groups, but there was no significant difference between the late and metastatic groups (Fig. 7C-D). In a comparison of patients with the same pathological stage (phase ⅢA) and the same pathological type (Invasive ductal carcinoma of the breast) (Supplementary material 5), we drew survival curves and found that the prognosis of patients with high-density CXCL10+monocyte expression was significantly better than that of patients with low-density CXCL10+monocyte expression (Fig. 7E) (Supplementary material 6) (Table 1), and the results were statistically significant, consistent with our previous speculation.
Table 1
Comparison of clinical characteristics of patients with high and low CXCL10+Monocyte density
Characteristic | Low density of CXCL10+Monocyte | high density of CXCL10+Monocyte | p |
N | 13 | 13 | |
Age, n (%) | | | 0.539 |
<=65 | 11 (84.62%) | 12 (92.31%) | |
> 65 | 2 (15.38%) | 1 (7.69%) | |
T stage, n (%) | | | 0.658 |
T1 | 4 (30.77%) | 3 (21.43%) | |
T2 | 9 (69.23%) | 10 (78.57%) | |
N stage, n (%) | | | 0.308 |
N1 | 13 (100%) | 12 (92.31%) | |
N2 | 0 (0%) | 1 (7.69%) | |
OS event, n (%) | | | < 0.05 |
Alive | 3 (23.08%) | 8 (53.85%) | |
Dead | 10 (76.92%) | 5 (46.15%) | |