1. single-cell analysis of tumor microenvironment of liver metastases of colorectal cancer.
To gain a better understanding of tumor microenvironment, and investigate how TME response to preoperative chemotherapy in liver metastases of CRC, we performed scRNA-seq of 15 samples from three sites (primary CRC, matched liver metastases and blood) of six CRC patients with liver metastatic disease (Table S1). While patients COL15, COL17 and COL18 have been treated with preoperative chemotherapy, the others are treatment naïve. All patients were classified as microsatellite-stable (MSS) with invasive adenocarcinomas and late-stage (IV) disease (Methods, Table S1).
Viable single cells were sorted and used for droplet-based single-cell RNA sequencing (scRNA-seq). After quality control (see Methods), we obtained transcriptome data for 111,292 single cells from primary CRC (n = 6), matched liver metastases (n = 6), and PBMCs (n = 3) for liver metastases of CRC (Table S2). Then we clustered single cells using shared nearest neighbor (SNN) clustering based on significant principal components (PCs), and visualized cell clusters using t-distributed stochastic neighbor embedding (t-SNE) (Fig. 1a, Methods). The major cell populations (including T cells /natural killer (NK) cells, B/plasma cells, cancer-associated fibroblasts (CAFs), endothelial and myeloid cells) were annotated with canonical marker genes (Fig. 1b), treatment state and tissue origin are mapped in Fig. 1c (also see Figure S1a), selected representative markers of each cell type is presented in Fig. 1d and Fig. 1e and Figure S1b.
To characterize the ecosystems of primary and metastatic tumor in CRC, and better understand how TME response to preoperative chemotherapy in liver metastases of CRC, we focus on major cell types of TME (T/ NK cells, B cells, Myeloid cells and CAFs), for each compartment, we re-centred, scaled, normalized and re-clustered the data. Ultimately, we obtained 28 myeloid clusters (6 dendritic cells, 18 tumor-associated macrophages (TAMs), 1 monocytes, 2 myeloid-derived suppressor-like cells (MDSCs-like) and 1 mast cells), 17 B cell clusters, 10 mesenchymal cell clusters (1 endothelial, 6 cancer-associated fibroblasts and 3 myofibroblasts) and T/NK cell clusters. Each cluster was composed of cells from different patients, and for each cluster, the distribution of tissue and treatment status were different (Figure S1c).
2. Preoperative chemotherapy promotes the activation of B cells in the primary CRC
First, we investigate how preoperative chemotherapy influenced the ability of B cells. Interestingly, B cells were prevalent in primary CRC (account for 0.1% of all stromal cells in primary CRC), but depleted significantly in liver metastases (account for 0.01% of all stromal cells in liver metastases) (Fig. 2a). Sub-clustering of B cells revealed 17 subpopulations (Fig. 2b). Among them, 15 clusters were mature B cells (with eight subsets from tumor lesion and five clusters from peripheral blood, Fig. 2b), which are characterized by highly expression of CD20 (MS4A1). While cluster 3 represents plasma cells which is characterized by highly abundant immunoglobulin (IGHG1, IGHG2, IGHG3, IGHG4 and IGHA2), cluster 16 represent plasmablast, characterized by upregulation of immunoglobulin and cell proliferation markers (e.g. MKI67, CDC20, CDKN3 and CCNB2) (Table S3).
Since few B cells infiltrated in the liver metastases, we mainly focused on B cells in the primary CRC. In primary tumor, B cells were separated into different compartments, the B cell populations from tumor treated with preoperative chemotherapy were quite distinct from that in treatment-naïve tumors (Fig. 2c and 2d). We found that cluster 0, 2, 6 and 11 were enriched specific in treatment-naïve tumors, whereas cluster 1, 4, 9, 10 and 12 were almost exclusively present in treated tumors (Fig. 2c, 2d).
On closer examination of the two groups, we find that they had distinct phenotypes (Table S3). Untreated tumors-derived B cells exhibited a naïve and inflammatory phenotype, with cluster 6 expressing IgD (IGHD), cluster 0 expressing immature B cell marker VPREB3 and cluster 2 expressing inflammatory transcription factor NF-kB (NFKBIA), lipid molecules (e.g. APOE and APOC1) and cytokines (e.g. AREG).Whereas B cells derived from treated tumors (cluster 1, 4, 9, 8 and 12) have an activated immune activation phenotype with upregulation of immunocostimulatory molecules (e.g.,CD82, CD83 and CLECL) and MHCs molecules (e.g.,HLA-DRA, HLA-DRB5) (Table S3). GO analysis also confirms that upregulated genes associated with treated tumors were enriched in antigen processing and presentation (Figure S2a).
Comparing with B cells in untreated tumors, higher expression of immunoglobulin, such as IGLC3, JCHAIN, IGHG1, IGHG3, IGHG4 and IGHA1, were observed in treated tumors (Fig. 2e, 2f, and Figure S2b), implying class switch recombination (CSR, also known as isotype switching) may occur after preoperative chemotherapy. Moreover, the expression of MHCs molecules (such as HLA-A, HLA-DQB1, HLA-DQA2 and HLA-DRB1) were also elevated in treated tumors (Fig. 2g). In line with these results, GO analysis shows that B cells in tumors treated with PC enriched in immune activation signaling pathway, including immune response-activating cell surface receptor, antigen processing and presentation signal pathway, and T cell co-stimulation (Fig. 2h). Whereas B cells in treatment-naïve tumors were enriched in process of responding to reactive oxygen species (ROS), a general feature of inflammation. Thus, those results support the phenotype of B cells in treatment-naïve tumors are associated with inflammation. On the contrary, in tumors treated with preoperative chemotherapy, preoperative chemotherapy promote the activation and the generation of class-switched antibody of B cells.
If B cells can contribute to anti-tumor processes, an effective immune response against tumor progress may be reflected by the presence of a gene expression signature of B cell activation. This lead us to hypothesize that gene expression signatures of activated B cells may be correlated with prognosis of colorectal cancer. To test this, we turn to The Cancer Genome Atlas (TCGA) colon adenocarcinoma (COAD) clinical data, and found the gene signatures of activated B cells were associated with a good prognosis in CRC patients marginally significantly (Fig. 2i, n = 266, including MSI-high, MSI-low and MSS subtype, p = 0.05, Cox regression). Interestingly, the correlation become more prominent in MSS tumors (Fig. 2j, n = 174, MSS subtype, p = 0.034, Cox regression), suggesting that activation of B cells could be more effective in MSS tumors.
Collectively, our data demonstrate that preoperative chemotherapy promotes B cells from a less activated and inflammatory state to a switched and more activated state in the primary tumor of CRC patients with liver metastases, and the activation of B cell could be a potential predictor of effective chemotherapy and good prognosis, especially in patients with MSS CRC. This was also reported in tumors treated with immunotherapy in recently studies (Cabrita et al., 2020; Helmink et al., 2020; Hollern et al., 2019; Petitprez et al., 2020).
3. Preoperative chemotherapy reprogramming TAMs from highly heterogeneity to immature and less activated phenotypes
Myeloid cells are the key components in tumor microenvironment, with an important role in tumor progression and metastasis (Qian et al., 2010).We identified 15,366 myeloid cells, sub-clustered into 28 clusters. Among these myeloid compartments, we designated 18 clusters of tumor-associated macrophages (TAMs), which displayed various features. One cluster of monocytes (M25: FCGR3A) and two clusters of myeloid-derived suppressor-like cells (MDSCs-like, M02 and M16). Six clusters were DCs, including CD1c + DCs (M07 and M10: CD1C), cross-presenting DCs (M21: CLEC9A), pDCs (M17 and M27: LILRA4) and LAMP3 + DCs (M22: LAMP3) (Zhang et al., 2019a). One cluster represents mast cells (M05: TPSAB1). The selected marker genes in myeloid cells were presented in Figure S3a. While TAMs and DCs were enriched in both primary tumor and liver metastases, monocytes and MDSCs-like cells were prevalent in blood and mast cells were mainly enriched in primary CRC (Fig. 3a and 3b).
M02 and M16 enriched the signature of myeloid-derived suppressor cells (MDSCs) (Condamine et al., 2016; Zhao et al., 2012). Consistent with results in Zhang et al., 2019a, S100A family genes are highly expressed, including S100A12, S100A9, S100A8, together with FCN1 and VCAN, but MHC I and MHC II molecules tend to be lowly expressed (Table S3).
Macrophages commonly function as phagocytic cells, which can be activated and display varying phenotype in response to different stimulations (Buckley, 2019; Ngambenjawong et al., 2017). For the diversity and plasticity of TAMs, their heterogeneity and impacts on tumor progression remains largely uncharacterized (Engblom et al., 2016). In this study, we identified a total of 18 clusters of tumor-associated macrophages (TAMs), with no clear delineation between the phenotypes between M1 and M2, the M1 and M2 gene signatures are positively correlated in our TAMs compartment (Figure S3b), indicating that TAMs were more complex than the classical M1/M2 model, consistent with previous studies (Azizi et al., 2018; Cassetta et al., 2019; Davidson et al., 2018; Müller et al., 2017; Zhang et al., 2019a; Zheng et al., 2011). Based on the transcription state and expressed genes of TAMs, we identified various signature genes and classed 18 clusters TAM into four major TAM subsets (Fig. 3c).
The first subset, including cluster 0, 12 and 13, were classified as MHC low TAMs. MHC low TAMs exhibited weak capacity of antigen presentation and immune activation, with lowly expressed MHC I and MHC II genes (e.g. HLA-A, HLA-C, HLA-DRA and HLA-DRB1), and lowly expressed immune-costimulatory genes (e.g. CD80 and CD86). Cluster 12, classified as IL1B + MHC low TAMs, characterized by upregulation of a large number of inflammatory and chemokine genes (e.g., IL1B, IL6, S100A9, S100A8, CXCL8, CXCL3 and CXCL1) which involved in recruiting and regulating of immune cells. Cluster 12 were presented in the primary CRC, in patients who are treatment-naïve (Fig. 3c, right panel, Fig. 3d). Cluster 0 (THBS1 + MHC low TAMs) were prevalent in liver (Fig. 3c, right panel), with highly expression of THBS1, MARCO and genes promoting proliferation of epithelia cells and angiogenesis, such as EREG, AREG and VEGFA, which could stimulate tumor growth and progression. Cluster 13 (MARCO + MHC low TAMs), mainly derived from liver resident Kupffer cells (KCs; the tissue resident macrophage of liver), may exert a tolerogenic or immune inhibitory functions in the liver metastases through MARCO, VSIG4 and CD163.
The second subset, including cluster 6, 9, 14, 18 and 19, were classified as lipid-associated TAMs (LAMs), which were characterized by upregulation of genes involved in lipid metabolism (e.g., APOC1, APOE and LPL), extracellular matrix (ECM) degradation (e.g., MMP7, MMP9 and SPRAC) and complement activation (e.g., C1QA, C1QB and C2). They also highly expressed TREM2 (encoding lipid receptor) and LGALS3 (associated with immune suppression), which were recently found to be associated with metabolic diseases (Jaitin et al., 2019). Among the subset, cluster 6 and cluster 9 showed higher expression levels of matrix metalloproteinases (MMPs) compared to other clusters. GO analysis of upregulated genes in each cluster demonstrated that TAMs were enriched in neutrophil activation and degranulation (Figure S3c). In the primary CRC, lipid-associated TAMs (LAMs) is mainly presented in treatment-naïve tumors, but is shared in treated and untreated in the liver metastases (Fig. 3c, right panel, Fig. 3d).
The third subset, including cluster 4 and 11, were enriched for genes involved in regulation of myeloid leukocyte differentiation, were prevalent in tumors treated with PC (Fig. 3c, right panel, Fig. 3d). With highly expression of FOS, cluster 4, mainly from liver metastases, may be monocyte-derived, which could differentiate into macrophage (Azizi et al., 2018). In addition to the canonical myeloid marker gene CD33, cluster 11 was also associated with highly expression of CD4, which is typically expressed by monocytes and involved in triggering cytokine expression and the differentiation of monocytes into functional mature macrophages (Azizi et al., 2018; Dumitru et al., 2012; Zhen et al., 2014). Moreover, MHC genes tend to be expressed at a higher level in cluster 11 cells than cluster 4, suggesting cluster 11 cells are more mature and activated than the latter. Thus, we classify this subset as monocyte-derived immature TAMs.
The fourth subset, including cluster 20 and 23, are immune regulatory TAMs characterized by upregulation of immune suppressive genes, such as CD274, CCL2, IL10 and TGFB2, and prevalent in primary CRC in treatment-naïve tumors (Fig. 3c, right panel). Strkingly, both of the two clusters are also associated with highly expression of MHC and co-stimulating genes, a signature of immune activation and anti-tumor activities. Together, we uncovered “double-agent” immune regulatory TAMs with co-expression of both immune activated and immune suppressive genes, which indicate complex interactions between anti-tumor and pro-tumor activities.
In addition to the four clusters described above, there are three more clusters, named as M01, M03 and M26, that can not be classified into the above-mentioned five subsets. Genes highly expressed in these three clusters are associated with different bioligical processes and cellular functions. M26 (MKI67 + TAMs) is marked with highly expression of genes in cell proliferation (e.g. MKI67, Table S3). In addition to EGF and MACRO, heat-shock genes are also highly expressed in M01. M03 (CXCL10 + TAMs) highly expressed genes involved in response to interferon-gamma (GBP1, STAT1, IFITM3 and PARP14) and cellular response to zinc ion (MT2A, MT1X and MT1F).
Preoperative chemotherapy reprogramming TAMs from highly heterogeneity to immature and less activated phenotypes
TAMs from treatment-naïve tumors and chemotherapy-treated tumors exhibited distinct phenotypes in both primary tumor and liver metastases. In the primary CRC, TAMs in untreated tumors showed higher heterogeneity, they present distinct phenotypes compared with TAMs in tumors treated with preoperative chemotherapy (Fig. 3d). Focusing on the treatment state, we observed some phenotypes shared between treated and untreated tumors, such as MKI67 + TAMs (M26). However, immature TAMs (M04, M11) and MHC low TAMs (M00, THBS1 + MHC low TAMs) were largely specific to treated tumors (Fig. 3e). Instead, clusters of more activated TAMs (MHC high TAMs, M20 and M23), MMPs + lipid-associated TAMs (M06 and M09), pro-inflammatory TAMs (IL1B + MHC low TAMs, M12), immune suppressive TAMs (CXCL10 + TAMs, M03) and HSPH1 + TAMs (M01) were specific to treatment-naïve tumors (Fig. 3d, and Fig. 3c, right panel).
Concordant with observations in primary tumors, TAMs populations with heterogeneous phenotypes were presented in both treated and untreated tumors in liver metastasis, including lipid-associated TAMs (M06, M09, M14, M18 and M19) and MKI67 + TAMs (M26). But the MHC low TAMs (M00, M13) and immature TAMs populations (M04, M11) were dominantly enriched in tumors treated with preoperative chemotherapy (Fig. 3d, 3e). While immune suppressive TAMs (M03, CXCL10 + TAMs) and HSP + TAMs (M01) were enriched in treatment-naïve tumors (Fig. 3d). Importantly, we found the gene signature of M11 is associated with good prognosis marginally significantly (p = 0.088, Cox regression, Fig. 3f) in TCGA (The Cancer Genome Atlas) COAD (colon adenocarcinoma) patients with MSS subset (n = 174, MSS subtype), but they are not associated with outcomes in MSI CRC tumors, suggesting the infiltration of M11 could be a potential predictor of good prognosis in MSS CRC.
TAMs from tumors treated with PC are distinct from those in treatment-naïve tumors (Figure S3d) at the transcriptome level. In the niche of treatment-naïve tumors, TAMs from the primary CRC were enriched for processes of neutrophil activation, responding to IFN-γ and regulation of fibroblast proliferation (Figure S3d), indicating its proinflammatory phenotype (Lee et al., 2017). TAMs from the metastases enriched for antigen processing and presentation, neutrophil activation and responding to IFN-γ. In the ecosystem of treated tumors, TAMs were characterized by unregulated genes involved in protein targeting to endoplasmic reticulum (ER) and RNA catabolic progress in the primary tumors, while in the metastatic sites, TAMs were enriched for myeloid leukocyte chemotaxis and migration processes, which could be related to monocyte-like TAMs aggregated in this lesion.
In general, our results showed that preoperative chemotherapy suppressed the diversity of TAMs, and only the immature TAMs and THBS1 + MHC low TAMs infiltrated in the primary tumors, while more TAMs aggregated in liver metastases and tended to be immature and less activated. Thus, chemotherapy reprogramming TAMs from highly heterogeneity to immature and less activated phenotypes
4. Preoperative chemotherapy decreases the abundance of ECM-remodeling CAFs and dysfunctional T cells, but promotes the accumulation of myofibroblast
In our study, 1383 CAFs were detected and classified into 9 clusters. Notably, CAFs were significantly more abundant in primary CRC than liver metastatic tumor (Fig. 4a). Using gene expression profile, we classified CAFs into three major subsets (Fig. 4b, c), including secretory CAFs (cluster 0, 1, 6 and 7), ECM-remodeling CAFs (cluster 2 and 8) and contractile CAFs (cluster 3, 4 and 5) (Fig. 4b). The subset of secretory CAFs highly expressed secretory genes such as various growth factor (e.g. IGF1, PDGFD, FGF7 and VEGFB) which mediate angiogenesis and cancer cell proliferation, some signal molecules (e.g. BMP4 and WNT2B) that were able to maintain cancer stem cell niche, complements (e.g. C1S and C3) and chemokines (e.g. CCL2, CXCL12 and CXCL14) that regulated tumor immunity and inflammation. The ECM-remodeling CAFs highly expressed abundant extracellular matrix (ECM) proteins (such as ECM collagens and fibronectin), strongly associated with a fibrotic matrix (Fig. 4c). They also expressed a large number of ECM proteases, which altered ECM and assisted angiogenesis and metastasis (Koliaraki et al., 2017). The contractile CAFs enriched for genes involved in regulation of cell contraction (Fig. 4c), suggesting some distinct phenotypes. Since CAFs have numerous potential cellular sources, we suggested the cluster 4 tended to be myofibroblast nature for their upregulation of myofibroblast markers (e.g. ACTA2 and TAGLN) and genes involved in myogenesis (e.g. MYH11, PLN and CNN1) (Elyada et al., 2019; Koliaraki et al., 2017). Cluster 5, highly expressing pericyte-associated markers (e.g. RGS5 and CSPG4), were largely pericytes source (Koliaraki et al., 2017). Cluster 3 upregulate genes involved in stress response (e.g., JUN, BAG3 and HSPA2) and only present in liver metastases, suggesting they may be triggered as adaptation to tumor microenvironments.
CAFs derived from treated and untreated tumor exhibited distinct phenotypes (Fig. 4b, right panel). ECM-remodeling CAFs were prevalent in the primary CRC in treatment-naïve tumors (Fig. 4d), whereas contractile CAFs were more prevalent in tumors treated with preoperative chemotherapy both in primary (cluster 4) and liver metastases (cluster 3 and cluster 5). Secretory CAFs were observed in the primary CRC, and were mainly enriched in treated tumors (except cluster 0) (Fig. 4d). Comparison the CAFs from treatment-naive tumors and treated tumors, CAFs from treatment-naive tumors were strongly enriched for processes of ECM organization and collagen metabolism (Fig. 4e), whereas CAFs in tumors treated with PC were significantly enriched pathways involved in regulating muscle cell differentiation, immune system (T cell activation) and epithelial cell proliferation (Fig. 4e). As we know, ECM remodeling is one important feature of CAFs common to progressive tumor and promotes metastasis (Bonnans et al., 2014). The observations here indicated that preoperative chemotherapy suppressed ECM remodeling in CAFs, but promotes accumulation of myofibroblast and diverse of secretory CAFs in metastases of CRC.
5. Preoperative chemotherapy suppress the CD8 + dysfunctional T cells.
Here, we also identified different phenotypes of T cells (Fig. 5a), including naïve T cells (TN), central memory T cells (TCM), intraepithelial lymphocytes (IELs), tissue-resident memory T cells (TRM) /effector memory T cells (TEM), recently activated effector memory T cells (TEMRA), dysfunctional or “exhausted” T cells (TEX), TH17-like cells, CXCL13 + TH1-like cells, MKI67 + T cells and regulatory T cells (Tregs). Within these sub-populations, TN, TCM and TEMRA were mainly enriched in blood, while IELs were most exclusively populated with cells form primary cancer, TEX cells and TRM cells were presented both in primary and liver metastasis. The treatment state and tissue origin are mapped in Figure S4a. The annotation was confirmed by expression of canonical markers (Figure S4b, S4c and Method).
In addition, we also identified two MKI67 + CD8 + T cell populations (T36 and T24) in the T cell populations (Fig. 5b). Closer examination of those two clusters revealed cluster 36 were prevalent in the liver metastases (Fig. 5b, 96% in liver and 4% in the primary CRC), whereas cluster 24 were enriched in both the primary CRC and the liver (53% in CRC and 45% in liver), indicating that primary CRC and liver metastases share cluster 24, and the liver metastases has their specific T cells characterized by high proliferation property (T36). To identify differences of those two subpopulations of T cells, we detected differentially expressed genes (DEGs) between T24 and T36 with p-value less than 1%, log2-fold change more than l and more than 10% cells expressing the gene. The results revealed heat-shock protein genes such as HSP90AA1, HSPA6, HSPA1A, HSPA1B and DNAJA4, and some molecular chaperones (e.g. BAG3 and HSPB1) were greatly upregulated in cluster 36 (Fig. 5b). Gene ontology (GO) analysis showed that cluster 36 was mainly enriched for protein folding, and response to stress (response to temperature stimulus, response to heat, response to topologically incorrect protein), suggesting they may be activated as adaptation to tumor microenvironments.
Many studies revealed that chemotherapy could impact the compartment of T cells. For example, chemotherapy could increase TIL (tumor-infiltrating lymphocyte) infiltration and decrease Tregs accumulation and proliferation in CRC patients (Maeda et al., 2011; Terme et al., 2013). Comparison of the cellular diversity of T cells in treatment- naïve tumors and tumors treated with PC, we found most subsets of CD4 + T cell were shared in chemotherapy- treated and untreated tumors (Fig. 5c, 5d), but the phenotypes of CD8 + T cells were significantly different (Fig. 5c, 5d).
In the treatment- naïve tumors, different phenotypes of CD8 + T cell were infiltrated in the primary tumor, including effector T cells and exhausted T cells, but in the metastatic sites, for the CD8 + T cells, only dysfunctional or exhausted T cells were accumulated in the liver (Fig. 5c). Most significantly, preoperative chemotherapy suppressed accumulation of dysfunctional T cell significantly both in the niches of primary CRC and liver metastases (Fig. 5e), which were validated by immunofluorescence (Figure S5). Consistent with previous studies (Maeda et al., 2011; Terme et al., 2013), chemotherapy decrease the accumulation of Tregs in the primary CRC, but in liver metastases, the abundance of Tregs were comparable between treated and untreated tumors (Fig. 5f). The differentiation trajectories of CD8 + and CD4 + T cells (Fig. 2e, 2f, Figure S4e, S4f) also confirm the CD8 + dysfunctional T cells were most prevalent in treatment-naïve tumors, while the Tregs is shared in treated and untreated tumors.
6. Cell-cell cross-talks within the tumor microenvironments in liver metastasis of CRC
Tumor microenvironment are a complex ecosystems, the crosstalks between each cell determine tumor biology and response to therapies (Tirosh and Suvà, 2019). In order to systematically map cellular interactions especially different immune cells and mesenchymal cells in the TME of the primary CRC and metastases, and investigate the potential cellular communication which contribute to the cancer progression, metastases and immune evasion, we used CellPhoneDB v2.0.6 to study the cross-talks between stromal cells in tumor microenvironment. To visualize the cross-talks between different cells types, a chord diagram was built using the circlize package (Gu et al., 2014) in R.
Firstly, we provided a landscape of cross-talks within major stromal cell types both in the niche of the primary CRC and the liver metastases. Then based on each sub cell types annotated above, we focused on sub-cell types, investigated intercellular communications in tumors treated with preoperative chemotherapy and untreated tumors. Selected LR paired were summarized in Fig. 6. The full list of results that were unique to different tumor microenvironment were available in Supplementary Tables (Table S4).
The cross-talks within major stromal cells reveals that TAMs had the most broad cross-talks with other cells, both in the niche of primary CRC and liver metastases (Figure S6a, S6b). They displayed a rich ligands and receptors (LR) profile, broadly communicated with mesenchymal compartment (CAFs and endothelial cells), and immune compartment (including T cells, NK cells, mast cells, DCs and TAMs) in both the primary and liver metastases (Figure S6a). Comparison of the interactions in the primary CRC and liver metastases, we found that TAMs in the primary CRC communicated more frequently with CAFs than that in liver metastases, whereas the cross-talk between TAMs and DCs was greatly increased in the metastatic niche when compared with that in primary CRC (Figure S6a, S6b).
Then on the basis of the annotation above, we focused on cell subpopulations, to identify key interactions between each cell subpopulations.The cross-talks between different phenotypes of TAMs were most abundant in the TME, compared with other cell types. Subsequently, the cross-talks between sub populations of DCs, CAFs, TAMs, dysfunctional T cells and endothelial cells were relatively more frequent. Whereas immune cell such as B cell, plasma cells, other T cell subtypes (including TEM, TRM, TEMRA and Tregs), NK cells and mast cells communicated less in TME (Figure S6c).
Interestingly, among the TAMs compartments, monocyte-derived immature TAMs, MKI67 + TAMs and M03-CXCL10 + TAMs were more activated, with more interactions with other cell types (Figure S6c), hinting that they mediate signals between diverse cell types. Moreover, compared with other T cell subtypes, we observed MKI67 + CD8 + T cells and TEX displayed a rich LR profile (Figure S6c), communicated densely with mesenchymal compartment (CAFs and endothelial cells), and immune compartment, indicating their immune-modulating functions. Both the MKI67 + CD8 + T cells and MKI67 + TAMs displayed a rich LR profile, hinting cells with highly proliferation may have more cell-to-cell communication between different cell types.
Cell-Cell communications in the niches of tumors treated with and without preoperative chemotherapy
Whether preoperative chemotherapy can reprogram the interactions within stromal cells is still unclear. To investigate the effects of chemotherapy on tumor microenvironment of liver metastases of CRC, we further compared cell interaction networks between tumors treated with PC and without PC both in primary tumor and metastatic sites.
Our data reveal the phenotypes of TAMs were highly heterogeneous, and the dysfunctional T cells and ECM-remodeling CAFs were expanded in the primary CRC of treatment-naïve tumors. While in the microenvironment of tumors treated with preoperative chemotherapy, in the primary CRC, chemotherapy promoted the activation of B cells, increased the abundance of immature and less activated TAMs, and expanded the myofibroblasts.
The LR map showed that communications related to immune regulation were more frequently and broadly in the tumors treated preoperative chemotherapy (Fig. 6a). In the niche of treatment-naïve primary tumors, we found that non-PC enriched TAMs (MHC high TAMs (M20 and M23), MHC low inflammatory TAMs (M12)) and MKI67 + TAMs (M26) expressed T cell immune checkpoint ligands CD86, interacted with dysfunctional T cells and Tregs, directly inhibited T cells function (Buchbinder and Desai, 2016; Santarpia and Karachaliou, 2015) through LR pair CD86- CTLA4 (Fig. 6a, left). While in the microenvironment of primary tumors treated with preoperative chemotherapy, the interactions involved in immunomodulation were more densely: PC-enriched TAMs expressed T cell immune checkpoint ligands CD86 and CD80, interacted with dysfunctional T cells and MKI67 + T cells through LR pair CD28- CD86, CD86- CTLA4 and CD80- CTLA4. Moreover, they mediated release of chemokine CCL20 and CXCL16, recruited CXCR3 + CD8 + T cells, CCR6 + T cells and CXCR6 + T cells through CCR6-CCL20, CXCR3-CCL20 and CXCR6-CXCL16 (Fig. 6a, right).
Moreover, preoperative chemotherapy promoted the activation of B cells. In the niche of primary tumors treated with PC, B cells, MKI67 + T cells (T24) and dysfunctional T cells expressed HLA-F, and/ or HLA-DPB1, interacted with PC-enriched TAMs through HLA-F -LILRB2 and HLA-DPB1 -NRG1 (Fig. 6a, right). Importantly, PC-enriched TAMs–M0 also expressed NRG1 and immunomodulatory gene LILRB2, broadly interacted with DCs, TEMA and other TAMs through HLA-F -LILRB2 and HLA-DPB1 -NRG1 (Fig. 6a, right), implying their role involved in immune regulation.
Compared with the untreated tumors, the Notch signaling was exclusively present in tumors treated with PC. The LR map in tumors treated PC showed that myofibroblasts high expressed JAG1, interacted with Notch receptors (NOTCH2, NOTCH3 and NOTCH4) on themselves and endothelial cells (Fig. 6a, right). Notch signaling could regulate myofibroblasts phenotype, tissue fibrosis (Ni et al., 2018), and macrophage differentiation and functions (Palaga et al., 2018). While in the treatment-naïve tumors, we identified pro-angiogenic interactions among ECM- remodeling CAFs, non-PC enriched TAMs and endothelial cells. MHC high TAMs (M20 and M23) and inflammatory TAMs (M12) highly expressed pro-angiogenic factor vascular endothelial growth factor A (VEGFA), activated and recruited ECM- remodeling CAFs and endothelial cells that might generate vascular networks in the microenvironment though NPR2-VEGFA (Hirano et al., 2001; Yeo et al., 2014) (Fig. 6a, left). Notably, among the TAMs compartment, highly expression of VEGFA are found in inflammatory MHC low TAMs (M12), which interacted with lipid-associated TAMs (M14, M18 and M19), MHC high TAMs (M20 and M23) and MKI67 + TAMs (M26) through NPR2-VEGFA (Fig. 6a, left).
Additionally, in the tumors treated with preoperative chemotherapy, we found endothelial cells densely communicated with other cells in the primary tumors. They product ACKR1 (DARC), broadly interacted with dysfunctional T cells, DCs, myofibroblast and TAMs through CCL5-ACKR1, CXCL8-ACKR1, CXCL1-ACKR1 and CCL17-ACKR1 (Fig. 6a, right). ACKR1 has a crucial role in regulation leukocytes recruitment (Pruenster et al., 2009), the high expression of ACKR1 can inhibits tumorigenesis, vascularity and metastasis (Addison et al., 2004). While in untreated tumors, endothelial cells expressed immunosuppressive genes LGALS9, communicated with MKI67 + T cells (T36), MHC high TAMs (M20 and M23) and dysfunctional T cells through LGALS9- HAVCR2 (Fig. 6a, left).
In liver metastases, chemotherapy promoted the abundance of DCs and myofibroblasts, but decreased the dysfunctional T cells.
Both in tumors treated with and without PC, LAMP3 + DCs producing CD86, CD274 (PDL1) and LGALS9, interacted with dysfunctional T cells, MKI67 + T cells and Tregs through CTLA4- CD86, PDCD1-CD274 and LGALS9- HAVCR2 (Fig. 6b). While in the untreated tumors, LAMP3 + DCs also expressed CCL19 and CXCL10, recruit CCR3 + Tregs, dysfunctional T cell and MKI67 + T cells (Fig. 6b, left).
In tumors treated with PC, in line with the observations in the primary tumors, Notch signaling was also identified in the liver metastases after tumors treated with chemotherapy (Fig. 6b, left). Importantly, TNF signaling was unique present in tumors treated with PC. MKI67 + T cells (T36) expressed tumor necrosis factor (TNF) and TNFSF14, broadly interacted with immune cells (DCs and TAMs) and non-immune cells (CAFs) though TNF- TNFRSF1A, TNF- TNFRSF1B and TNFSF14- TNFRSF14 (Fig. 6b, right), which can positively regulate T cell response and contribute to the function of effector T cells in previous studies (Croft, 2009; Zhang et al., 2018).
Comparison of tumors after treatment, in treatment-naive tumors, cross-presenting DCs expressed DPP4, interacted with TAMs through CXCL2 –DPP4, CXCL10 –DPP4 and CCL3L3- DPP4 (Fig. 6b, left). DPP4 (also known as CD26) has been showed to be positive correlated with distance metastasis in colorectal cancer (Lam et al., 2014), and the high expression of DPP4 had significantly worse overall survival in CRC (Lam et al., 2014). In addition, Cross-presenting DCs also interacted with TAMs through LGALS9- HAVCR2, implying their roles in immunosuppression.
In summary, our LR interaction map highlights ACKR1, Notch signaling and molecules mediated immune regulation may contribute to the reprogramming of TME after tumors treated with preoperative chemotherapy.