https://orcid.org/0009-0001-1226-4607
Figures
Abstract
Colorectal cancer (CRC) is a major cause of cancer mortality and a serious health problem worldwide. Mononuclear phagocytes are the main immune cells in the tumor microenvironment of CRC with remarkable plasticity, and current studies show that macrophages are closely related to tumor progression, invasion and dissemination. To understand the immunological function of mononuclear phagocytes comprehensively and deeply, we use single-cell RNA sequencing and classify mononuclear phagocytes in CRC into 6 different subsets, and characterize the heterogeneity of each subset. We find that tissue inhibitor of metalloproteinases (TIMPs) involved in the differentiation of proinflammatory and anti-inflammatory mononuclear phagocytes. Trajectory of circulating monocytes differentiation into tumor-associated macrophages (TAMs) and the dynamic changes at levels of transcription factor (TF) regulons during differentiation were revealed. We also find that C5 subset, characterized by activation of lipid metabolism, is in the terminal state of differentiation, and that the abundance of C5 subset is negatively correlated with CRC patients’ prognosis. Our findings advance the understanding of circulating monocytes’ differentiation into macrophages, identify a new subset associated with CRC prognosis, and reveal a set of TF regulons regulating mononuclear phagocytes differentiation, which are expected to be potential therapeutic targets for reversing immunosuppressive tumor microenvironment.
Author summary
Colorectal cancer (CRC) is a serious health problem worldwide. Mononuclear phagocytes are the main immune cells in the tumor microenvironment of CRC with remarkable plasticity. In this study, single-cell RNA sequencing technology was used to study the subsets and characteristics of mononuclear phagocytes in CRC. We classify mononuclear phagocytes in CRC into 6 different subsets and demonstrate the heterogeneity in infiltrating mononuclear phagocytes of CRC. We find that tissue inhibitor of metalloproteinases (TIMPs) involved in the differentiation of proinflammatory and anti-inflammatory mononuclear phagocytes and reveal the trajectory and characteristics of circulating monocytes differentiation into tumor-associated macrophages (TAMs). We also find a new immunosuppressive subset of mononuclear phagocytes, which is associated with the prognosis of CRC patients. Our findings identify the heterogeneity in infiltrating mononuclear phagocytes of CRC, find a new subset associated with CRC prognosis. Also, we reveal the differentiation of mononuclear phagocytes would lead to an immunosuppressive tumor microenvironment and provide potential targets for reversing the process.
Citation: Ji T, Fu H, Wang L, Chen J, Tian S, Wang G, et al. (2024) Single-cell RNA profiling reveals classification and characteristics of mononuclear phagocytes in colorectal cancer. PLoS Genet 20(2): e1011176. https://doi.org/10.1371/journal.pgen.1011176
Editor: David J. Kwiatkowski, Brigham and Women’s Hospital, UNITED STATES
Received: May 18, 2023; Accepted: February 8, 2024; Published: February 26, 2024
Copyright: © 2024 Ji et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: The single-cell data generated during the current study are available in GEO database under accession code GSE217774.
Funding: This work was supported by by the National Natural Science Foundation of China (No. 81801578 to TJ, 81974382 to ZW, 82173315 to ZW and 82072743 to GW), the National Key R&D Program of China (No. 2022YFC2408100 to LW) and the Major Scientific and Technological Innovation Projects in Hubei Province (No. 2022BCA013 to ZW). The funders had no role in study design, data collection and analysis, decision to publish or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Colorectal cancer (CRC) is one of the most common malignant tumors worldwide [1]. Although much progress has been made in cancer biology and treatment, the prognosis of CRC patients still remains poor. The process of tumor development is not only related to gene mutations and functional changes of tumor cells, but also affected by tumor microenvironment [2]. As a major component of tumor microenvironment, tumor-associated macrophages (TAMs) are remarkably plastic. TAMs reportedly promote tumor progression, invasion, and dissemination through various mechanisms [3].
TAMs, mainly originated from circulating monocytes, are recruited to tumor sites through growth factors and chemokines secreted by tumor cells. In response to various signals in tumor microenvironment, circulating monocytes differentiate into TAMs, and different types of signals within tumor environment (such as a plethora of tumor cell-derived soluble molecules and metabolites) can result in different TAMs gene transcriptional expression profiles and functions. The current classification divides macrophages into classically activated M1 macrophages and alternatively activated M2 macrophages. However, this classification cannot fully cover all functional subtypes of macrophages, and is not suitable for functional classification of TAMs. Some researchers have divided M2 macrophages into M2a, M2b, M2c and M2d subtypes according to external stimulating signals and their main involvement in immune responses [4–8]. Even so, these classifications are insufficient to cover functional heterogeneity of macrophages. To understand immunological functions of TAMs comprehensively and deeply, the study on how to more accurately cluster TAMs was highly needed.
Single-cell RNA sequencing (scRNA-seq) is a powerful new way to investigate transcriptional profiles of individual cells, which has facilitated the exploration of molecular heterogeneity of TAMs. Using single-cell sequencing, some studies found that a higher level of SPP1+ macrophages was correlated with shorter progression-free survival (PFS) in CRC patients [9,10], indicating that scRNA-seq can distinguish different functional subsets in TAMs from gene expression levels at a single-cell resolution. Therefore, the analysis of different subsets of TAMs in CRC using scRNA-seq would deepen our understanding of tumor microenvironment and point to new directions of cancer therapy.
Monocytes and macrophages are two major components of the mononuclear phagocyte system (MPS) [11]. During MPS differentiation, monocytes differentiate into macrophages, and macrophages exhibit different functions according to different physiological or pathological conditions. Study on MPS would reveal dynamic changes from monocytes to pro-tumorigenic TAMs. Here, our study using scRNA-seq revealed MPS’ heterogeneity and identified an immunosuppressive subset of MPS, which was associated with poor prognosis of CRC patients. We also found possible targets for reversing immunosuppressive tumor microenvironment. This work provides new insights into MPS in CRC and might facilitate the development of immunotherapy for CRC patients.
Results
Single-cell RNA sequencing reveals heterogeneity of mononuclear phagocytes in CRC
To explore the characteristics of monocyte differentiation into TAM in CRC, we collected the tumor tissues of three CRC patients (WUH1, WUH2 and WUH3), and selected the marker gene CD14 for magnetic bead sorting (MACS). The cells sorted were used to perform scRNA-seq (Fig 1A). After quality control (S1A Fig) and removal of the batch effect (S1B Fig), 5153 cells were identified using SingleR package [12]. All cells were confirmed to be MPS, with marker gene of mononuclear phagocytes (CD14) highly expressed and marker genes of intestinal dendritic cells (CX3CR1) and neutrophils (MPO and CEACAM8, also named CD66b) rarely expressed (S1C Fig) [13–15].Then these cells were clustered into six clusters: cluster 1 (C1) to cluster 6 (C6) (Fig 1B).
(A). Schematic illustration of the workflow of the sample preparation and scRNA-seq. (B). UMAP plot of mononuclear phagocytes in CRC. Cells are clustered into six clusters and colored by clusters. (C). Heatmap showed the marker genes of six clusters. Gene expression is color-coded according to a scale based on Z score distribution from −0.5 (blue) to 2 (red). (D). Bar plot showed proportion of different clusters in different samples. (E). Radar charts showed GSVA score of biological process GO term in six clusters.
Different clusters showed unique gene expression profile (Figs 1C and S2A), reflecting MPS heterogeneity. The existences of six clusters were further validated by immunofluorescence staining in CRC (S2B Fig). C1 was characterized by the high expression level of the S100 family (S100A8, S100A9 and S100A12). Since S100A8 and S100A9 are mainly derived from immune cells and participate in inflammatory process [16], C1 could be involved in the inflammatory process of CRC. C2 had high MALT, SLC25A37, and SLC2A3 expression, and C3 had high PLAC8 expression. C4, C5 and C6 had high expression of complement component C1q (C1QA, C1QB, and C1QC). Given that complement component C1q can stimulate phagocytic function of human monocytes and macrophages [17], these clusters might have stronger phagocytic function. C5 also highly expressed APOE, RNASE1, and APOC1, while C6 highly expressed cell cycle genes (TOP2A, MKI67, UBE2C and CDK1), suggesting that C6 is a proliferative cluster. In addition, MPS also exhibited inter-patient heterogeneity (Fig 1D). WHU1, WHU2, and WHU3 had a higher proportion of C5, C1 and C2, respectively. These results indicate that MPS displays a high molecular heterogeneity in CRC patients.
Functional analysis showed different functional characteristics in each cluster (Fig 1E). C1, C2 and C3 were enriched in pathogen-associated molecular patterns (PAMP) receptor signaling pathway, indicating function of innate immunity [18]. While C4, C5 and C6 showed an enrichment in antigen processing and presentation, indicating their role in adaptive immunity [19]. Besides, C1 was enriched in monocyte aggregation, C5 had ability to clear lipoprotein particle, and C6 showed a stronger cell cycle DNA replication.
Clustering of mononuclear phagocytes in CRC is different from that in blood and normal colon
In healthy intestinal tissue, MPS consists of blood-derived monocytes and tissue-resident macrophages [20]. To explore the origins of the identified MPS subsets in CRC, we aligned our data with data of blood and normal colon tissue from healthy people and CRC patients. The integrated data showed similar clusters of MPS (Fig 2A, B), while proportion of each cluster varied in different sample and tissue (Fig 2C and 2D). The main MPS cluster in blood was C1 and blood-derived monocyte marker (S100A8) [20] was highly expressed in blood data (Fig 2E), indicating C1 as blood-derived monocytes. In normal colon, C1 and C5 were main MPS clusters and mucosal macrophages marker (SELENOP) [20] was highly expressed in C5 (Fig 2E), suggesting that C5 was tissue-resident macrophage. Compared with healthy people, normal colon in CRC patient showed a higher proportion of C1. MPS in cancer tissue showed completely different clustering, with proportion of C2 increased and C3, C6 specifically appeared in cancer tissue. These results indicated tumor microenvironment of CRC would result in unique polarization of MPS.
(A). UMAP plot of mononuclear phagocytes from blood, normal colon and CRC from healthy and CRC samples. (B). Dot plot showed marker genes of six clusters. (C and D). Split UMAP plot (C) of mononuclear phagocytes and bar chart (D) showed proportion of six clusters in different sample. (E). Split UMAP plots showed expression levels of markers of intestinal monocytes (S100A8) and mucosal macrophages (SELENOP) in mononuclear phagocytes.
Tissue inhibitor of metalloproteinases (TIMPs) participate in the differentiation of M1-like and M2-like mononuclear phagocytes
The classic classification groups macrophages into M1 or M2 types. Further studies found, according to different stimuli, the M2 macrophages can be subdivided to M2a, M2b, M2c and M2d macrophages [21]. However, markers of these subdivisions [6,22,23] could not distinguish MPS in CRC (Fig 3A). To identify more specific markers that could distinguish MPS in human CRC, we calculated the M1 and M2 polarization scores using related gene sets (Fig 3B) and observed M1-like and M2-like cells present in MPS. C1, C2, C3 were M1-like cells, and C4, C5, C6 were M2-like cells (Fig 3C). We also analyzed a public scRNA-seq dataset of CRC (GSE132465) for validation and found the similar results (S3A–S3C Fig).
(A). UMAP plots showed expression levels of M1 markers (TLR2, TLR4), M2a markers (MRC1, CD209), M2b markers (CD86, CCL1), M2c markers (CD163, MERTK) and M2d markers (VEGFA, IL10). (B). Scatterplot showed M1 polarization and M2 polarization score of each cluster, the colored triangles coded by clusters present mean score of each cluster. (C). UMAP plot of mononuclear phagocytes in CRC, color-coded by M1-like and M2-like cell. (D). Scatterplots showed expression levels of TIMP1 (top) and TIMP2 (bottom). (E). UMAP plots and violin plots showed expression levels of TIMP1 and TIMP2 in M1-like and M2-like cell. (F). Spatial plots showed areas of tumor and stroma based on tumor score in tissue sections (top), and distribution signature of TIMP1 and TIMP2 in CD14-expressed area (bottom). (G). Violin plots showed expression level of TIMP1 and TIMP2 in different CD14-expressed areas. * represent p value < 0.05, ** represent p value < 0.01, *** represent p value < 0.001.
To explore the characteristics of M1-like and M2-like cells, differential gene expression analysis was performed (S3D Fig). Tissue inhibitor of metalloproteinases (TIMPs) were differentially expressed in M1, M2-like cells. TIMP1 and TIMP2 were highly expressed in M1-like and M2-like cells, respectively (Fig 3D and 3E), which was consistently verified when using the public scRNA-seq dataset (S3E Fig). To further verify this finding, we obtained human blood-derived monocytes and induced them into M1 or M2 macrophages. Q-PCR analysis showed that, compared with M2 macrophages, M1 macrophages expressed higher TIMP1 and lower TIMP2 (S3F Fig).
To decipher the spatial distribution of mononuclear phagocytes with high expression of TIMP1 or TIMP2 in the tumor microenvironment, spatial transcriptome data of CRC was analyzed. Tumor area was identified by tumor score and CD14 was used to label mononuclear phagocytes. As the result shown, MPS with high expression of TIMP1 were distributed more within the tumor, while MPS with high expression of TIMP2 were distributed more in stroma (Figs 3F, 3G and S4A, S4B). Besides, spatial distribution of different clusters varied as well (S4C Fig). Despite being considered as secreted endogenous inhibitors of metalloproteinases, TIMP1 was reported to promote cancer progression [24], but its roles in mononuclear phagocytes remained unclear. Our results not only indicate the complexity of mononuclear phagocytes’ classification in human CRC, but also suggest that TIMPs’ involvement in mononuclear phagocytes’ differentiation.
Trajectory analysis reveals the differentiation characteristics of mononuclear phagocytes into immunosuppressive macrophages
Differentiation of mononuclear phagocytes in CRC is a dynamic process [11]. Pseudotime analysis was performed to investigate this ongoing process (Fig 4A), and showed that different clusters were at different stages of differentiation (Figs 4B and S5A). RNA velocity analysis of mononuclear phagocytes showed that C1 was the initial state of differentiation; C3 and C4 were the transitional clusters; C5 was in the terminal state of differentiation (Fig 4C).
(A). Pseudotime plot showed the dynamics of mononuclear phagocytes. (B). Pseudotime plot of mononuclear phagocytes showed stages of each cluster. (C). Developmental trajectory of mononuclear phagocytes inferred by RNA velocity. Arrows indicate developmental trend of each cell. (D). Heatmap showed downregulated inflammatory genes and upregulated immunosuppression-related genes in the differentiation of mononuclear phagocytes.
During the differentiation of mononuclear phagocytes, the expression of inflammatory genes (TNF, IL1A, IL1B, IL23A) and immunosuppression-related genes (CD86, LILRB4, CD274 [25] PDCD1LG2) was gradually decreased and increased, respectively (Figs 4D and S5B), indicating that the immune function of mononuclear phagocytes decreases towards forming an immunosuppressive tumor microenvironment. Here we reveal the differentiation characteristics of mononuclear phagocytes that different clusters were at different stages of differentiation and formation of the immunosuppressive tumor microenvironment.
Changes of transcription factors during of mononuclear phagocytes differentiation
To discover potential targets for inhibiting or reversing the formation of the immunosuppressive tumor microenvironment, we combined SCENIC [26] to analyze the dynamic changes of transcription factors (TFs) and TF regulons (consisting of the TF and their target genes) during differentiation of mononuclear phagocytes. The results showed the relative expression of TFs genes and the normalized activity of TF regulons in each cluster (Fig 5A and 5B). We found that the TF regulons ELF2, RCOR1, CREB1, ETV6, ETV3, CHD1 were downregulated (Fig 5C), and the TF regulons MAF, USF2, NFIC, ETV5, NFATC2, TCF4 were upregulated (Fig 5D).
(A and B). Heatmaps showed the relative expression (z-score) of top TFs (A) and TF regulons (B) in each cluster. Gene expression level is color-coded according to a scale based on Z score distribution from −1.5 (blue) to 2 (red). (C and D). Scatterplots showed downregulated (C) and upregulated (D) TF regulons in the differentiation of mononuclear phagocytes. Correlation coefficient of pseudotime and AUC score of TF regulons was calculated using the Pearson correlation coefficient presented at top right corner of each scatterplot.
Our results reveal the expression characteristics of TFs and TF regulons during mononuclear phagocytes differentiation, providing potential targets for reversing immunosuppressive tumor microenvironments.
Terminally differentiated cluster C5 is an immunosuppressive cluster
Trajectory analysis revealed that the differentiation of mononuclear phagocytes might form an immunosuppressive tumor microenvironment and C5 was in the terminal state of the differentiation, suggesting that C5 plays an important role in the formation of an immunosuppressive tumor microenvironment. To explore the role of C5 in tumor treatment, we analyzed publicly available scRNA-seq dataset from CRC patients treated with or without preoperative chemotherapy (PC) (GSE178318) and CRC patients before and after chemotherapy (GSE232525). The public scRNA-seq dataset showed similar clusters of mononuclear phagocytes with our data (S6A–S6D Fig). PC significantly reduced C5 abundance in CRC (Fig 6A and 6B). CRC before and after chemotherapy showed similar results (Fig 6C and 6D). These patients who have undergone treatment showed significant tumor shrinkage, suggesting that C5 is a potential target cluster for CRC treatment.
(A). UMAP plots of mononuclear phagocytes in CRC with treatment-naïve (left) and treated with PC (right). Cells are colored by clusters. (B). Bar plot showed proportion of different clusters in CRC with treatment-naïve and treated with PC. (C). UMAP plots of mononuclear phagocytes in CRC before chemotherapy (left) and after chemotherapy (right). Cells are colored by clusters. (D). Bar plot showed proportion of different clusters in CRC before and after chemotherapy. (E). Violin plots showed expression levels of LILRB4 (top) and CD274 (bottom) in each cluster. (F). Dot plot showed metabolism activity of different clusters. (G). Pairwise Pearson’s correlation for infiltration patterns of ten major cell types in TCGA CRC cohort.
We further explored the characteristics of C5. Consistent with trajectory analysis, our scRNA-seq dataset showed that C5 highly expressed the immunosuppression-related genes (LILRB4, CD274) (Fig 6E), indicating that C5 is an immunosuppressive cluster. ScMetabolism [27] was applied to explore the metabolic characteristics of C5 and found a significant enhancement in lipid metabolism (Fig 6F). In addition, paired Pearson correlation analysis showed that in CRC tumor microenvironment, the fibroblasts and C5 mononuclear phagocytes were the most highly correlated populations in TCGA CRC cohort (Fig 6G).
Some studies also identified SPP+ macrophages as a pro-tumorigenic/ pro-metastatic cluster in CRC [9,10], but doubts arose as to whether C5 were the same as SPP+ macrophages. The comparison of gene signature between C5 and SPP+ macrophages revealed a clear difference between these two cell clusters (S7A and S7B Fig), verifying C5 as a unique cluster. Together, these results indicate that C5 is an immunosuppressive cluster with a significantly activated lipid metabolism.
C5 is associated with poor prognosis for CRC patients
Immunosuppressive tumor microenvironment can promote cancer progression and affect patient prognosis. To explore the relationship between MPS subsets with CRC prognosis, we correlated our scRNA-seq data to public datasets. Abundance of each cluster in tumor tissue was evaluated with CYBERSORTx, a computational framework for inferring cell type abundance and cell-type-specific transcriptomes from RNA profiles of intact tissue [28]. Results showed only C5 is associated with worse prognosis in CRC (Fig 7A). The public scRNA-seq dataset of CRC showed that APOE and RNASE1, two marker genes of C5, were highly co-expressed in mononuclear phagocytes only (Fig 7B). Survival analyses performed on the TCGA CRC cohort also showed patients with high expression level of APOE and RNASE1 were tightly associated with a worse prognosis (Fig 7B and 7C). These results indicate that C5 is negatively associated with the prognosis of CRC patients.
(A). Kaplan–Meier survival curves of TCGA COAD and READ patients grouped by abundance of six mononuclear phagocytes subsets. (B). Violin plots showed expression levels of APOE (top) and RNASE1 (bottom) in different types of cells. (C). Kaplan–Meier survival curves of TCGA COAD and READ patients grouped by the expression levels of APOE and RNASE1. (D). Kaplan–Meier survival curves of TCGA BLCA, GBM, HNSC, KIRC, LIHC, LUSC, PAAD and STAD patients grouped by C5 abundance. P values were calculated using the log-rank test.
To expand these observations to more tumors, we explored scRNA-seq datasets of bladder cancer (BLCA), glioblastomas (GBM), head and neck cancer (HNSC), kidney clear cell cancer (KIRC), liver cancer (LIHC), lung cancer (LUSC), pancreatic cancer (PAAD) and gastric cancer (STAD). The results showed that MPS in these tumors had similar clusters with C5 present in all tumors (S8A–S8H Fig). Survival analysis showed C5 was associated with a worse prognosis of patients with BLCA, GBM, HNSC, PAAD and STAD, with no significance in LUSC patients, and was associated with better prognosis in KIRC and LIHC (Fig 7D).
Discussion
Many studies have shown that macrophages play an important role in cancer development. CRC scRNA-seq dataset shows that mononuclear phagocytes are highly enriched in cancer tissues. TAMs in CRC are mostly derived from and are maintained by circulating monocytes. A better understanding of how circulating monocytes differentiate into macrophages and the function of various subsets of mononuclear macrophages in CRC may open new avenue for targeting mononuclear macrophages for CRC therapy. Here, scRNA-seq was used to analyze the MPS in CRC. We revealed the trajectory and characteristics during circulating monocytes’ differentiation into TAMs, identified the heterogeneity in infiltrating mononuclear phagocytes of CRC, and found a new cluster of mononuclear phagocytes as an immunosuppressive cluster that was associated with CRC patients’ prognosis. Also, we revealed the differentiation of mononuclear phagocytes would lead to an immunosuppressive tumor microenvironment and provide potential targets for reversing the process.
Consistent with the studies in other tumors [29], our study found the presence of M1-like and M2-like cells in mononuclear phagocytes of CRC, but the traditional M1 markers and M2 markers were unable to effectively distinguish M1-like cells and M2-like cells. Failure of conventional classification underscored the complexity of tumor microenvironment and required more in-depth research. Differential gene expression analysis of M1-like and M2-like cells showed TIMP1 and TIMP2 was highly expressed in M1-like and M2-like cells, respectively. TIMP1 and TIMP2 are the inhibitors of MMPs. As MMPs were reportedly associated with poor prognosis of cancer patients [30], TIMPs were first thought to be disease-protective molecules. Although TIMP1 can be antitumorigenic and antimetastatic functions in cell lines or mice [31–33], clinical findings showed that TIMP1 was consistently elevated in cancerous tissues as well as in blood, and clearly correlated with poor prognosis of most cancer patients [24,34]. This was consistent with the spatial transcriptome data we analyzed, which showed that MPS with high expression of TIMP1 were distributed more within the tumor, while MPS with high expression of TIMP2 were distributed more in stroma. However, clinical results of TIMP2 varied with tumor types, which might be caused by different interactions and non-anti-proteolytic functions of TIMPs. Current studies on TIMPs are limited to tumor cells, and clinical results of TIMPs do not often take the complex components of TME into account. Thus, the relationship between TIMPs and mononuclear phagocytes remained undetermined. Our study found different expression patterns of TIMP1 and TIMP2 in mononuclear phagocytes, suggesting that the TIMPs may participate in differentiation of mononuclear phagocytes. However, how TIMPs in mononuclear phagocytes affect CRC progression remains unclear and awaits further investigation.
Our trajectory analysis revealed different clusters were at different differentiation stages and would form an immunosuppressive tumor microenvironment during the differentiation of mononuclear phagocytes. To find out potential targets that might reverse the immunosuppressive tumor microenvironment, we analyzed the dynamic changes of TF regulons during the process of differentiation. Several members of ETS family of TF regulons showed dynamic changes. TF regulons ELF2, ETV3, ETV6 were downregulated and ETV5 was upregulated during the differentiation. ETV3 was reportedly a transcriptional repressor that block terminal macrophage differentiation in a Ras-dependent manner [35], while the mechanistic relationship between ETS and mononuclear phagocytes was not clear. Our results suggest the important role of ETS in mononuclear phagocytes, which required further validation. Recently, therapeutic strategies based on ETS proteins have made much progress. Small molecule inhibitors of ETS proteins exhibited antitumor effects, possibly opening up novel avenues for CRC therapies [36–38]. Among the TF regulons that are downregulated during macrophage differentiation, CREB1 induces the upregulation of anti-proliferative Ets repressor factor [39], while loss of RCOR1 could enhance monocytic cell survival and self-renewal [40]. Consistent with MAF as a critical controller for immunosuppressive macrophage polarization and function in cancer [41,42], MAF was found as a transcription factor regulon upregulated during macrophage differentiation. Together, these TF regulons could be the potential targets that regulate differentiation of mononuclear phagocytes and further reverse the immunosuppressive tumor microenvironment.
In exploring the origins of the identified MPS subsets in CRC, we found that C1 was blood-derived monocyte and mucosal macrophages marker (SELENOP) was highly expressed in C5, suggesting that C5 was tissue-resident macrophage. Immunosuppressive cluster C5, which highly expresses C1QA, APOE and RNASE1, was at the terminal state of differentiation and played an important role in the formation of an immunosuppressive tumor microenvironment [43,44]. With the unique gene signature, C5 was also associated with poor prognosis for CRC patients. C5 had a strong gene enrichment in extracellular matrix disassembly and a significantly activated in lipid metabolism, which could collectively promote angiogenesis and tumor-promoting microenvironment. Thus, blockade of C5 formation could be beneficial in treating CRC. Cortese N, et al. have reported that GPNMB+ TAM is negatively associated with CRC prognosis. GPNMB+ TAM more often accumulated at the tumor margin, therefore more exposed to tumor immunosuppressive signals, favor tumor growth and metastasis [45]. In our study, C5 was more accumulated at the tumor margin. Whether the C5 subset can promote CRC metastasis remains to be explored.
We were also aware of some limitations of this study. Currently, the prognostic data of public datasets were based on bulk RNA sequencing, and lacked prognostic data in scRNA-seq dataset. These limitations result in deficiency of clinical results for cell subsets. Although some algorithms such as CYBERSORTx were developed to impute cell fractions of bulk RNA samples based on scRNA-seq dataset, it still could not make a clear distinction between two cell subsets with similar gene signatures and could lead to inaccurate results. Therefore, the prognostic data for scRNA-seq dataset is in an urgent need and it would help guide treatment decisions.
In conclusion, our study revealed the characteristics of the MPS in CRC and advanced the understanding of colorectal tumor microenvironment. Further study on TIMPs and TF regulons in mononuclear phagocytes is required towards deepening our understanding of MPS.
Methods
Ethics statement
All experimental procedures were in accordance with the Declaration of Helsinki. Ethical approval was obtained from the Institutional Review Board of Tongji Medical College, Huazhong University of Science and Technology (permission number: S002).
Human specimens
All our samples were obtained from the Union Hospital of Tongji Medical College, Huazhong University of Science and Technology, Wuhan, China. Three patients (WHU1, WHU2, WHU3) were pathologically diagnosed CRC and specimens (1*1*1 cm) were took directly from the tumor and collected for scRNA-seq. Clinical information of patients were included in Table 1. Tumor tissues of five CRC patients were collected for immunofluorescence staining.
Immunofluorescence staining
Immunofluorescence staining of CRC was performed according to the general protocol. After antigen retrieval and autofluorescence quenching, tissue sections were incubated with BSA for 30 min and then incubated with antibody of marker genes of each cluster (Table 2) overnight. Tissue sections were washed with PBS three times and incubated with secondary fluorescent antibodies for 50 min. After washed with PBS three times, nuclei were counterstained with DAPI. The fluorescent photographs were visualized and taken using fluorescence microscope.
Isolation, culture and differentiation of human blood-derived monocytes
To obtain human blood-derived monocytes, we collected peripheral blood from healthy people. PBMC were isolated through Ficoll density gradient centrifugation (Absin, abs930-200mL). Then monocytes were isolated from PBMC using Human Monocyte Isolation Kit (Stem Cell, Catalog #19359). Monocytes were counted and cultured in RIPM medium with M-CSF (25ng/mL) (Absin, UA040016) for 4 days. Then LPS (10 ng/mL) (Sigma, L5293) and IFN-γ (10 ng/mL) (Absin, UA040053) was added to induce monocytes to differentiate into M1 macrophages. And IL-4 (10 ng/mL) (Absin, UA040026) was added to induce monocytes to differentiate into M2 macrophages. M1 and M2 macrophages were collected after 2 days of culture.
Quantitative RT-PCR
Total RNA was extracted from cultured cells by the RNA-easy Isolation Reagent (Vazyme, R701-02). Then total RNA was used to synthesize cDNA with the M-MLV (H-) Reverse Transcriptase System (Vazyme, R021-01). The cDNA was subjected to RT-qPCR by the AceQ qPCR SYBR Green Master Mix (Vazyme, Q111-02). The analysis of target gene expression was followed by 2–ΔΔCT method. GraphPad Prism 8.3.0 was used for analysis, and all summary data are shown as the mean ± SEM. Student’s t tests were assessed to compare the two different groups. The specific primers were listed in Table 3.
Sorted single-cell suspension preparation
Tumor tissues were immediately cut into approximately 1 mm3 small pieces in RPMI medium after collected from CRC patients. Then each sample was enzymatically digested using Collagenase IV (1 mg/mL), DNase I (100 μg/mL) and hyaluronidase (100 μg/mL) on a 37°C shaker for 1 h. Subsequently, suspension was filtered using a 70 μm cell strainer and centrifuged at 200 g for 5 min. The supernatant was discarded, and then the cells were washed with PBS twice. After the supernatant was removed, the cell pellet was suspended in 3 mL red blood cell lysis buffer and incubated at 4°C for 5 min to lyse red blood cells. With 3 mL RPMI medium added, suspension was centrifuged at 200 g for 5 min. After washing twice with PBS, supernatant was removed and the cell pellet was resuspended in buffer (PBS with 0.5% BSA and 2 mM EDTA).
Single-cell suspension was sorted using magnetic bead sorting according to the manufacturer’s protocol. Cell count of suspension was counted and then Centrifuge cell suspension at 200 g for 5 min. After removing supernatant completely, cell pellet was resuspended to 107 total cells in 100 μL of buffer and then incubated with biotinylated Anti-Human-CD14 at 4°C for 5 min. Next, suspension was washed with buffer twice and centrifuged at 200 g for 5 min. After supernatant moved, cell pellet was resuspended in 80 μL of buffer per 107 total cells and incubated with Anti-Biotin MicroBeads at 4°C for 15 min. Then suspension was washed with buffer twice and centrifuged at 200 g for 5 min. Cell pellet was resuspend to 108 cells in 500 μL of buffer after supernatant moved. Labeled suspension was applied onto MACS column that rinsing with the appropriate amount of buffer. After unlabeled cells pass through, column was washed with buffer three times and placed on a suitable collection tube. Add the 5 mL buffer onto the column and immediately flush out the magnetically labeled cells by firmly pushing the plunger into the column. Finally, the sorted single-cell suspension was counted and diluted in an appropriate concentration.
Single-cell RNA sequencing
scRNA-seq was performed according to the manufacturer’s protocol. Single-cell suspensions were mixed with Chromium Single cell 3′ Reagent v3 kits and then were loaded on Chromium Controller to generate single-cell gel beads in emulsions (GEMs). Cell lysis and reverse transcription reactions were performed in GEMs to generate barcoded full-length cDNA. Then GEMs were disrupted and cDNA was amplificated. The amplificated cDNA was fragmented, end-repaired, A-tailed, and ligated to an index adaptor to generate barcoded scRNA-seq libraries. Subsequently, scRNA-seq libraries were sequenced on the Illumina NovaSeq platform and raw sequencing data were generated for downstream analyses.
Statistics
scRNA-seq data processing and quality control.
Raw data was processed using Cell Ranger toolkit. With cellranger count function, raw data was mapped to the reference genome GRCh38 and generated the gene-cell UMI matrix for downstream analyses.
The UMI matrix was converted into a Seurat object using the R package Seurat2 (version 4.1.0) [46]. To filter out possible empty droplets, low-quality cells, and possible multiplets, cells with <200 or >6000 features and with <500 UMI counts were excluded. Also, to filter out cells of low quality, cells with >20% of their transcripts coming from mitochondrial genes were excluded from analysis. Cell-type annotation was performed based on R package SingleR (version 1.8.1) [12]. Monocytes and macrophages were selected as the research object. Finally, 5153 cells were remained for downstream analyses.
As the dissociation and digestion during single-cell suspension preparation would lead to high levels of immediate early genes (Fos, Jun and other activating protein 1 complex genes) and heat-shock proteins (HSPs) [47], these dissociation and stress associated genes were excluded from further analysis [48].
To remove potential batch effect, data from three patients were integrated using a Seurat approach. The UMI matrix was normalized and separated by patient sample. With FindIntegrationAnchors function of R package Seurat, Top 2000 high variable genes were used to create anchors. Then, we use IntegrateData function to these matrixes were integrated into a new matrix without batch effect.
The principal component analysis (PCA) was performed to reduce dimensionality of the integrated matrix. Top 30 principal components determined by ElbowPlot function were used in the downstream analyses. With resolution = 0.2 and Louvain clustering, six clusters were identified using FindClusters function. And then clustering results were visualized with two-dimensional UMAP plots.
Bioinformatics analysis of single-cell sequence data.
Our single-cell sequence data, public scRNA-seq datasets of blood from healthy people (GSE167363), colon from healthy people (GSE214695), blood and normal colon from CRC patients (GSE161277), CRC samples (GSE132465), CRC patients treated with or without PC (GSE178318), CRC patients before and after chemotherapy (GSE232525), BLCA patients (GSE190888), GBM patients (GSE162631), HNSC patients (GSE234933), KIRC patients (GSE121636), LIHC patients (GSE242889), LUSC patients (GSE123904), PAAD patients (GSE214295), STAD patients (GSE234129) and spatial transcriptomics data from GSE226997 in Gene Expression Omnibus (GEO) were used for bioinformatics Analysis.
Functional enrichment analysis.
We performed functional enrichment analysis to investigate biological characteristic of different clusters. Gene Set Variation Analysis (GSVA) of biological process GO term was performed using R package GSVA (version 1.42.0).
Polarization score of MPS.
To evaluate M1/M2 polarization ability of mononuclear phagocytes. We applied GSVA to calculate M1/M2 polarization score. The reference gene set of M1/M2 polarization was described by Sun et al [49] and all the parameters were set as default.
Trajectory analysis.
Pseudotime analysis and RNA velocity analysis were used to infer developmental trajectory of mononuclear phagocytes. Pseudotime analysis was performed based on R package Monocle2 (version 2.22.0) [50]. The UMI matrix was converted into a CellDataSet object and genes with mean expression < 0.1 were filtered out. Then significantly changed genes identified by the differentialGeneTest function were used to order mononuclear phagocytes in pseudotime analysis. Cell trajectories were constructed and visualized with plot_cell_trajectory function. Key genes during the differentiation progress of mononuclear phagocytes were detected with Monocle2 plot_pseudotime_heatmap function. RNA velocity analysis was performed based on python package velocyto (version 0.17) and python package scvelo (version 0.2.5) [51]. In this process, spliced reads and unspliced reads of single-cell sequence data were calculated based on aligned bam files. Then RNA velocity values were estimated with dynamical model and the velocities were visualized in two-dimensional plot.
Regulatory network analysis.
To explore gene regulatory network of mononuclear phagocytes, we performed single-cell gene regulatory network analysis using R package SCENIC (version 1.3.0) [26]. Following SCENIC tutorial, candidate regulatory modules were inferred from co-expression patterns between genes, and co-expression regulatory modules were refined using TF regulons motif information. Then, AUC score of each TF regulons in each cell was calculated with AUCell function.
Metabolic analysis.
To explore metabolism activity of clusters, we applied scMetabolism (version 0.2.1), a new algorithm used to quantify metabolism activity at the single-cell resolution [27], on single-cell sequence data. Metabolism score of each cell was calculated with sc.metabolism.Seurat function.
Correlation analysis.
Ten major cell types of TCGA CRC cohorts were imputed with CYBERSORTx, then Pearson correlation analysis was performed to assess the relationship among the proportions of ten cell types.
Survival analysis.
Transcriptome data and clinical data from TCGA were used to evaluate the prognostic performance of genes and cell clusters. The TCGA data were obtained from UCSC Xena (http://xena.ucsc.edu/). Transcriptome data was normalized to TPM data. The optimal cutpoint was determined by R package survminer (version 0.4.9). Then samples were grouped into high expression and low expression groups according to the cut-off expression of the gene. Kaplan-Meier survival curves was used to evaluate the correlation between genes and survival.
To evaluate the role of cell clusters in predicting the prognosis of cancers, we apply CYBERSORTx to impute abundance of cell clusters [28]. Then samples were grouped into high abundance and low abundance groups based on cut off value determined by abundance. Kaplan-Meier survival curves was also used to evaluate the correlation between abundance and survival.
Supporting information
S1 Fig. Six clusters of MPS are Identified.
(A). Quality control of single-cell RNA sequence data. (B). UMAP plots of mononuclear phagocytes before integration (left) and after integration (right). Cells are colored by samples. (C). UMAP plots showed expression levels of marker genes of mononuclear phagocytes (CD14), intestinal dendritic cells (CX3CR1) and neutrophils (MPO and CEACAM8, also named CD66b).
https://doi.org/10.1371/journal.pgen.1011176.s001
(TIF)
S2 Fig. Marker genes of mononuclear phagocytes subsets are verified in CRC tissues.
(A). UMAP plots showed expression levels of marker genes of each cluster. (B). Immunofluorescence staining of CD14 and marker genes of each cluster in CRC tissues confirmed the existence of six clusters. White arrows point cells of each cluster. Scale bars, 20 μm.
https://doi.org/10.1371/journal.pgen.1011176.s002
(TIF)
S3 Fig. M1-like and M2-like cell of MPS is verified in public scRNA-seq dataset and human peripheral blood-derived macrophages.
(A). UMAP plot of mononuclear phagocytes in GSE132465 CRC patients. Cells are colored by clusters. (B). Scatterplot showed M1 polarization and M2 polarization score of each cluster in GSE132465 CRC patients. The mean score of each cluster is represented by a color-coded triangle. (C). UMAP plot of mononuclear phagocytes in GSE132465 CRC patients, color-coded by M1-like and M2-like cell. (D). Heatmap displaying differential expression pattern of genes in M1-like and M2-like cells. (E). UMAP plots showed expression levels of TIMP1 and TIMP2 in GSE132465 CRC patients. (F). Expression levels of TIMP1 and TIMP2 in human peripheral blood-derived M1 and M2 cells (N = 3), *P < 0.05. The statistical significance was obtained by Student’s t-tests.
https://doi.org/10.1371/journal.pgen.1011176.s003
(TIF)
S4 Fig. Spatial distribution of MPS in CRC tissue.
(A). Spatial plots of another CRC sample showed areas of tumor and stroma based on tumor score in tissue sections (top), and distribution signature of TIMP1 and TIMP2 in CD14-expressed area (bottom). (B). Violin plots showed expression level of TIMP1 and TIMP2 in different CD14-expressed areas. (C). Spatial feature plots showed distribution signature of six clusters in tissue sections, respectively. ns represent no significance, * represent p value < 0.05, ** represent p value < 0.01, *** represent p value < 0.001.
https://doi.org/10.1371/journal.pgen.1011176.s004
(TIF)
S5 Fig. Dynamic changes of clusters and genes during the differentiation of MPS.
(A). Split pseudotime plot of mononuclear phagocytes showed stages of each cluster. (B). Scatterplots showed downregulated inflammatory genes and upregulated immunosuppression-related genes in the differentiation of mononuclear phagocytes. Cells are colored by clusters.
https://doi.org/10.1371/journal.pgen.1011176.s005
(TIF)
S6 Fig. Similar clusters are identified in MPS of public scRNA-seq dataset.
(A). UMAP plot of mononuclear phagocytes in GSE178318 CRC patients. Cells are colored by clusters. (B). Heatmap showed similar expression characteristics of marker genes between GSE178318 and our data. (C). UMAP plot of mononuclear phagocytes in GSE232525 CRC patients. Cells are colored by clusters. (D). Heatmap showed similar expression characteristics of marker genes between GSE232525 and our data.
https://doi.org/10.1371/journal.pgen.1011176.s006
(TIF)
S7 Fig. C5 is a unique cluster different from SPP+ macrophages.
(A). UMAP plot showed expression level of SPP1 and RNASE1 in MPS. (B). Split UMAP plots showed expression levels of SPP1 in each cluster.
https://doi.org/10.1371/journal.pgen.1011176.s007
(TIF)
S8 Fig. Identification of mononuclear phagocytes subsets in eight types of tumor.
(A-H). UMAP plot (left) and dot plot (right) showed mononuclear phagocytes subsets in BLCA (A), GBM (B), HNSC (C), KIRC (D), LIHC (E), LUAD (F), PAAD (G) and STAD (H).
https://doi.org/10.1371/journal.pgen.1011176.s008
(TIF)
References
- 1. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, et al. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin. 2021;71 (3), 209–249. pmid:33538338
- 2. Hinshaw DC, Shevde LA. The Tumor Microenvironment Innately Modulates Cancer Progression. Cancer Res. 2019;79 (18), 4557–4566. pmid:31350295
- 3. Qian BZ, Pollard JW. Macrophage diversity enhances tumor progression and metastasis. Cell. 2010;141 (1), 39–51. pmid:20371344
- 4. Duan J, Liu X, Wang H, Guo SW. The M2a macrophage subset may be critically involved in the fibrogenesis of endometriosis in mice. Reprod Biomed Online. 2018;37 (3), 254–268. pmid:30314882
- 5. Anderson CF, Mosser DM. A novel phenotype for an activated macrophage: the type 2 activated macrophage. J Leukoc Biol. 2002;72 (1), 101–6. pmid:12101268
- 6. Wang LX, Zhang SX, Wu HJ, Rong XL, Guo J. M2b macrophage polarization and its roles in diseases. J Leukoc Biol. 2019;106 (2), 345–358. pmid:30576000
- 7. Avdic S, Cao JZ, McSharry BP, Clancy LE, Brown R, Steain M, et al. Human cytomegalovirus interleukin-10 polarizes monocytes toward a deactivated M2c phenotype to repress host immune responses. J Virol. 2013;87 (18), 10273–82. pmid:23864618
- 8. Duluc D, Delneste Y, Tan F,Moles MP,Grimaud L,Lenoir J, et al. Tumor-associated leukemia inhibitory factor and IL-6 skew monocyte differentiation into tumor-associated macrophage-like cells. Blood. 2007;110 (13), 4319–30. pmid:17848619
- 9. Zhang L, Li Z,Skrzypczynska KM, Fang Q, Zhang W, O’Brien SA, et al. Single-Cell Analyses Inform Mechanisms of Myeloid-Targeted Therapies in Colon Cancer. Cell. 2020;181 (2), 442–459 e29. pmid:32302573
- 10. Qi J, Sun H, Zhang Y, Wang Z, Xun Z, Li Z, et al. Single-cell and spatial analysis reveal interaction of FAP(+) fibroblasts and SPP1(+) macrophages in colorectal cancer. Nat Commun. 2022;13 (1), 1742. pmid:35365629
- 11. Hume DA, Irvine KM, Pridans C. The Mononuclear Phagocyte System: The Relationship between Monocytes and Macrophages. Trends Immunol. 2019;40 (2), 98–112. pmid:30579704
- 12. Aran D, Looney AP, Liu L, Wu E, Fong V, Hsu A, et al. Reference-based analysis of lung single-cell sequencing reveals a transitional profibrotic macrophage. Nat Immunol. 2019;20 (2), 163–172. pmid:30643263
- 13. Geissmann F, Manz MG, Jung S, Sieweke MH, Merad M, Ley K. Development of monocytes, macrophages, and dendritic cells. Science. 2010;327 (5966), 656–61. pmid:20133564
- 14. Arnold IC, Mathisen S, Schulthess J, Danne C, Hegazy AN, Powrie F. CD11c(+) monocyte/macrophages promote chronic Helicobacter hepaticus-induced intestinal inflammation through the production of IL-23. Mucosal Immunol. 2016;9 (2), 352–63. pmid:26242598
- 15. Tian S, Chu Y, Hu J, Ding X, Liu Z, Fu D, et al. Tumour-associated neutrophils secrete AGR2 to promote colorectal cancer metastasis via its receptor CD98hc-xCT. Gut. 2022;71 (12), 2489–2501. pmid:35086885
- 16. Wang S, Song R, Wang Z, Jing Z, Wang SMa J. S100A8/A9 in Inflammation. Front Immunol. 2018;9, 1298. pmid:29942307
- 17. Benoit ME, Clarke EV, Morgado P, Fraser DA, Tenner AJ. Complement protein C1q directs macrophage polarization and limits inflammasome activity during the uptake of apoptotic cells. J Immunol. 2012;188 (11), 5682–93. pmid:22523386
- 18. Nguyen TA, Pang KC, Masters SL. Intercellular communication for innate immunity. Mol Immunol. 2017;86, 16–22. pmid:27769523
- 19. Pishesha N, Harmand TJ, Ploegh HL. A guide to antigen processing and presentation. Nat Rev Immunol. 2022;22 (12), 751–764. pmid:35418563
- 20. Hegarty LM, Jones GR, Bain CC. Macrophages in intestinal homeostasis and inflammatory bowel disease. Nat Rev Gastroenterol Hepatol. 2023;20 (8), 538–553. pmid:37069320
- 21. Roszer T. Understanding the Mysterious M2 Macrophage through Activation Markers and Effector Mechanisms. Mediators Inflamm. 2015;2015, 816460. pmid:26089604
- 22. Zhang Q, Sioud M. Tumor-Associated Macrophage Subsets: Shaping Polarization and Targeting. Int J Mol Sci. 2023;24 (8). pmid:37108657
- 23. Lai W, Xian C, Chen M, Luo D, Zheng J, Zhao S, et al. Single-cell and bulk transcriptomics reveals M2d macrophages as a potential therapeutic strategy for mucosal healing in ulcerative colitis. Int Immunopharmacol. 2023;121, 110509. pmid:37369160
- 24. Eckfeld C, Haussler D, Schoeps B, Hermann CD, Kruger A. Functional disparities within the TIMP family in cancer: hints from molecular divergence. Cancer Metastasis Rev. 2019;38 (3), 469–481. pmid:31529339
- 25. Pardoll DM. The blockade of immune checkpoints in cancer immunotherapy. Nat Rev Cancer. 2012;12 (4), 252–64. pmid:22437870
- 26. Aibar S, Gonzalez-Blas CB, Moerman T, Huynh-Thu VA, Imrichova H, Hulselmans G, et al. SCENIC: single-cell regulatory network inference and clustering. Nat Methods. 2017;14 (11), 1083–1086. pmid:28991892
- 27. Wu Y, Yang S, Ma J, Chen Z, Song G, Rao D, et al. Spatiotemporal Immune Landscape of Colorectal Cancer Liver Metastasis at Single-Cell Level. Cancer Discov. 2022;12 (1), 134–153. pmid:34417225
- 28. Newman AM, Steen CB, Liu CL, Gentles AJ, Chaudhuri AA, Scherer F, et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol. 2019;37 (7), 773–782. pmid:31061481
- 29. Azizi E, Carr AJ, Plitas G, Cornish AE, Konopacki C, Prabhakaran S, et al. Single-Cell Map of Diverse Immune Phenotypes in the Breast Tumor Microenvironment. Cell. 2018;174 (5), 1293–1308 e36. pmid:29961579
- 30. Kessenbrock K, Plaks V, Werb Z. Matrix metalloproteinases: regulators of the tumor microenvironment. Cell. 2010;141 (1), 52–67. pmid:20371345
- 31. Albini A, Melchiori A, Santi L, Liotta LA, Brown PD, Stetler-Stevenson WG. Tumor cell invasion inhibited by TIMP-2. J Natl Cancer Inst. 1991;83 (11), 775–9. pmid:1645772
- 32. Khokha R. Suppression of the tumorigenic and metastatic abilities of murine B16-F10 melanoma cells in vivo by the overexpression of the tissue inhibitor of the metalloproteinases-1. J Natl Cancer Inst. 1994;86 (4), 299–304. pmid:8158685
- 33. Baker AH, George SJ, Zaltsman AB, Murphy G, Newby AC. Inhibition of invasion and induction of apoptotic cell death of cancer cell lines by overexpression of TIMP-3. Br J Cancer. 1999;79 (9–10), 1347–55. pmid:10188875
- 34. Jackson HW, Defamie V, Waterhouse P, Khokha R. TIMPs: versatile extracellular regulators in cancer. Nat Rev Cancer. 2017;17 (1), 38–53. pmid:27932800
- 35. Klappacher GW, Lunyak VV, Sykes DB, Sawka-Verhelle D, Sage J, Brard G, et al. An induced Ets repressor complex regulates growth arrest during terminal macrophage differentiation. Cell. 2002;109 (2), 169–80. pmid:12007404
- 36. Pop MS, Stransky N, Garvie CW, Theurillat JP, Hartman EC, Lewis TA, et al. A small molecule that binds and inhibits the ETV1 transcription factor oncoprotein. Mol Cancer Ther. 2014;13 (6), 1492–502. pmid:24737027
- 37. Erkizan HV, Kong Y, Merchant M, Schlottmann S, Barber-Rotenberg JS, Yuan L, et al. A small molecule blocking oncogenic protein EWS-FLI1 interaction with RNA helicase A inhibits growth of Ewing’s sarcoma. Nat Med. 2009;15 (7), 750–6. pmid:19584866
- 38. Mohamed AA, Xavier CP, Sukumar G, Tan SH, Ravindranath L, Seraj N, et al. Identification of a Small Molecule That Selectively Inhibits ERG-Positive Cancer Cell Growth. Cancer Res. 2018;78 (13), 3659–3671. pmid:29712692
- 39. Sawka-Verhelle D, Escoubet-Lozach L, Fong AL, Hester KD, Herzig S, Lebrun P, et al. PE-1/METS, an antiproliferative Ets repressor factor, is induced by CREB-1/CREM-1 during macrophage differentiation. J Biol Chem. 2004;279 (17), 17772–84. pmid:14754893
- 40. Yao H, Goldman DC, Fan G, Mandel G, Fleming WH. The Corepressor Rcor1 Is Essential for Normal Myeloerythroid Lineage Differentiation. Stem Cells. 2015;33 (11), 3304–14. pmid:26119982
- 41. Chen Z, Zhou L, Liu L, Hou Y, Xiong M, Yang Y, et al. Single-cell RNA sequencing highlights the role of inflammatory cancer-associated fibroblasts in bladder urothelial carcinoma. Nat Commun. 2020;11 (1), 5077. pmid:33033240
- 42. Liu M, Tong Z, Ding C, Luo F, Wu S, Wu C, et al. Transcription factor c-Maf is a checkpoint that programs macrophages in lung cancer. J Clin Invest. 2020;130 (4), 2081–2096. pmid:31945018
- 43. Donadon M, Torzilli G, Cortese N, Soldani C, Di Tommaso L, Franceschini B, et al. Macrophage morphology correlates with single-cell diversity and prognosis in colorectal liver metastasis. Journal of Experimental Medicine. 2020;217 (11). pmid:32785653
- 44. Murray PJ. Cancer metastasis linked to macrophage size, shape, and metabolism. Journal of Experimental Medicine. 2020;217 (11). pmid:32997733
- 45. Cortese N, Carriero R, Barbagallo M, Putignano AR, Costa G, Giavazzi F, et al. High-Resolution Analysis of Mononuclear Phagocytes Reveals GPNMB as a Prognostic Marker in Human Colorectal Liver Metastasis. Cancer Immunol Res. 2023;11 (4). pmid:36652202
- 46. Hao Y, Hao S, Andersen-Nissen E, Mauck WM 3rd, Zheng S, Butler A, et al. Integrated analysis of multimodal single-cell data. Cell. 2021;184 (13), 3573–3587 e29. pmid:34062119
- 47. van den Brink SC, Sage F, Vertesy A, Spanjaard B, Peterson-Maduro J, Baron CS, et al. Single-cell sequencing reveals dissociation-induced gene expression in tissue subpopulations. Nat Methods. 2017;14 (10), 935–936. pmid:28960196
- 48. Caushi JX, Zhang J, Ji Z, Vaghasia A, Zhang B, Hsiue EH, et al. Transcriptional programs of neoantigen-specific TIL in anti-PD-1-treated lung cancers. Nature. 2021;596 (7870), 126–132. pmid:34290408
- 49. Sun Y, Wu L, Zhong Y, Zhou K, Hou Y, Wang Z, et al. Single-cell landscape of the ecosystem in early-relapse hepatocellular carcinoma. Cell. 2021;184 (2), 404–421 e16. pmid:33357445
- 50. Qiu X, Mao Q, Tang Y, Wang L, Chawla R, Pliner HA, et al. Reversed graph embedding resolves complex single-cell trajectories. Nat Methods. 2017;14 (10), 979–982. pmid:28825705
- 51. Bergen V, Lange M, Peidli S, Wolf FA, Theis FJ. Generalizing RNA velocity to transient cell states through dynamical modeling. Nat Biotechnol. 2020;38 (12), 1408–1414. pmid:32747759