Single-cell analysis reveals cellular reprogramming in advanced colon cancer following FOLFOX-bevacizumab treatment

Introduction The combination of FOLFOX and bevacizumab (FOLFOX-Bev) is a promising treatment for advanced colorectal cancer (CRC). However, the response of the tumor microenvironment to FOLFOX-Bev is still largely unexplored. Methods We conducted single-cell transcriptomic analysis of CRC samples derived from a patient before and after treatment to gain insights into the cellular changes associated with FOLFOX-Bev treatment. Results We found that cancer cells with high proliferative, metastatic, and pro-angiogenic properties respond better to FOLFOX-Bev treatment. Moreover, FOLFOX-Bev enhances CD8+ T cell cytotoxicity, thereby boosting the anti-tumor immune response. Conversely, FOLFOX-Bev impairs the functionality of tumor-associated macrophages, plasma cells, and cancer-associated fibroblasts, leading to a decrease in VEGFB-mediated angiogenesis. Furthermore, FOLFOX-Bev treatment reset intercellular communication, which could potentially affect the function of non-cancer cells. Discussion Our findings provide valuable insights into the molecular mechanisms underlying the response of advanced CRC to FOLFOX-Bev treatment and highlight potential targets for improving the efficacy of this treatment strategy.


Introduction
Colorectal cancer (CRC) is one of the leading causes of cancerrelated deaths worldwide (1), with metastasis being the primary cause of mortality among CRC patients (2). Despite advancements in treatment, the prognosis for advanced CRC remains poor, with a 5-year survival rate of only 12% (3). FOLFOX, which is a combination of three chemotherapeutic agents (oxaliplatin, leucovorin, and fluorouracil), has emerged as a standard therapeutic option for patients with metastatic CRC. Bevacizumab, a monoclonal antibody targeting vascular endothelial growth factor (VEGF), was approved by the FDA as a first-line therapy for metastatic CRC in February 2004 (4). Despite the clinical benefits seen with bevacizumab, not all patients with metastatic CRC respond equally to the treatment, and the underlying mechanism of this variability remains poorly understood.
Recently, it has become increasingly evident that CRC is a highly heterogeneous disease, with the tumor microenvironment (TME) playing a critical role in tumor progression and treatment response. To improve the efficacy of chemotherapy combined with bevacizumab for the treatment of metastatic CRC, it is essential to comprehend the intricate cellular types and their functional characteristics in the TME. Single-cell transcriptomic analyses greatly contribute to our understandings in CRC, including intra-tumor cell heterogeneity (5,6), cancer cell development driven by oncogenic mitogen-activated protein kinase (7), immunosuppressive landscape fostered by oncogenic mutations (8), immune phenotypic linkage between CRC and liver metastasis (9), and alterations of the TME in CRC with chemotherapy (10) or CSF1R blockade (11). However, the response of the TME to chemotherapy combined with bevacizumab is still largely unexplored. Therefore, there is an urgent need to investigate the cellular and molecular mechanisms underlying the response of the TME to FOLFOX-Bev, which could optimize the treatment strategy for metastatic CRC.
The aim of this study was to investigate the impact of FOLFOX-Bev treatment on the TME in advanced CRC using single-cell RNA sequencing (scRNA-seq). Our results revealed important insights into the cellular and molecular changes occurring within the TME following FOLFOX-Bev treatment. Specifically, we observed that cancer cells with high proliferative, metastatic, and pro-angiogenic potential exhibited increased sensitivity to FOLFOX-Bev treatment. Furthermore, we found that FOLFOX-Bev treatment enhanced the cytotoxicity of CD8 + T cells while impairing the functions of tumorassociated macrophages (TAMs), plasma cells, and cancerassociated fibroblasts (CAFs). Additionally, we identified significant inhibition of pathways affecting intercellular communication following FOLFOX-Bev treatment. These findings contribute to our understanding of the TME's response to FOLFOX-Bev treatment in metastatic CRC and have significant implications for the development of personalized treatment strategies.

Tumor samples
This study was conducted with the approval of the institutional review boards of Guangdong Provincial People's Hospital, and informed consent was obtained from the patient with advanced CRC, whose tumor biopsy samples were collected before and after FOLFOX-Bev therapy for scRNA-seq.

Preparation of single-cell suspensions and scRNA-seq
The preparation of single-cell suspensions and scRNA-seq was conducted according to previously published methods (12,13). In brief, fresh tumor biopsy tissue was dissociated, erythrocytes were lysed, and cell viability was tested. The resulting single-cell suspension was then subjected to scRNA-seq using the 10x Genomics Single Cell protocol.

scRNA-seq data preprocessing
The preliminary sequencing results were processed into fastq files using Cell Ranger (v3.0.2). The fastq files were aligned to the human reference genome hg19, and filtering and counting of barcodes and UMI were performed. The Seurat package (version 4.0.4; http://satijalab.org/seurat/) was used to analyze the resulting gene expression matrix. Cells with low quality (gene numbers < 200 or > 6,000, and with a percentage of mitochondrial genes > 50%) were filtered out, and genes detected in less than 10 cells were excluded. The data from the two samples underwent normalization using the default "LogNormalize" method. To address batch effects, we utilized the FindIntegrationAnchors function, employing canonical correlation analysis to identify a set of anchors. These anchors were subsequently used with the IntegrateData function (14) to integrate the objects, thereby correcting for batch effects. Ultimately, a total of 7,343 cells were included for further analysis.

Cell clustering, visualization and annotation
We obtained gene expression data from 7,343 cells and performed principal component analysis to identify principal components. The clustering analysis was conducted using the Louvain algorithm, which was implemented through the "FindClusters" function of the Seurat package. The identified clusters and sub-clusters were visualized using uniform manifold approximation and projection (UMAP) or t-distributed stochastic neighbor embedding (t-SNE) analysis. We annotated cell clusters based on canonical cellular markers from the literature (10,15,16).

Identification of cancer cells
We identified cancer cells using epithelial cell marker genes (EPCAM, KRT19, KRT18) and inferred copy number variations (CNVs) profiles as previously described in our studies (12,13). We performed InferCNV within each sample using immune cells, including T cells, B cells, and myeloid cells, as the reference group and cancer cells as the observation group.

Pathway analysis
To investigate the biological state differences of cancer cells between the sensitive and insensitive subgroups, we performed gene set enrichment analysis (GSEA) based on Hallmark gene sets using the fgsea package. We identified differentially expressed genes upregulated in cancer-associated fibroblast (CAF) metabolism cells using the "FindMarkers" function and screened them based on the condition of p_value < 0.05 and avg_log2FC > 0.5. We then performed Gene Ontology (GO) analysis using the R toolkit ClusterProfiler (17).

Cell cycle analysis
We performed cell cycle analysis of cancer cells using the "CellCycleScoring" function, which scored the probability of S phase and G2M phase of each cell based on marker genes directly related to S phase and G2M phase provided by the Seurat package. The corresponding score indicated the degree of high expression of the S phase (or G2M phase) gene set in this cell. When both S phase score and G2M score were less than 0, the cells were in the G1 phase.

Trajectory analysis
We performed pseudotime trajectory analysis using the Monocle package (v2.20.0) (19)(20)(21). Raw counts from the Seurat RNA assays were used as inputs to create an object for pseudotime trajectory analysis. Data dimensionality was reduced using the "reduceDimension" function with the parameter reduction_method =

Statistics
Statistical analyses were conducted using the "stat_ compare_means" function of the ggpubr package with two-tailed Student's t-test and Wilcoxon Rank Sum Test.

Results
scRNA-seq profiling of TME in advanced colon cancer before and after FOLFOX-Bev treatment To investigate the effects of FOLFOX-Bev treatment on the tumor microenvironment (TME) of advanced colon cancer, we collected biopsy samples from a patient with advanced colon cancer and multiple metastases ( Figure 1A). Before treatment, positron emission tomography (PET) confirmed the diagnosis of sigmoid colon cancer with adjacent bladder wall involvement and multiple metastases to perienteric lymph nodes and liver, as well as multiple small nodules in both lungs ( Figure 1B), indicating high metastatic potential. The patient received four cycles of FOLFOX chemotherapy and four cycles of bevacizumab, resulting in a partial response ( Figure 1C). We obtained biopsy samples from the primary tumor before treatment (Naive) and after treatment (Treat) for further scRNA-seq analysis.
We observed that all cell types were present in both Naive and Treat samples ( Figure 1F), although the proportion of cell clusters varied ( Figure 1G). Notably, the abundance of fibroblasts, and endothelial cells was reduced in the Treat sample, while the infiltration levels of T cells and B cells were relatively increased, suggesting that FOLFOX-Bev treatment may alter the TME of advanced colon cancers.
Cancer cells with strong proliferation, metastasis and pro-angiogenesis are more sensitive to FOLFOX-Bev Next, we classified cancer cells into two subgroups based on their sensitivity to FOLFOX-Bev treatment: a sensitive subgroup (Sensitive) and an insensitive subgroup (Non-sensitive) ( Figure 2A). We then examined the expression of the proangiogenic factor VEGF ( Figure 2B) and found that in the Naive sample, cancer cells mainly secreted VEGFB, which was predominantly expressed by the sensitive subset. However, its expression decreased in the Treat sample, while the insensitive subgroup upregulated the expression of VEGFA, indicating that FOLFOX-Bev treatment could influence the expression of VEGF. Furthermore, we observed that the upregulated genes of the sensitive subgroup compared with the insensitive subgroup were significantly enriched in pathways related to proliferation (mitotic spindle, G2M checkpoint, MYC targets V1, E2F targets, MYC targets V2), metastasis (epithelial-mesenchymal transition), and angiogenesis ( Figure 2C). This suggests that cancer cells with strong proliferation, metastasis, and pro-angiogenesis are more sensitive to FOLFOX-Bev treatment.
We then conducted cell cycle analysis and found that most of the cancer cells from the non-sensitive subgroup and the Treat sample were in the G1 phase, with a minor population in the G2M and S phases, compared to those from the sensitive subgroup and the Naive sample ( Figures 2D, E). This suggests that FOLFOX-Bev treatment could induce G1 cell cycle arrest in cancer cells.
In addition, we observed that the upregulated genes of the non-sensitive subgroup were significantly enriched in energy metabolism pathways, including glycolysis, fatty acid metabolism, and oxidative phosphorylation ( Figure 2C). We scored the energy metabolism signatures of cancer cells ( Figure 2F), and found that glycolysis, fatty acid metabolism, and oxidative phosphorylation were significantly enhanced in the Treat sample, indicating that FOLFOX-Bev treatment could influence the energy metabolism of cancer cells. Overall, our results suggest that FOLFOX-Bev treatment could affect the expression of VEGF, cell cycle progression, and energy metabolism of cancer cells, and that cancer cells with high proliferation, metastasis, and pro-angiogenesis are more sensitive to this treatment.

Impaired function of tumor-associated macrophages (TAMs) after FOLFOX-Bev treatment
We conducted an unsupervised clustering analysis of myeloid cells and identified four distinct clusters, including two tumorassociated macrophages (TAMs) clusters, a neutrophil cluster, and a mast cell cluster ( Figure 3A). TAM M0 cluster only expressed macrophage marker genes (CD14 and CD68) without expressing any polarization genes. In contrast, TAM M1M2 cluster expressed both M1 polarization genes (CD86, IL1B and TNF) and M2 polarization genes (CD163, MS4A4A and TGFB1). Neutrophils highly expressed neutrophil-specific genes (S100A8, S100A9, and G0S2), while mast cells highly expressed mast cell-specific genes (KIT, GATA2, and TPSAB1) ( Figure 3B). We observed that the infiltration level of TAMs, particularly the TAM M0 cluster, was significantly increased in the Treat samples ( Figure 3C). To further examine the effect of FOLFOX-Bev treatment on the functional characteristics of TAMs, we calculated the M1/M2 polarization and pro-/anti-inflammatory scores of TAMs. Our analysis showed that the M1/M2 polarization and pro-inflammatory scores of TAMs decreased markedly after FOLFOX-Bev treatment ( Figure 3D). Moreover, the angiogenic activity of TAMs was also suppressed by FOLFOX-Bev treatment, mainly indicated by a reduction in the expression of VEGFB (Figures 3E, F). Thus, TAMs appeared to be in a relatively quiescent state after FOLFOX-Bev treatment, with diminished functional activity.
We then evaluated the distribution of T cell subpopulations and observed increased infiltration of CD8 + T cells and decreased infiltration of immunosuppressive Treg FOXP3 cells in the Treat sample ( Figure 4C). CD8 + T cells are crucial for anti-tumor immunity, and therefore, we further examined the effect of FOLFOX-Bev treatment on their functional signatures. We extracted CD8 + T cells for further analysis and found that the activation and cytotoxic scores of CD8 + T cells were significantly elevated in the Treat sample compared to the Naïve sample ( Figure 4D). These findings suggest that FOLFOX-Bev treatment might enhance anti-tumor immunity mediated by CD8 + T cells.

Decreased angiogenic function of plasma cells after FOLFOX-Bev treatment
The B cells were re-clustered into four distinct subsets ( Figure 5A), with B Naive cells exhibiting high expression of CD19, MS4A1 and BANK1. Two clusters of plasma cells were identified based on their high expression levels of MZB1, IGLL5 and IGJ, one of these clusters exhibited expression of VEGFA and VEGFB, thus called plasma cell VEGF. Cycling cells were identified as a proliferative cluster with high expression of MKI67, TOP2A and PCNA ( Figures 5B, S3A, B). Interestingly, we observed increased infiltration of B Naive cells and reduced infiltration of plasma cells that release VEGF factors in the Treat sample, indicating a potential reduction in the angiogenic function of plasma cells mediated by VEGF factors following FOLFOX-Bev treatment ( Figure 5C).
To further understand the dynamic cell transitions of B cell subsets, we used monocle for unsupervised analysis to construct the trajectory. We observed that B Naive cells and the proliferative cluster (Cycling cells) were located at the beginning of the trajectory, whereas the two plasma cell clusters were more concentrated at the end of the trajectory ( Figure 5D). Moreover, three sets of differentially expressed genes were identified along the trajectory. The first set, consisting of naive B cell markers (CD22, CD19, MS4A1 and BANK1), decreased along the trajectory, while the second set, consisting of proliferative genes (MKI67, TOP2A, and PCNA), were highly expressed in the middle of the trajectory. The last set, consisting of plasma cell markers (IGLL5, IGJ and SDC1), increased towards the end of the trajectory ( Figure 5E). These analyses reveal the dynamic cell transitions of B cells, from naive to terminally differentiated plasma cells, providing insight into the potential role of B cells in the immune response to cancer and the effects of FOLFOX-Bev treatment on this response.

Impaired secretory function of cancerassociated fibroblasts (CAFs) after FOLFOX-Bev treatment
Next, we re-clustered cancer-associated fibroblasts (CAFs) into four distinct subtypes based on their gene expression profiles ( Figure 6A). CAF ECM cells exhibited high expression of extracellular matrix (ECM) molecules such as COL1A1, COL3A1,   Figures 6B, S4A). CAF metabolism cells were significantly enriched in pathways related to RNA metabolism ( Figure S4B). Unexpectedly, our analysis revealed a significant decrease in the number of CAF secretory cells in the Treat sample ( Figure 6C), indicating that FOLFOX-Bev treatment might have inhibited the secretory function of CAFs, including the reduction of VEGFB, CXCL1 and CCL13 secretion ( Figure 6D), which could potentially impair the angiogenic function of CAFs.

Altered intercellular communication in advanced colon cancer after FOLFOX-Bev treatment
To further investigate the potential effects of FOLFOX-Bev treatment on intercellular communication within the TME, we utilized CellChat for cell communication analysis. Notably, our analysis revealed significant inhibition of the VEGF pathway after treatment ( Figure 7A), particularly the VEGFB-VEGFR1 pathway which was the primary VEGF pathway in the Naive sample and was markedly suppressed in the Treat sample. Conversely, we observed transcriptional upregulation of the VEGFA-mediated pathway in the Treat sample ( Figures 7B,  S5A), indicating that FOLFOX-Bev treatment might have a complex impact on the VEGF pathway.
Further analysis of the TME pathways revealed that PERIOSTIN, ncWNT, TWEAK, CSF, ANGPT, and APELIN pathways were exclusively present in the Naive sample ( Figure  S5B). Moreover, our analysis of the distribution of these pathways among various cell populations ( Figure 7C) showed that PERIOSTIN, ncWNT, and TWEAK pathways primarily targeted fibroblasts and were involved in promoting fibroblast activation, proliferation, and inducing tissue fibrosis and metastatic niche formation (23)(24)(25). The CSF pathway primarily targeted TAMs and may play a role in controlling their proliferation, differentiation, survival, and function (26,27). Additionally, we found that ncWNT, ANGPT, and APELIN pathways could target endothelial cells, regulating endothelial differentiation and angiogenesis (28)(29)(30)(31)(32).
In conclusion, our findings suggest that FOLFOX-Bev treatment not only targets cancer cells but also affects intercellular communication in the TME, leading to inhibition of pathways related to the survival, proliferation, differentiation, and function of non-cancer cells.
To address potential limitations in our study, including small sample size, single-center design, and patient selection bias, we sought to validate our findings by performing additional analyses using an independent publicly available dataset (GSE178318). In this dataset, we focused on untreated (Naïve, n=3) and FOLFOX-Bev-treated (Treat, n=1) colon cancer orthotopic samples. Initially, we annotated the cell types in this dataset using the top differentially expressed genes (DEGs) of the respective cell clusters (Figures S6A,  B). Subsequently, we classified cancer cells into two subgroups based on their sensitivity to FOLFOX-Bev treatment: a sensitive subgroup (Sensitive) and an insensitive subgroup (Non-sensitive) ( Figure S6C). Consistent with our original study, we observed that the upregulated genes in the sensitive subgroup, compared to the insensitive subgroup, were significantly enriched in pathways related to proliferation (E2F targets, G2M checkpoint, mitotic spindle) and metastasis (epithelial-mesenchymal transition) ( Figure S6D). This finding supports the notion that cancer cells exhibiting higher proliferation and metastatic potential are more responsive to FOLFOX-Bev treatment.
Furthermore, we conducted cell cycle analysis and found that FOLFOX-Bev treatment induced G1 cell cycle arrest in cancer cells, as evidenced by an increased proportion of cells in the G1 phase ( Figures S6E, F). This observation aligns with our original study and strengthens the evidence for the ability of FOLFOX-Bev treatment to influence cell cycle progression in cancer cells. Regarding the tumor-associated macrophages (TAMs), our analysis confirmed that TAMs displayed a relatively quiescent state with reduced functional activity following FOLFOX-Bev treatment, as indicated by scoring the functional characteristics of TAMs ( Figures S6G, H). This validation further supports our original finding of TAM impairment upon treatment. Additionally, we observed a decrease in the number of cancer-associated fibroblast (CAF) secretory cells and a reduction in the expression of vascular endothelial growth factors (VEGFs) in the Treat sample, suggesting that FOLFOX-Bev treatment may hinder the secretory function of CAFs and potentially impede angiogenesis ( Figures S6I-K). These results align with our original study and provide further evidence of the impact of FOLFOX-Bev treatment on CAF-mediated angiogenesis

Discussion
Despite advances in the treatment of advanced CRC, there is limited knowledge about the changes in the TME that occur in response to chemotherapy combined with bevacizumab. In this study, we used scRNA-seq to characterize the cellular landscape of advanced CRC before and after FOLFOX-Bev therapy, encompassing cancer cells, immune cells, and stromal cells. Our findings shed light on the response of specific cell subsets to FOLFOX-Bev therapy, highlighting the potential of this treatment regimen in treating advanced CRC. Importantly, we identified the cancer cell subset with high proliferative, metastatic, and proangiogenic properties showed sensitivity to FOLFOX-Bev therapy. Furthermore, we demonstrated that FOLFOX-Bev treatment impaired VEGFB-mediated angiogenesis. These findings have significant implications for the clinical management of advanced CRC and offer insight into the potential use of predictive biomarkers to improve patient outcomes.
Bevacizumab has demonstrated significant preclinical and clinical benefits against metastatic CRC, especially when combined with chemotherapy (33). However, clinical trials investigating the efficacy of adding bevacizumab to FOLFOX in early-stage colon cancer have produced conflicting results. While Allegra et al. reported that the addition of bevacizumab did not significantly prolong disease-free survival in stages II and III colon cancer, despite significant but transient effects being observed (34). Gramont et al. found no benefit for disease-free survival and even a potential detrimental effect for overall survival when bevacizumab was added to adjuvant chemotherapy in stage III colon cancer (35). The reason for the different responses to bevacizumab in metastatic and early-stage CRC remains unclear. In this study, we used scRNA-seq to identify the functional characteristics of FOLFOX-Bev-sensitive cancer cells. We found that these cells exhibited enhanced epithelial-mesenchymal transition and proangiogenic capacity in addition to strong proliferation, which are key factors for the metastasis of CRC. Our findings suggest that CRC with a strong metastatic tendency may be more sensitive to FOLFOX-Bev, which could partly explain why FOLFOX-Bev is more effective in treating metastatic CRC.
The VEGF family plays a critical role in regulating both normal and pathological angiogenesis, which supports the growth, invasion, and metastasis of CRC. The expression of the predominant VEGF subtypes is correlated with the progression of CRC (36,37). For example, high expression of VEGFB and VEGFC has been found to significantly correlate with lymph node metastasis and lymphatic invasion, while VEGFA is higher in CRC with liver metastasis than in those without (Kawakami et al., 2016) (36). In this study, we observed that VEGFB, secreted by both cancer cells and non-cancer cells, was the predominant VEGF subtype in advanced CRC with multiple metastases to the liver, lung, and perienteric lymph nodes. Bevacizumab is a humanized monoclonal antibody that binds to the VEGFA subtype. Surprisingly, we found that FOLFOX-Bev treatment significantly inhibited the VEGFB-VEGFR1 pathway while enhancing the VEGFA pathway based on transcriptional level. We propose several hypotheses to explain this result. Firstly, FOLFOX-Bev treatment may kill sensitive cancer cells that primarily secrete VEGFB, leading to compensatory upregulation of VEGFA in the remaining insensitive cancer cells, which is competitively bound by bevacizumab, thereby dually inhibiting the VEGF pathway. Additionally, we found that M2 polarization of TAMs was significantly inhibited and remained in a relatively quiescent state after FOLFOX-Bev treatment, resulting in the inhibition of their proangiogenic function. Furthermore, we observed that CAF secretory cells, which mainly secrete VEGF, were significantly reduced after FOLFOX-Bev treatment, potentially impairing their angiogenic function. Lastly, we identified a distinct subpopulation of VEGFsecreting plasma cells that were significantly reduced after FOLFOX-Bev treatment. Although the primary function of plasma cells is to synthesize, store, and secrete antibodies, our study and previous research (38) have shown that they can also produce VEGF.
Interestingly, we made an intriguing observation that Treg FOXP3 cells exhibited high expression of cytotoxic T-lymphocyte protein 4 (CTLA-4), a potent molecule known to inhibit T cell activation. CTLA-4 has emerged as a critical immune checkpoint receptor with implications for immune evasion by cancer cells in various cancer types, including CRC (39). Studies by Afshin et al. demonstrated an increasing trend of CTLA-4 expression in CRC tissues and cell lines, and treatment with capecitabine effectively suppressed the expression of CTLA-4 in CRC cells, suggesting a potential bridge between immunotherapy approaches and chemotherapy (40). In addition, multiple research groups have suggested a significant role for CTLA-4 in CRC progression (41,42), immunosuppression (39), and therapeutic interventions (43). These findings highlight the potential of CTLA-4 as an important target for colon cancer treatment.
There are several important limitations that should be taken into account: i) the findings were based on biopsy samples obtained from a limited number of advanced colon cancer patients before and after FOLFOX-Bev therapy. Although we have obtained some further confirmation through analysis of a published dataset, it is essential to conduct future research with larger and more diverse cohorts to enhance the generalizability of the findings. ii) In terms of cell annotation, we acknowledge the limitations associated with relying solely on canonical cellular markers, as they may not fully capture the complexity and heterogeneity of cell populations. Therefore, we recognize the need for complementary approaches, such as single-cell proteomics or spatial transcriptomics, to provide a more comprehensive and accurate cell characterization. iii) the cell annotation in this study was based on canonical cellular markers from the literature. While these markers provide a useful framework for cell classification, they may not encompass the full complexity and diversity of the cell populations under investigation. Therefore, future studies should aim to incorporate additional markers and employ more comprehensive approaches to capture the full range of cellular heterogeneity. iv) in examining functional differences between cell populations, we utilized predefined functional gene sets derived from previous published papers. However, it is crucial to acknowledge that these gene sets may not encompass the full complexity and intricacies of the biological processes under investigation. Future studies should consider incorporating a broader range of functional gene sets or employ more comprehensive and unbiased methods to explore the functional landscape of the cell populations.
In conclusion, our findings highlight the complex impact of FOLFOX-Bev treatment on the VEGF pathway in advanced CRC. By comprehensively characterizing the cellular composition and transcriptional profiles of the TME in advanced CRC with FOLFOX-Bev therapy, our study provides valuable insights that could aid in the identification of potentially advantageous populations. Moving forward, it is crucial to further investigate and compare the transcriptional profiles among patients exhibiting different responses to FOLFOX-Bev treatment. This approach will enable the identification of clear biomarkers that can guide the selection of patients who are more likely to benefit from FOLFOX-Bev therapy.

Data availability statement
The data presented in the study are deposited in the NCBI Gene Expression Omnibus database (https://www.ncbi.nlm.nih.gov/), accession number GSE232525 and GSE178318.

Ethics statement
The studies involving human participants were reviewed and approved by the institutional review boards of Guangdong Provincial People's Hospital. The patients/participants provided their written informed consent to participate in this study. Written informed consent was obtained from the individual(s) for the publication of any potentially identifiable images or data included in this article.

Author contributions
MY performed the scRNA-seq analyses, wrote the manuscript and revised the manuscript. CY, DM and ZL interpreted the data and revised the manuscript. WZ and DY interpreted the data, wrote the manuscript, and provided supervision.