Single-cell and spatial transcriptome profiling reveal CTHRC1+ fibroblasts promote EMT through WNT5A signaling in colorectal cancer

Abstract

Background

Cancer-associated fibroblasts (CAFs), known for facilitating the progression and metastasis of colorectal cancer (CRC), have become a promising therapeutic target. However, the significant heterogeneity of CAFs and their intricate crosstalk with tumor cells present substantial challenges in the development of precise and effective therapeutic strategies.

Methods

Single-cell RNA sequencing (scRNA-seq) technology was used to identify various cell subtypes. Spatial transcriptomics (ST) was employed to map the spatial niches and colocalization patterns of these cell subtypes. Cell-cell interactions among these subtypes were analysed via CellChat and NicheNet software. Tumor cell invasion, migration, and proliferation were assessed through wound healing assays, transwell assays, colony formation assays, and xenograft mouse models.

Results

We identified a significant spatial colocalization between CTHRC1+ CAFs and a distinct subtype of malignant epithelial cells, both residing within the EMT-active spatial niche. Our results demonstrate that CTHRC1+ CAFs, as a major source of WNT5A, promote epithelial-mesenchymal transition (EMT) and enhance tumor cell invasiveness by upregulating MSLN expression in adjacent malignant epithelial cells. This signaling axis contributes significantly to CRC progression and metastasis.

Conclusions

Targeting the CTHRC1+ CAF-WNT5A-MSLN signaling axis presents a promising therapeutic strategy for advanced CRC patients. Our study provides new insights into the role of CAFs in CRC progression and offers potential avenues for developing targeted therapies to disrupt this pathway.

Background

Colorectal cancer (CRC) stands as the third most frequently diagnosed cancer worldwide, accounting for approximately 10% of all cancer cases [1]. The incidence and mortality rates have consistently been high in recent years. The five-year survival rate for patients with metastatic CRC is less than 20%, making it a leading cause of cancer-related deaths globally [2]. Nowadays, systemic therapy is the primary approach for metastatic CRC. However, owing to individual heterogeneity among CRC patients, systemic treatment regimens vary from person to person, resulting in significant differences in treatment efficacy [3]. The tumor microenvironment (TME) significantly embodies the heterogeneity among CRC patients. Currently, the cells within the TME and the cytokines they secrete are recognized as crucial in cancer progression and present promising therapeutic targets [4].

Cancer-associated fibroblasts (CAFs) are pivotal components of the TME in solid tumors, exhibiting significant heterogeneity and playing a pivotal role in modulating the behavior of malignant epithelial cells (EPCs) [5]. Given the diverse roles of CAFs in tumor cell dynamics, targeting the interaction between CAFs and tumor cells has garnered increasing attention as a promising strategy to impede CRC progression and improve patient prognosis [6]. Understanding the interaction between CAFs and malignant EPCs is vital for the development of therapies designed to disrupt this crosstalk.

Before the emergence of single-cell RNAsequencing (scRNA-seq) technology, the substantial heterogeneity of CAFs made it difficult to explore their complexity and plasticity in an unbiased manner [7]. The emergence of scRNA-seq has made this exploration possible. However, a major limitation of scRNA-seq is the absence of spatial details [8]. Recently developed spatial transcriptomics (ST) technology overcomes this limitation by integrating spatial gene expression with histological information. Therefore, combining scRNA-seq and ST technologies can provide a spatially resolved gene expression map of CAFs and tumor cells at single-cell resolution [9].

This study primarily aims to leverage multi-omics strategies to investigate the interaction between CAFs and metastatic EPCs in CRC. Utilizing scRNA-seq technology, we examined the function of a previously under-studied population of CTHRC1+ fibroblasts in CRC. Our results indicate that tumor-specific CTHRC1+ fibroblasts can remodel the stroma within the TME, mediate the exclusion of CD8+ T cells, and induce epithelial-mesenchymal transition (EMT) in tumor EPCs. Spatial transcriptomics analysis further identified colocalization of CTHRC1+ fibroblasts with specific malignant EPCs, with strong crosstalk promoting EMT and tumor cell metastasis. We hope to identify novel therapeutic targets on the basis of these cellular interactions to optimize the efficacy of targeted therapies for CRC patients.

Methods

Data sources

The single-cell data used in our study were derived from the GSE132465 and GSE144735 single-cell datasets, which included 16 normal samples, 6 tumor border samples, and 29 tumor samples [10]. Detailed clinical information can be found in Table S1. Two spatial transcriptomic slice data (10X Genomics) were obtained from the Cancer Diversity Asia database [11]. The RNA-seq data and corresponding clinical characteristics from TCGA-COAD and TCGA-READ, which included 584 CRC tissues, were retrieved from UCSC Xena.

Single-cell data preparation

This study employed the R package Seurat (v4.4.0) [12] to analyse scRNA-seq data from CRC samples in the GSE132465 and GSE144735 datasets. Initially, the datasets were integrated, and low-quality cells were excluded using criteria of fewer than 200 or more than 6000 genes expressed per cell, or a mitochondrial gene proportion exceeding 20%. After quality control, the data were normalized and log-transformed, and the top 3000 highly variable genes were selected for subsequent analysis. Principal component analysis was then conducted via the “RunPCA” function to reduce data dimensionality, and batch effects between samples were corrected using the Harmony algorithm [13]. To determine the cell types, the cells were annotated with reference to the CellMarker 2.0 database [14] and related literature.

InferCNV analysis to identify epithelial cell states

To distinguish between non-malignant and malignant EPCs, we employed the R package infercnv (v1.16.0) [15] to analyse copy number variation (CNV) in scRNA-seq data from epithelial cells. In this process, CNV scores from epithelial cells in normal samples served as spiked-in controls. We conducted hierarchical clustering analysis on the CNV profiles of epithelial cells from both normal and tumor samples. The cells in clusters that are predominantly composed of epithelial cells from normal samples are considered non-malignant cells. In contrast, cells in clusters in which almost all exhibit chromosomal deletions and amplifications are identified as malignant cells [16].

Spatial transcriptomics data analysis

We employed the R package Seurat to analyse the spatial transcriptomics data from the spatial spot expression matrix filtered by SpaceRanger. Initially, we performed quality control on the ST data by retaining genes expressed in at least three cells and filtering out spatial spots with fewer than 300 measured genes and fewer than 600 UMIs. Each individual slice matrix was normalized and standardized using “SCTransform”, setting the number of variable features to 3,000. Cell2location (v0.1.3) [17] was run in a python (v3.9) environment to determine the cell type composition of each spot. For each spatial spot, we calculated the proportions of various cell types via cell type-specific abundance estimates.

Spatial niche definition

The spatial niche refers to regions adjacent in space with similar compositions of cell types, reflecting local microenvironments within tissues [18]. To identify clusters of spots characterized by similar compositions of cell types across different samples, we converted the estimated cell type proportions for each spatial transcriptomic spot and slice into isometric log-ratio coordinates and clustered the spots into groups. We created shared nearest neighbor graphs with various neighbor counts [10, 20, 30] using the “buildSNNGraph” function from the R package scran (v1.28.2) [19] and performed Louvain clustering. The clustering resolution was then optimized to maximize the average silhouette score. Finally, the median proportion of each cell type within each niche was used to quantify the cellular abundance in each niche.

MISTy for capturing the spatial map of cellular dependencies

We used the R package mistyR (v1.10.0) [20] to assess the importance of different groups of malignant epithelial cell types in explaining the abundance of each group of fibroblast types. The spatial cell type abundances predicted by cell2location are modeled using three different spatial contexts in a multiview framework: (1) the intrinsic view, which measures the relationship between deconvolution estimates within a single point; (2) the juxta view, which aggregates deconvolution estimates from neighboring cells within a maximum distance of 5 points; and (3) the para view, which applies weighted deconvolution estimates for distant neighbors within an effective radius of 15 points. By constructing and analysing these views, we can evaluate the spatial dependencies between malignant epithelial cells and fibroblasts across different spatial environments. The importance of each view in spatial slices is interpreted in terms of cell type dependencies in various spatial contexts, including colocalization or mutual exclusion.

Cell-cell interaction analysis

We utilized CellChat (v2.1.2) [21] and CellPhoneDB (v5.0.0) [22] to calculate the intensity of cell-cell interactions among major cell types in spatial niches. To enhance the reliability and robustness of the CellChat analysis, we filtered out communication interactions involving fewer than 10 cells. The “netVisual” function was employed to visualize ligand-receptor-mediated cell interactions. We subsequently used the R packages liana (v0.1.13) [23] and nichenetr (v2.0.5) [24] to determine interactions between ligands, receptors, and target genes across two target cell types. The top 50 ligands predicted by LIANA were used as potential ligands from “sender cells” in NicheNet analysis, integrating 80 differentially expressed target genes from “receiver cells” for ligand, receptor, and target gene activity analysis.

Tissue sample collection and cell culture

All samples were collected from the First Affiliated Hospital of Nanjing Medical University, with approval from the Medical Ethics Committee and informed consent obtained from all participating individuals. The Chinese Academy of Sciences’ Cell Bank of Type Culture Collection in Shanghai, China, provided the cell lines NCM460, SW480, RKO, DLD-1, LoVo, SW620, and HCT116. All the cell lines were cultured at 37 °C in a humidified atmosphere containing 5% CO2. The cells were validated by STR analysis in 2022 and were periodically screened for mycoplasma contamination via PCR.

Multiplex immunohistochemistry staining (mIHC)

Tissue sections embedded in paraffin were first deparaffinized using xylene and then rehydrated with graded alcohols. Antigen retrieval was performed in an EDTA buffer, followed by the blocking of endogenous peroxidase activity with hydrogen peroxide. To minimize non-specific binding, 3% BSA was applied as a blocking agent. The sections were then incubated sequentially with primary antibodies (anti-KLK10, mouse, 1:200, Santa Cruz, Cat# sc-100551; anti-CTHRC1, rabbit, 1:50, Proteintech, Cat# 16534-1-AP; anti-EPCAM, rabbit, 1:500, Proteintech, Cat# 21050-1-AP) and species-specific secondary antibodies. Tyramide signal amplification was conducted after each incubation, with antigen retrieval and antibody stripping performed between rounds. Finally, DAPI was applied to counterstain the nuclei, and the sections were mounted using antifade mounting medium.

Immunohistochemistry (IHC) analysis

We collected tumor and adjacent normal tissue samples from 80 patient pairs after written informed consent was obtained from the CRC patients. Detailed clinical information is presented in Table S1. The IHC staining was performed by experienced pathologists. To quantify MSLN expression levels, we calculated the H-Score using the IHC Profiler plugin in ImageJ.

Isolation, and culture of CAFs from human primary colorectal cancer tissues

The human colorectal cancer tissues were first washed and then incubated in PBS containing 3% penicillin/streptomycin/amphozone. The tissue was then cut into small pieces approximately 1 mm³ in size using ophthalmic scissors. These tissue fragments were subjected to enzymatic digestion by repeatedly pipetting with collagenase IV (Sigma, USA), hyaluronidase (BioFroxx, Germany), and DNase I (BioFroxx, Germany). The digestion process was carried out at 37 °C with shaking for 60 min. After digestion, the resulting mixture was filtered through a 40 μm mesh (BD Falcon, USA) to remove undigested debris, followed by centrifugation. The digested tissue was then resuspended in complete DMEM (Gibco, USA) containing 20% fetal bovine serum (FBS, Gibco, USA), and incubated at 37 °C with 5% CO2 for subsequent expansion. A specific group of CAFs in good condition was selected for further passaging and subsequent experiments. CAFs are verified by immunofluorescence staining of α-SMA.

Co-culture of CAF and CRC cells and Preparation of Conditioned Medium

CAFs (1 × 10⁵ cells) were seeded into the lower chamber of 0.4 μm pore size transwell chambers designed for 6-well plates (Corning, USA), while CRC cells (1 × 10⁵ cells) were seeded into the upper chamber. After a 48-h incubation, the transwell chambers were then transferred to new 6-well plates, washed, and tumor cells were harvested. The samples were then treated with lysis buffer and stored at −80 °C for further analysis. To prepare conditioned media, the culture medium was replaced with fresh, serum-free DMEM once the CAFs reached the desired confluence. The cells were cultured for an additional 48 h, after which the CAF-conditioned supernatants were collected. These supernatants were centrifuged to remove any remaining cell debris and stored at −80 °C until needed for subsequent experiments.

Enzymelinked immunosorbent assay (ELISA)

The WNT5A levels in the supernatant of CAFs were measured using human WNT5A ELISA kits (#MM-50715H2, MEIMIAN, China) according to the manufacturer’s instructions.

RNA extraction and quantitative real‑time PCR (qRT‑PCR)

The detailed process was performed in line with previous research [25]. The primers used in this study are listed in Table S2.

Cell transfection

The lentiviral vectors carrying shRNA targeting MSLN and CTHRC1 were produced by Tsingke (Beijing, China). The MSLN and CTHRC1 overexpression plasmids was obtained from Obio (Shanghai, China). Transfection of shRNAs and plasmids was performed with Lipo3000 (Invitrogen, Carlsbad, USA) according to the manufacturer’s protocol. The efficiency of transfection was validated using qRT-PCR and Western blotting assay. The sequences of the shRNAs are provided in Table S2.

Western blot assay

The cells were lysed in RIPA buffer (NCM Biotech, China) containing protease inhibitors (NCM Biotech, China). Protein concentrations were quantified using a BCA protein assay kit (Beyotime, China). The primary antibodies used in this study included CTHRC1 (1:2000, Proteintech Cat# 16534-1-AP), WNT5A (1:1000, Proteintech Cat# 55184-1-AP), MSLN (1:1000, Proteintech Cat# 66404-1-Ig), E-cadherin (1:2000, Proteintech Cat# 60335-1-Ig), N-cadherin (1:2000, Proteintech Cat# 22018-1-AP), Vimentin (1:1000, Proteintech Cat# 22031-1-AP), and GAPDH (1:50,000, Proteintech Cat# 60004-1-Ig).Western blot analysis was conducted as described previously [26].

Cell proliferation assay

Cell proliferation was assessed using the Cell Counting Kit-8 (Beyotime, China) and colony formation assays, according to previously established methods [25].

Wound healing assay

The treated cells were plated in six-well plates at a density of 5×$\:{10}^{5}$ cells per well and cultured until they reached 90% confluence. A scratch was made in the center of each well using a sterile 200 µL pipette tip. Following the scratch procedure, the corresponding CAF-conditioned medium was added, and the culture was continued. Cell migration was assessed by microscopy at 0 and 48 h.

Transwell assay

To evaluate the ability to migrate and invade, transwell assays were conducted using transwell chambers (8 μm pore size; Corning) with and without Matrigel. In the upper chamber, a total of 2×$\:{10}^{4}\:$CRC cells were resuspended in 200 µL of serum-free medium. The corresponding CAFs (4×$\:{10}^{4}$/well) were plated in the lower chamber. Due to experimental needs in some wells, the lower chamber was supplemented with 350 ng/mL recombinant human WNT5A (R&D Systems, 645-WN, USA) and 60 µmol/L Box5 (MCE, HY-123071, Shanghai). Cells were incubated at 37 °C for 36 h. After incubation, the migrated cells were fixed, stained, and counted in three randomly selected fields of view to evaluate migration and invasion.

Animal model

MSLN shRNA and overexpression plasmids were separately transfected into HCT116 and DLD-1 cells in vitro. Subsequently, 1×$\:{10}^{6}\:$treated cells were injected under the skin of five-week-old male BALB/c nude mice. Tumor volumes were monitored and measured every 5 days. The mice were euthanized 25 days post-injection, and the tumors were removed for final volume and weight assessments.

Statistical analysis

Statistical analysis was performed using R version 4.4.1, Python 3.9, and GraphPad Prism 9.0. Differences between the two groups were evaluated using either Student’s t-test or the Wilcoxon rank-sum test. Survival data were analyzed using the log-rank test. All experiments were conducted in triplicate or more. Data are presented as mean ± standard deviation (SD). A p-value of less than 0.05 was considered statistically significant for all statistical comparisons.

Results

Single-cell transcriptomic profiling of CRC patients

The graphic overview of the study design is shown in Fig. 1A. To elucidate the heterogeneity of the TME in CRC patients, we integrated the single-cell datasets GSE132465 and GSE144735. These datasets encompass tumor core tissues, tumor edge tissues, and normal tissues from CRC patients [10]. The integration yielded a high-quality expression matrix with 91,103 cells for further analysis. Based on cell markers from the CellMarker 2.0 database, we classified the 91,103 single cells into seven major types (Fig. 1B), including 24,708 epithelial cells marked by EPCAM and KRT8; 13,586 stromal cells marked by COL1A1 and DCN; 9,465 myeloid cells marked by S100A9 and LYZ; 28,874 T cells marked by CD3D and CD3E; 5,609 B cells marked by CD79A and MS4A1; 8,417 plasma cells marked by CD79A and MZB1; and 444 mast cells marked by KIT and CPA3 (Fig. 1C, D). Furthermore, compared to normal tissues, we observed a significant elevation in the proportion of epithelial and myeloid cells, coupled with a marked reduction in stromal and plasma cells, both in the core and edge regions of the tumors. (Fig. 1E). This may indicate the remodeling of the TME to some extent.

Fig. 1

Single-cell transcriptomic profiling of CRC patients. (A) Flowchart of the overall research design. (B) t-SNE plots of scRNA-seq data for 91,481 cells in this study. (C) Stacked-violin plots showing average expression of known markers in indicated cell types. (D) Feature plots illustrating the expression of canonical marker genes for each cell type. (E) Bar plots showing the proportions of the 7 major cell types across different tissue types, including adjacent normal tissues, tumor border, and tumor core

TME in CRC shapes a tumor-specific fibroblast population

Fibroblasts are key components of the TME and play significant roles in promoting tumor growth, invasion, and metastasis through various mechanisms [27]. Investigating the role of specific fibroblasts in the TME is crucial for understanding cancer progression and developing new therapeutic strategies. To study fibroblast heterogeneity, we performed subcluster analysis on 13,586 stromal cells and identified 11 subtypes, including 7 types of fibroblasts (Fig. 2A). Using tSNE dimensionality reduction based on tissue origin, we found that fibroblasts with high expression of the CTHRC1 gene (Fig. S1A, B) existed only in tumor tissues, not in normal tissues (Fig. 2B). The proportion of CTHRC1+ fibroblasts was relatively low in early-stage (stage I-II) patients but significantly increased in advanced stages (Fig. 2C). By integrating clinical prognosis data from CRC patients with the scissor algorithm [28], we found that CTHRC1+ fibroblasts aggregated a large number of cells associated with poor patient prognosis (Fig. 2D). Utilizing the ssGSEA algorithm, we quantified the infiltration status of CTHRC1+ fibroblasts based on the expression of the top 30 marker genes, classifying TCGA CRC patients into high and low CTHRC1+ fibroblast infiltration groups. Patients with a high infiltration of CTHRC1+ fibroblasts were found to have an unfavorable prognosis (Fig. 2E). Differential analysis of the 11 stromal cell subtypes revealed that CTHRC1+ fibroblasts were particularly enriched with genes related to extracellular matrix remodelling, such as POSTN, MMP11, and TGFB1 (Fig. 2F). GSEA functional enrichment analysis indicated that CTHRC1+ fibroblasts were more active in promoting EMT and angiogenesis compared to other fibroblasts (Fig. 2G, Fig. S1C). We then employed the R package monocle (v2.28.0) to analyze the polarization trajectory of the fibroblasts. Pseudotime analysis revealed that proliferative fibroblasts and CTHRC1+ fibroblasts were mainly concentrated in the early stages of fibroblast differentiation. One branch differentiated into F3+ fibroblasts with myofibroblast properties, whereas the other branch differentiated into DPT+ fibroblasts (Fig. 2H). To verify these findings, we used the python package CytoTRACE (v0.3.3), another tool for inferring cell developmental trajectories. The results demonstrated that CTHRC1+ fibroblasts exhibit increased stemness and developmental potential (Fig. S1D, E). These findings indicate that CTHRC1+ fibroblasts might contribute to extracellular matrix remodelling and are closely associated with EMT in tumor cells, predicting poor prognosis for CRC patients.

Fig. 2

Tumor-specific fibroblast population shaped by TME in CRC. (A)t-SNE visualization of stromal cell composition, with colors indicating different cell types. (B) t-SNE plots of stromal cells faceted by tissue types. (C) Stacked bar plots showing the percentage of fibroblast subtypes across different CRC stages. (D) t-SNE plot showing cells selected by Scissor. Red dots indicate Scissor+ (worse survival) cells, while gray dots represent background cells. (E) Survival curve of overall survival for CTHRC1+ fibroblasts in TCGA datasets. (F) Heatmap of characteristic gene expression and pathway enrichment analysis for stromal cell subtypes. (G) GSEA of EMT pathway between CTHRC+ fibroblasts and other fibroblasts. (H) The arrangement of fibroblasts subtypes along pseudotimes within the two-dimensional state space as determined by Monocle2

Identification of a malignant epithelial cell subtype associated with EMT

We used tSNE dimensionality reduction to visualize the distribution of normal tumor EPCs, tumor margin EPCs, and tumor core EPCs and observed clear boundaries between tumor EPCs and normal EPCs (Fig. 3A). To differentiate malignant cells from non-malignant cells within epithelial subtypes, we first performed inferCNV analysis on each epithelial cell to infer large-scale chromosomal CNV and provided each cell with a CNV score. Normal sample epithelial cells presented low CNV scores (Fig. 3B). Hierarchical clustering was subsequently conducted on the CNV profiles of epithelial cells from both normal and tumor samples. This analysis resulted in the identification of 15,064 malignant and 9,644 normal epithelial cells. Compared with normal epithelial cells, malignant epithelial cells had significantly higher CNV scores (Fig. 3C). CNV scores and states for all epithelial cells are listed in Table S3. The identified malignant EPCs were further divided into 8 subtypes (Fig. 3D). Among the malignant EPC subtypes, malignant EPC_1 significantly changed with advancing CRC stage (Fig. 3E) and was associated with tumor progression, distant metastasis, and MSI status (Fig. S2A), suggesting a correlation with poor prognosis. We found that the genes KLK10 and CDA, which are specifically expressed by malignant EPC_1, were associated with poor prognosis in CRC patients (Fig. 3F, Fig. S2B, C). To understand the roles of each malignant EPC subtype in the TME, we scored each subtype via hallmark gene sets. The malignant EPC_1 demonstrated elevated activity in both EMT and angiogenesis (Fig. 3G). The malignant EPC_6 also had high EMT scores but was associated with better patient prognosis (Fig. S2D). Marker genes for each malignant epithelial subtype are listed in Table S4. We identified the malignant EPC_6 as goblet cells, and their elevated EMT score could be linked to chronic inflammation and hypoxic conditions. Moreover, pySCENIC analysis [29] revealed that the malignant EPC_1 exhibited increased activity of transcription factors associated with EMT, such as FOXC1, FOXD1, and ISL1 (Fig. S3E). These results suggest that malignant EPC_1 contributes to CRC metastasis and is linked to unfavourable patient prognoses.

Fig. 3

A malignant epithelial cell subtype linked to higher CRC stages. (A) t-SNE visualization of epithelial cell composition, with colors indicating different tissue types. (B) Feature plot illustrating the copy number variation (CNV) score for epithelial cells. (C) Comparison of CNV scores between malignant and normal epithelial cells. (D)UMAP visualization of malignant epithelial cell composition, with colors representing different cell types. (E) Stacked bar plots illustrating the distribution of malignant epithelial cell subtypes across different stages of CRC. (F) Survival curve of KLK10, a marker of malignant EPC_1, in public datasets. (G) Dotplot showing scores of 50 hallmark gene sets for 8 malignant epithelial cell subtypes

Characterization of the spatial niche associated with EMT

We hypothesized that there might be some connection between malignant EPC_1 and CTHRC1+ fibroblasts. To explore this, we conducted spatial transcriptomics analysis using data from Cancer Diversity Asia. The unsupervised clustering analysis classified the spatial tissue sections into several distinct regions, including tumor, fibroblast, normal epithelial, smooth muscle, and lamina propria areas (Fig. S3A). To determine the cellular composition of each spatial spot, we utilized the cell2location deconvolution algorithm, which maps scRNA-seq-defined cell states onto spatial tissue sections. Unsupervised clustering based on spot composition identified 10 distinct cellular niches (Fig. 4A, Fig. S3B). The spatial niche types of each spot are displayed in Table S5. Spatial niches 1–3 were predominantly characterized by immune and stromal cells, with minimal tumor cell presence. The epithelial cells were located primarily in niches 4–5 and 7–10. Notably, spatial niche_4 contained a significant concentration of malignant EPC_1 and CTHRC1+ fibroblasts, suggesting a colocalization relationship and potential close communication between them. In spatial slices, we also observed that CTHRC1+ fibroblasts clustered together with malignant EPC_1 (Fig. 4B, Fig. S3C). Next, we assessed the spatial dependency between EPCs and fibroblasts from three perspectives using MISTy analysis. Among the malignant EPCs, malignant EPC_1 had the strongest predictive association with the surrounding abundance of CTHRC1+ fibroblasts, suggesting possible interdependence (Fig. 4C, Fig. S3D). The mIHC confirmed the proximity between malignant EPC_1 and CTHRC1+ fibroblasts (Fig. 4D, Fig. S3E). Niche_4 also demonstrated more active EMT processes compared to other spatial niches (Fig. 4E). These findings led us to speculate that CTHRC1+ fibroblasts might promote the EMT of malignant EPC_1. Additionally, spatial niche_4 exhibited higher infiltration of TREM2 + macrophages, SPP1+ macrophages, regulatory T cells, and type 17 helper T cells compared to other epithelial-containing niches. Spatial niche_4 had minimal infiltration of CD8+ T cells, CD4+ effector T cells, and B cells. These results indicate that the TME in niche_4 may play a role in facilitating tumor cell evasion from immune surveillance and attack. Previous research has demonstrated that CAFs can contribute to immunotherapy resistance by excluding T cell infiltration into tumors [30]. We subsequently investigated the association between CTHRC1+ fibroblasts and immune evasion in CRC patients. Using the Tumor Immune Dysfunction and Exclusion (TIDE) [31] Database, we assessed immune exclusion and dysfunction in these patients. Our results showed a strong positive correlation between the presence of CTHRC1+ fibroblasts and the exclusion of CD8+ T cells from the tumors (Fig. S3F). Our comprehensive results describe the colocalization events of CTHRC1+ fibroblasts, malignant EPC_1, and other cell types, revealing active EMT processes and both T cell exclusion and inhibitory activities within spatial niche_4.

Fig. 4

CTHRC1+ fibroblasts colocalize with malignant EPC_1 in spatial niches characterized by high EMT activity. (A) Scaled median cell-type compositions within each spatial niche. (B) Deconvoluted abundance of malignant EPC_1 and CTHRC+ fibroblasts in a CRC sample. (C) Heatmap displaying the significance of malignant epithelial cell subtypes in predicting the abundances of fibroblast subtypes within a spot. (D) Representative IF staining of human CRC tissues. DAPI (blue), CTHRC1 (red), KLK10 (green) and EPCAM (yellow) are shown in merged channels. Bar, 100 μm. The experiment was conducted on samples from three independent patients. (E) Bar chart showing upregulated and downregulated pathways in spatial Niche_4, with upregulated hallmark pathways on the right and downregulated hallmark pathways on the left

Interactions between CTHRC1+ fibroblasts and malignant EPC_1

Given that juxtacrine and paracrine signals operate only within a limited spatial range [32], we restricted our cell communication analysis to the cell types within spatial niche_4. We used the CellPhoneDB and CellChat packages to establish communication networks among the predominant cell types in this niche. We observed that CTHRC1+ fibroblasts in niche_4 had the highest number and strength of interactions with malignant EPC_1, and the signal flow primarily directed from CTHRC1+ fibroblasts to malignant EPC_1 (Fig. 5A-C). Further analysis revealed two distinct communication patterns within niche_4 based on ligand-receptor interaction strengths (Fig. S4A). Pattern 1, which included malignant EPC_1, CTHRC1+ fibroblasts, and F3+ fibroblasts, was driven predominantly by pathways such as the WNT, ncWNT, Notch, FN1, and COLLAGEN pathways (Fig. S4A). This suggests that fibroblasts are actively involved in remodelling the extracellular matrix, regulating cell polarity and movement, and promoting tumor cell invasion and migration. Additionally, LIANA combined with NicheNet was used to predict interactions between ligands, receptors, and target genes between malignant EPC1and CTHRC1+ fibroblasts. The heatmap highlighted the top 20 most active ligands predicted by LIANA for CTHRC1+ fibroblasts that act on malignant EPC_1 (Fig. 5D). Ligands such as APOE, THBS2, COMP, and WNT5A were highly active. We subsequently investigated the potential receptors for the top 20 ligands with the highest activity (Fig. 5E). It was found that CTHRC1+ fibroblasts interact directly with malignant EPC_1 through adhesion ligand-receptor pairs such as COL1A1-ITGB1 and LAMA1-ITGA3. Additionally, CTHRC1+ fibroblasts enhanced the migration and invasion of malignant EPC_1 through WNT5A-LRP6/FZD3 and THBS1-ITGAV interactions. To identify the downstream target genes influenced by CTHRC1+ fibroblasts in malignant EPC_1, we intersected the potential targets of the top 20 ligands (Fig. 5F), niche_4-specific genes (Table S6), and genes specifically expressed in malignant EPC_1 (Table S4). We found that MSLN was the only shared downstream target gene (Fig. 5G). The expression of MSLN in the malignant EPC_1 of stage III-IV CRC patients was notably higher compared to that in stage I-II CRC patients (Fig. 5H). Importantly, MSLN expression in malignant EPCs was regulated primarily by WNT5A produced by CTHRC1+ fibroblasts (Fig. 5F). At the single-cell level, we then found that the WNT5A gene was predominantly expressed in CTHRC1+ fibroblasts (Fig. S4B). Among various cell types in the spatial niche_4, CTHRC1+ fibroblasts were also identified as the main source of WNT and ncWNT signaling pathways (Fig. S4C, D). These results suggest that MSLN expression in malignant EPCs is regulated by WNT5A secreted by CTHRC1+ fibroblasts and may be crucial for tumor differentiation, migration, and metastasis.

Fig. 5

Extensive interactions exist between CTHRC1+ fibroblasts and malignant EPC_1. (A) Heatmap illustrating the predicted number of potential ligand-receptor pairs between cell subtypes in spatial niche_4 by CellphoneDB. Heatmap showing the strength (B) and number (C) of potential cell interactions between cell subtypes in spatial niche_4 predicted by CellChat. (D) Heatmap displaying the top 20 most potent ligands from CTHRC1+ fibroblasts targeting malignant EPC_1, as predicted by LIANA and CellChat software. (E) Ligand-receptor pairs showing interaction between CTHRC1+ fibroblasts and malignant EPC_1, as predicted by NicheNet software. (F) Heatmap showing regulatory potential of top 20 ranked ligand and the downstream target genes in malignant EPC_1, as predicted by NicheNet software. (G) Venn diagram of potential target genes for the top 20 ligands, Spatial Niche_4-specific genes, and genes specifically expressed in malignant EPC_1. (H) Violin plot of expression of MSLN in malignant EPC_1 at the single-cell RNA level

WNT5A derived from CTHRC1+ CAF can upregulate MSLN expression in CRC cells, indicating poor prognosis for patients

To verify whether CTHRC1+ CAF can affect MSLN expression in CRC cells by secreting WNT5A, we first explored the relationship between CTHRC1 and WNT5A. Previous bioinformatics analysis revealed that WNT5A is highly expressed in CTHRC1+ CAF (Fig. S4B). Further analysis using the TCGA database demonstrated a strong correlation between CTHRC1 and WNT5A expression in CRC patients (Fig. S5B). Prior studies have also indicated that CTHRC1 can upregulate WNT5A expression in cervical carcinoma and gastrointestinal stromal tumors [33, 34]. To further validate this mechanism, we isolated primary CAFs (Fig. S5A) and overexpressed or knocked down CTHRC1 in these cells (Fig. S5C), followed by the observation of changes in WNT5A levels. We found that overexpression of CTHRC1 in CAFs caused a marked increase in both the expression and secretion of WNT5A, whereas knocking down CTHRC1 in CAFs resulted in the opposite effect (Fig. S5D). This suggests that CTHRC1+ CAF can secrete more WNT5A compared to other CAF types. Next, we co-cultured CAFs with CRC cells to further demonstrate the impact of CTHRC1+ CAF on CRC cells via WNT5A. The results showed that overexpression of CTHRC1 in CAFs caused a marked increase in MSLN expression in CRC cells. Notably, when using the WNT5A inhibitor Box5, the MSLN expression in CRC cells was significantly reduced. Similar results were observed by knocking down CTHRC1 in CAFs (Fig. S5E).

To elucidate the role of MSLN in tumors, we utilized the online tool TIMER [33] to explore its expression across various cancers, and find significantly higher expression in colon, rectal, gastric, and urothelial bladder cancers (Fig. 6A). High MSLN expression in the TCGA database was correlated with poor prognosis in CRC patients (Fig. 6B). Moreover, MSLN expression was greater in stage III-IV CRC patients compared to those in stage I-II (Fig. 6C). To confirm the importance of MSLN in CRC, we performed qRT-PCR and IHC on tumor and paired normal tissue samples from 80 CRC patients, as well as Western blot analysis on samples from 12 patients. MSLN mRNA levels were notably elevated in CRC tumor samples compared to normal tissues (Fig. 6D). MSLN protein levels were significantly higher in tumor tissues than in matched normal tissues (Fig. 6E). IHC analysis of a tissue microarray containing samples from 80 CRC patients confirmed the significant upregulation of MSLN in CRC tissues and its association with poor survival outcomes (Fig. 6F-H). By integrating the tumor staging information, we found that MSLN protein levels were elevated in stage III-IV patients compared to stage I-II CRC patients (Fig. S6A), indicating its association with CRC progression. In summary, our study indicates that CTHRC1+ CAF upregulates MSLN expression in CRC tumor cells by secreting the WNT5A, thereby influencing patient prognosis.

Fig. 6

High MSLN expression predicts advanced clinical stage in CRC. (A) Boxplots illustrating the expression of MSLN in cancer and adjacent normal tissues across multiple cancer types. (B) Kaplan–Meier curves showing overall survival analysis of CRC patients with low and high MSLN expression in TCGA datasets. (C) Comparison of MSLN expression levels between stage I-II and stage III-IV CRC patients in TCGA datasets. (D) MSLN mRNA levels measured by qRT-PCR in 80 CRC tissues and matched adjacent normal tissues. (E) MSLN protein levels detected by Western blot in 12 paired CRC tissues and adjacent normal tissues. (F, G) The protein level of MSLN detected by IHC in 80 CRC tissues and paired adjacent normal tissues. (H) Kaplan–Meier analysis of the overall survival of 80 CRC patients based on MSLN protein levels

CTHRC1+ CAF regulates MSLN expression to impact invasion and EMT activity in CRC cells

To investigate the mechanism through which MSLN influences CRC, we first assessed MSLN expression across seven different CRC cell lines. MSLN was significantly upregulated in highly invasive LoVo, SW620, and HCT116 cells (Fig. S6B). Next, we established stable shRNA-mediated MSLN knockdown in HCT116 cells, and stable MSLN overexpression in DLD-1 cells. qRT-PCR and Western blot analysis were used to validate the efficiency of MSLN knockdown and overexpression (Fig. S6C). The CCK-8 and colony formation assays revealed that cells with high MSLN expression exhibited stronger proliferative abilities (Fig. S6E-G). Xenograft tumor models indicated that MSLN overexpression promoted tumor growth, resulting in larger tumor volumes and weights compared to controls, while MSLN knockdown inhibited tumor progression (Fig. S6H). Wound healing and transwell assays revealed that MSLN overexpression enhanced the migratory and invasive abilities of DLD-1 cells (Fig. 7A, C), whereas MSLN knockdown had the opposite effect on HCT116 cells (Fig. 7B, D). GSEA analysis revealed that samples with high MSLN expression had higher EMT, myogenesis, and hypoxia scores than those with low MSLN expression (Fig. S6D). Additionally, MSLN knockdown resulted in elevated E-cadherin levels, along with a reduction in N-cadherin and Vimentin expression, whereas MSLN overexpression enhanced EMT activity in DLD-1 cells (Fig. 7I).

Fig. 7

CTHRC1+ CAF regulates MSLN expression via WNT5A to impact migration, invasion and EMT activity in CRC cells (A, C) Overexpression of MSLN enhanced migration and invasion in DLD-1 cells. (B, D) knockdown of MSLN weakened migration and invasion in HCT116 cells. (E, G) Overexpression of MSLN or addition of rhWNT5A in DLD-1 cells reversed the effect of CAFs-shCTHRC1 on DLD-1 cell migration and invasion. (F, H) Knockdown of MSLN or treatment with the WNT5A inhibitor Box5 in HCT116 cells reversed the effect of CAFs-CTHRC1 on HCT116 cell migration and invasion. (I) Western blot analysis demonstrated that MSLN influences the expression of key proteins involved in EMT. (J) Western blot analysis showed that knockdown of MSLN or treatment with the WNT5A inhibitor Box5 in HCT116 cells reversed the effect of CAFs-CTHRC1 on EMT activity in HCT116 cells. In contrast, overexpression of MSLN or addition of rhWNT5A in DLD-1 cells rescued the effect of CAFs-shCTHRC1 on EMT activity in DLD-1 cells. Scale bar: 200 μm. All data are presented as the means ± SDs of three independent experiments. **p < 0.01, and ***p < 0.001

Next, we investigated whether CTHRC1+ fibroblasts can influence the invasion ability and EMT activity of CRC cells by co-culturing CRC cells with CAFs and performing transwell and wound healing assays. The co-culture schematic representations for the transwell and wound healing assays are shown in Fig. S5F and G. We found that overexpression of MSLN or addition of rhWNT5A in DLD-1 cells reversed the effects of CAFs-shCTHRC1 on migration, invasion, and EMT activity in DLD-1 cells (Fig. 7E, G, J). In contrast, knockdown of MSLN or treatment with the WNT5A inhibitor Box5 in HCT116 cells could reverse the enhancing effects of CAFs-CTHRC1 on the migration, invasion, and EMT activity of HCT116 cells (Fig. 7F, H, J). Our results show that CTHRC1+ CAF regulates WNT5A, which subsequently modulates MSLN expression, thereby enhancing migration, invasion, and EMT activity in CRC cells, ultimately driving tumor progression.

Discussion

Recent studies have underscored the essential role of the TME in regulating tumor progression, metastasis, and the response to therapy [35]. CAFs are one of the most abundant components in the TME of solid tumors, and are characterized by significant heterogeneity and plasticity. They engage in extensive interactions with both immune cells and cancer cells, making them promising therapeutic targets in oncology [36, 37]. In the present study, we conducted a combined analysis of scRNA-seq and ST, revealing that WNT5A, secreted by CTHRC1+ CAFs, acts on adjacent EPCs, leading to the upregulation of MSLN expression. This, in turn, promotes EMT and enhances the invasiveness of CRC cells. These findings significantly enhance our understanding of the intricate interactions between CAFs and CRC cells within the TME, emphasizing the biological and therapeutic significance of the novel CTHRC1+ CAF-WNT5A-MSLN signaling axis.

CAFs can promote tumor proliferation, invasion, and metastasis and are generally associated with poor prognosis [38]. However, owing to the elusive heterogeneity and functional plasticity of CAFs, most cancer treatments do not target them directly [39]. Some studies indicate that the indiscriminate depletion of CAF may accelerate tumor progression and reduce patient survival [40,41,42]. Therefore, it is necessary to carry out extensive research on CAFs and their interactions with the TME [43]. This study found a specific enrichment of CTHRC1 + fibroblasts in CRC tumor tissues. CTHRC1, a secreted protein, has been linked to distant liver metastasis in CRC patients [44, 45], but its role in fibroblasts within CRC has rarely been studied. Our research revealed that the CTHRC1+ fibroblast subpopulation is particularly enriched in CRC patients with higher clinical stages. These CTHRC1+ fibroblasts are actively involved in matrix formation and remodelling, and exhibit strong activity in promoting EMT and angiogenesis. Notably, ECM remodeling hinders immune cell infiltration and impairs the effectiveness of immunotherapies [46]. We also found that CTHRC1+ fibroblasts could limit the infiltration of CD8+ T cells, thereby promoting tumor immune evasion and progression. Investigating the role of CTHRC1+ fibroblasts within the TME may yield important perspectives on the mechanisms underlying distant metastasis and the establishment of an immunosuppressive microenvironment in CRC.

Our research identified a strong spatial colocalization between CTHRC1+ fibroblasts and malignant EPC_1, both of which are situated within the EMT-active spatial niche_4. Malignant EPC_1 also exhibited high EMT scores, suggesting that CTHRC1+ fibroblasts may facilitate in the EMT process in this subtype. During the analysis of intercellular crosstalk in spatial niche_4, we observed that interactions between CTHRC1+ fibroblasts and malignant EPC_1 were the most frequent, with the signal flow directed primarily from CTHRC1+ fibroblasts to malignant EPC_1. CTHRC1 + fibroblasts can send highly active ligands such as WNT5A, COMP, and THBS2 to malignant EPC_1. Earlier investigations have shown that these ligands can induce EMT and enhance tumor metastasis [47,48,49]. Notably, WNT5A is associated with cancer progression and plays dual roles as both an oncogene and a tumor suppressor in CRC [50]. This discrepancy may be attributed to the different cellular sources of WNT5A. Previous research has shown that WNT5A is involved in ncWNT signaling, coordinating cell polarity and cytoskeletal rearrangement, and is expressed at relatively high levels in highly invasive cancer subtypes, such as prostate cancer, metastatic melanoma, and ovarian cancer [50,51,52]. Consistent with these findings, our study underscores the pro-tumorigenic role of WNT5A in the TME. We also confirm that CTHRC1+ CAFs are a major source of WNT5A, providing a new perspective on how CAFs influence tumor progression within the TME.

Our subsequent research revealed a previously unreported connection between WNT5A signaling and the upregulation of MSLN in CRC cells. MSLN is highly expressed in ovarian cancer, pancreatic cancer, mesothelioma, and other cancers [53, 54]. Previous research suggests that MSLN promotes cell metastasis by inducing the expression of matrix metalloproteinases [55]. MSLN is a marker of poor prognosis in stage II/III CRC patients and can guide preoperative chemotherapy for stage III patients [56, 57]. Most research on MSLN has focused on other solid tumors such as ovarian, breast, and lung cancers, with relatively few investigations and targeted therapies focused on MSLN in CRC. Our study suggests that targeting MSLN can block the downstream effects of WNT5A secreted by CTHRC1+ fibroblasts, which promotes EMT in malignant EPC_1. These findings suggest that MSLN holds considerable potential as a therapeutic target for advanced CRC patients.

Conclusions

In summary, this study employed multi-omics approaches to explore the interactions between CAFs and malignant EPCs in CRC. During this process, we identified a novel CTHRC1+ CAF-WNT5A-MSLN signaling axis that plays a pivotal role in promoting EMT and enhancing tumor invasiveness in CRC. This axis underscores the critical contribution of CTHRC1+ CAFs in shaping the TME and driving CRC progression. The identification of MSLN as a downstream effector of WNT5A signaling provides new insights into how CAF-derived signals regulate tumor behavior and facilitate metastasis. Ultimately, our research advances the understanding of CRC pathogenesis and opens new avenues for developing targeted therapeutic strategies aimed at the CTHRC1+ CAF-WNT5A-MSLN axis, offering potential approaches to attenuate cancer cell invasiveness and metastatic potential in CRC patients.