Identification of MC signature genes and decreased MC density in CRC

In this study, we analyzed three large CRC single-cell datasets (GSE164522, GSE178341, and 5-cohorts) comprising 341 samples and a total of 953,493 cells, identifying 8,875 MCs (Fig. 1a). We annotated major clusters based on defining marker genes and identified various cell types, including T cells (CD3D, CD3E), natural killer (NK) cells (KLRF1, GNLY), B cells (CD19, MS4A1), plasma cells (MZB1, IGHA1), MCs (TPSAB1, TPSB2, CPA3), myeloid cells (LYZ, CD68), endothelial cells (PECAM1, VWF), epithelial cells (EPCAM, AGR2), and fibroblasts (DCN, COL1A2) (Fig. 1a and b, Table S3).

Fig. 1

Identification of MC signature genes and decreased MC density in CRC.

a. UMAP plots displaying the major cell types in the GSE164522 (n = 52 samples), GSE178341 (n = 100 samples), and 5-cohorts (n = 189 samples) datasets. b. Dot plots of marker genes for each major cell type in the 5-cohorts dataset. c. Venn diagrams (center) of differential expressed genes (DEGs) in MCs in the GSE164522, GSE178341, and 5-cohorts datasets, with the intersection showing the five MC signature genes (right). The criteria for screening DEGs are depicted within the dotted lines (left). d. Heatmap showing the expression of MC signature genes in each cancer (Normal vs. Tumor). Histogram shows the number of genes with statistical significance (upper). Red represents an increase in tumor expression, green represents a decrease in tumor expression, and only p-values < 0.05 are displayed. e. Heatmap showing the expression of the MC signature genes in 10 bulk RNA-seq cohorts of CRC (Normal vs. Tumor). Blue represents a decrease in tumor expression, dataset source (top), sample size (bottom). f. The expression of TPSAB1 in human CRC tissue and paired NC tissue by Western blotting. g. Immunofluorescence staining of human CRC tissue and paired NC tissue. TPSAB1 (pink), CMA1 (green), CPA3 (green), DAPI (blue), Bar, 200 μm. CRC, colorectal cancer; DEGs, differential expressed genes; FC, fold change; FDR, False Discovery Rate; MC, mast cell; NC, normal colorectum or adjacent colorectum

To identify highly expressed genes that could serve as MC markers, we established strict criteria for differential gene screening, including a log2 fold change > 1 a proportion of expressed cells in MCs (PCT1) > 0.7, a proportion of expressed cells in other cell types (PCT2) < 0.3, and a P-value < 0.01. Five genes (TPSAB1, TPSB2, CPA3, HPGDS, and MS4A2) were identified as common DEGs across all three single-cell datasets and defined as MC signature genes (Fig. 1c).

The MC signature genes were used as MC markers to assess the density of MCs in bulk RNA-seq samples of both tumors and normal tissues. All five MC signature genes that achieved statistically significant (P < 0.05) were considered credible. Results showed that MC signature genes were significantly increased only in KICH, KIRC, and THCA but significantly decreased in BLCA, CESC, COAD, ESCA, LUAD, LUSC, READ, STAD, and UCEC (abbreviations of cancers are represented in Table S4), suggesting a reduction in MC density in the majority of tumors (Fig. 1d). Subsequent analyses in additional CRC cohorts verified a significant decrease in MC density in CRC (Fig. 1e). This was further corroborated by our Western Blot (Fig. 1f) and immunofluorescence (Fig. 1g) results derived from paired CRC and NC samples, which also pointed to a significant decrease in MC density within CRC.

In summary, we identified five reliable MC signature genes as markers, and our findings consistently demonstrated a reduction in MC density in CRC.

Activation of MCs in CRC

Activated MCs refer to MCs that have been stimulated by external factors, and they release biologically active substances, including histamine, cytokines, proteases, and lipid mediators such as leukotrienes and prostaglandins [6, 9]. Our research on the activation of MCs in CRC began with an unexpected finding. Using the CIBERSORT algorithm, we compared immune cell ratios between CRC and normal tissues in TCGA-CRC data and observed a significant increase in the proportion of activated MCs and a decrease in the proportion of resting MCs in CRC (Fig. 2a). This finding was further confirmed by data from nine additional CRC cohorts.

Fig. 2

Activation of MCs in CRC.

a. CIBERSORT-based analysis-generated heatmap showing the difference in the proportion of activated and resting MCs between NC and CRC in 10 bulk RNA-seq cohorts of CRC. Red indicates a higher proportion in CRC, while blue indicates a lower proportion in CRC. b. DEGs between NC MCs and CRC MCs from the GSE164522, GSE178341, and 5-cohort datasets (top). Venn diagram of upregulated DEGs in CRC MCs (bottom). Screening criteria are indicated by the dotted lines. c. Gene Ontology Biological Process (GO BP) enrichment analysis of upregulated DEGs in CRC MCs. d. Heatmap showing the expression of cytokine and growth factor, protease and histamine, lipid mediator, and various receptor-related genes across different major cell types (5-cohorts). e. Heatmap showing the expression of cytokine and growth factor, protease and histamine, lipid mediator, and various receptor-related genes in MCs (NC vs. CRC) (5-cohorts). The tissue type is indicated by the color above the heatmap. f. Heatmap showing the expression of cytokine and growth factor, protease and histamine, lipid mediator, and various receptor-related genes in MCs (NC vs. CRC) (GSE178341)

To gain a deeper understanding of the changes in MCs during the tumorgenesis of CRC, we compared gene expression differences between MCs in CRC tissue and normal tissue in three single-cell datasets (GSE164522, GSE178341, and 5-cohorts) (Fig. 2b and Table S5). Genes that met the criteria of log2 fold change > 0.25, a proportion of expressed cells in CRC MCs (PCT1) > 0.25, and a P-value < 0.05 were defined as DEGs. The number of DEGs enriched in MCs in CRC was significantly higher in all three datasets compared to those enriched in normal tissue, indicating a widespread gene expression increase in MCs during the progression of CRC. Notably, TMEM176B and CD52 were among the top 5 significant genes in all three datasets, with CD52 being reported as a marker of neoplastic MCs in patients with advanced systemic mastocytosis [31]. To study the functional changes in MCs during the progression from normal tissue to CRC, we used DEGs (n = 268) that were increased in two or more datasets for enrichment analysis. The results of the GO analysis showed that pathways related to cell activation were the most significantly enriched in MCs in CRC (Fig. 2c).

To better understand the activation features of MCs during CRC progression, this study analyzed the expression of receptor and mediator genes related to MCs in the 5-cohorts dataset (Fig. 2d). The results showed that MCs were the main population that expressed receptors for IL-33 (IL1RL1) and KIT, with the highest expression levels of MRGRPX2, CSF2RB, and AHR. Notably, MCs were also the only cell population that expressed all three subunits of the high-affinity IgE receptor FcεR1 (i.e., FCER1G, FCER1A, and MS4A2). In addition, MCs showed high expression of signature proteases, including TPSAB1, TPSB2, CPA3, CMA1, and CTSG, with CTSD and CTSW being enriched in MCs but not limited to them. Moreover, MCs displayed high expression of genes involved in histamine biosynthesis (HDC), leukotriene biosynthesis (LTC4S, ALOX5OP, and ALOX5), and prostaglandin biosynthesis (HPGDS, PTGS1, PTGS2). MCs were the only cell population that expressed mRNA encoding diverse cytokines, chemokines, and growth factors, including IL4, IL5, IL9, IL13, CCL1, LIF, CSF1, and AREG. MCs also showed high expression of IL18, VEGFA, and TGFB1, although expression of these genes was not restricted to MCs.

The expression of MC receptors and mediators in CRC MCs and normal MCs was compared in the 5-cohorts (Fig. 2e) and GSE178341 (Fig. 2f) datasets. The results showed that most genes were significantly upregulated in CRC MCs, indicating that CRC MCs have a more activated MC phenotype compared to normal MCs. The exception was CMA1, which encodes chymase and was significantly more highly expressed in normal tissue.

In summary, the evidences indicate activation of MCs in CRC.

Heterogeneity of MCs in CRC

The advent of single-cell analysis has enabled the characterization of MC activation during CRC from the perspective of MC heterogeneity. In the GSE178341 cohort, analysis of 4,155 MCs led to the identification of 12 clusters corresponding to 4 distinct MC subtypes (Fig. S1a). The MC11 and MC12 clusters, enriched in MZB1 and CD3D respectively, were considered B cell doublets and T cell doublets, respectively, while the MC09 and MC10 clusters, enriched in interferon-related genes (IFITs) and mitochondrial genes (MTs), respectively, were grouped as “Other MCs”. The MC08 cluster, enriched in genes related to proliferation such as MKI67, was named as “proliferation MCs”. Based on the overall expression levels of MC receptor and mediator genes, the relatively high MC01-04 clusters were named “activated MCs”, while the MC05-07 clusters were named “resting MCs” (Fig. 3a and Fig. S1a). The significant enrichment of the MC activation signature in activated MCs compared to resting MCs further supported the naming of these MC clusters (Fig. 3b). The proportion of activated MCs in CRC was significantly higher, while the proportion of resting MCs was significantly higher in normal tissue (Fig. 3a and c), indicating that the activation of MCs in CRC is due to a higher proportion of activated MCs.

Fig. 3

Heterogeneity of MCs in CRC.

a. UMAP plots of 4,155 MCs colored by cluster (left) and tissue type (center) in the GSE178341 dataset. Bar charts show the proportion of MC subtypes in different tissues (right). b. UMAP plots of MC activated signature. c. Log ratio of average fraction per MC clusters in tumor to normal tissue (top). Mann-Whitney U test, *: p < 0.05, **: p < 0.01, ***: p < 0.001. Dot plots of cytokine and growth factor, protease and histamine, lipid mediator and various receptor-related gene expression in MC clusters (bottom). d. Volcano plot of differentially expressed genes between resting and activated MCs. e. Differential pathway enriched in resting and activated MCs by GSVA, showing the top 10 significant enriched hallmark terms. f. Differentiation trajectory of MCs with each color coded for MC subsets (left), tissue types (center), and pseudotime (right). g. Pseudotime trajectory of CMA1, CPA3, TPSAB1, and KIT expression levels

In the comparison of activated MCs and resting MCs, it was found that most MC receptor and mediator-related genes, including the five major MC signatures (TPSAB1, TPSB2, CPA3, HPGDS, and HS4A2), were enriched in activated MCs, while CMA1 was enriched in resting MCs (Fig. 3d and Table S6). GSVA analysis revealed that the TNFA signature via NF-κB was most significantly enriched in activated MCs, whereas the angiogenesis-related pathway was enriched in resting MCs (Fig. 3e). This finding was supported by the result of the angiogenesis score, which confirmed the higher angiogenesis feature of resting MCs compared to activated MCs (Fig. S1b). Additionally, the enrichment of MHC-I and MHC-II related genes in activated MCs indicated that activated MCs have a stronger antigen-presenting function (Fig. S1c). The heterogeneity observed in MCs within CRC was also corroborated in the 5-cohorts dataset (Fig. S2a-d).

In addition, our pseudo-temporal analysis using Monocle 2 further supported a transition from resting to activated MCs (Fig. 3f), during which we also observed a decline in CMA1 expression (Fig. 3g).

High MC signature associated with favorable outcome in tumors

The impact of MCs on cancer prognosis remains controversial [9,10,11]. Based on the Kaplan-Meier model of TCGA pan-cancer data, we found that high expression of five MC signature genes was associated with better disease-free survival (Fig. 4a). In addition, we used eight prognostic indicators based on univariate Cox regression and Kaplan-Meier models for OS, DSS, DFI, and PFI to evaluate the impact of MC signature on prognosis in different cancers. A reliable result was considered if statistical significance (p < 0.05) was reached in at least four indicators. The results showed that the MC signature had a significant protective effect in 10 cancer types, including ACC, CESC, CHOL, HNSC, KIRC, KIRP, LIHC, LUAD, PRAD, and SARC, but was only associated with poor prognosis in STAD (Fig. 4b). Kaplan-Meier analysis of the TCGA-CRC (log-rank, p = 0.019) and GSE39582 (log-rank, p = 0.029) also indicated that a high MC signature was associated with better overall survival in CRC cohorts (Fig. 4c and d, Fig. S3a, and Fig. S3b).

Fig. 4

MC signature predicts better prognosis.

a. Kaplan-Meier disease-free survival curves grouped by MC signature genes (TPSAB1, TPSB2, CPA3, HPGDS, and MS4A2) in pan-cancer. b. Summary of the correlation between the expression of MC signature and overall survival (OS), disease-specific survival (DSS), disease-free interval (DFI), and progression-free interval (PFI) based on univariate Cox regression and Kaplan-Meier models. Red indicates factors that are detrimental to the prognosis of cancer patients, while green represents protective factors. Only p-values < 0.05 are displayed. c. The Kaplan-Meier overall survival curve of the MC signature in TCGA-CRC is shown, with the High-MC group and Low-MC group including patients with CRC who had MC signature expression in the top 30% and bottom 30%, respectively. d. Kaplan-Meier overall survival curve of MC signature in GSE39582

Additionally, we employed CIBERSORTx to delve into the influence of MC phenotypes on prognosis and observed a heightened proportion of activated MCs in CRC samples compared to normal tissues, concomitant with a reduction in resting MCs (Fig. S3c). However, no statistically significant differences were discovered concerning the proportions of activated and resting MCs calculated via CIBERSORTx in relation to CRC prognosis (Fig. S3d and Fig. S3e).

These findings demonstrate the important role of MCs in the prognosis of CRC patients and their potential as a protective prognostic biomarker.

KITLG/KIT signaling in MC activation and CRC inhibition

In this section, we aimed to identify the potential causes of MC activation in the CRC TME. Using the GSE178341 dataset, we performed CellphoneDB analysis and found that activated MCs had a higher number of interactions with other cell types compared to resting MCs. The prediction results also showed that activated MCs had the highest number of interacting receptors with myeloid cells, endothelial cells, fibroblasts, and epithelial cells (Fig. 5a).

Fig. 5

KITLG/KIT signaling in MC activation and CRC inhibition.

a. Heatmap generated by CellphoneDB analysis showing the potential ligand-receptor interactions between resting and activated MCs and other major cell types in CRC (GSE178341). Numbers indicate the number of potential ligand-receptor pairs. b. Dot plots of interactions between resting and activated MCs and other major cell types along the IL33-IL1RL1 and KITLG-KIT axes. c. Dot plots displaying the expression of IL1RL1, IL33, KIT, and KITLG in different major cell types (5-cohorts) (left), and UMAP plots displaying the expression of IL33 and IL1RL1 (right upper), and KITLG and KIT (right bottom). d. Dot plots showing the expression of IL1RL1, IL33, KIT, and KITLG in different major cell types (GSE178341). e. Bar plots comparing the positive rate of KITLG expression between normal tissue and CRC in fibroblasts (left) and endothelial cells (center). Dot plots showing the expression of KITLG in different endothelial cell subsets (right). f. Correlation analysis of KITLG and KIT expression in TCGA-CRC (Spearman test). g. Correlation analysis of KITLG and KIT expression in GSE39582 (Spearman test). h. qRT-PCR analysis shows KIT mRNA expression level in P815 after manipulating different concentrations of KITLG protein. i. Western blot analysis shows KIT protein expression level in P815 after manipulating different concentrations of KITLG protein. j. CCK-8 assay comparing the proliferative capacity of CRC cells when exposed to medium with only different KITLG concentrations (left), compared to medium from p815 coculture with varying KITLG concentrations (right). Optical density (OD) was monitored daily for a 5-day period. k. Transwell analysis showing the impact on CRC cell migration and invasion when exposed to medium with only different KITLG concentrations (left), compared to medium from p815 coculture with varying KITLG concentrations (right). All data are shown as the mean ± SD. **: p < 0.01, ***: p < 0.001

Previous studies have shown that KITLG (SCF, stem cell factor) [32,33,34] and IL33 [27, 35, 36] are key molecules that can promote MC activation by binding to the surface receptor KIT and IL1RL1, respectively. To investigate the relationship between IL33-IL1RL1 and KITLG -KIT axes with MC activation, we conducted a CellphoneDB analysis and found that the average expression of the KITLG -KIT axis interacting with upstream cells was significantly higher in activated MCs compared to resting MCs (Fig. 5b). Furthermore, single-cell analysis identified the upstream cell types expressing IL1RL1 and KITLG. IL1RL1 was found to be almost exclusively expressed in endothelial cells (Fig. 5c and d), while KITLG was mainly highly expressed in endothelial cells and fibroblasts, followed by epithelial cells.

We compared the expression of KITLG in CRC and normal tissues and found that the expression of KITLG in fibroblast cells in CRC was significantly higher than that in normal tissues (p = 5.45E-12) (Fig. 5e) (Table S7). Similarly, the expression of KITLG in endothelial cells in CRC was also significantly higher than that in normal tissues (p = 1.97E-20). Further sub-group analysis of endothelial cells revealed that KITLG expression is enriched in tip cells, a subset marked by ESM1 and PGF and significantly elevated in CRC (Fig. 5e) (Table S8) [37]. Consistent with a recent study that reported higher KITLG expression in ACTA2 + vascular smooth muscle cells (VSMCs) compared to other stromal cell types [38], our findings also indicate that KITLG expression is enriched in mural cells, which are marked by RGS5 and ACTA2 (Fig. 5e). Moreover, a significant positive correlation between KITLG and KIT expression was observed in both the TCGA-CRC (Fig. 5f) and GSE39582 cohorts (Fig. 5g). These findings suggest that the increased expression of KITLG in endothelial cells and fibroblasts in CRC may be a significant cause of MC activation.

We further explored the effects of KITLG/KIT pathway activation in MCs on the proliferation and migration of CRC cells through in vitro experiments. Initially, PCR and Western Blot analyses were performed to assess KIT expression in the P815 cell line after 48 h of co-culturing with KITLG. Both assays confirmed an elevated KIT expression in the P815 cells (Fig. 5h and i), validating the successful activation of the KITLG/KIT signaling pathway in MCs. Subsequent CCK8 (Fig. 5j) and Transwell assays (Fig. 5k) revealed that as the concentration of KITLG increased, the presence of KITLG-activated P815 cells led to a significant reduction in both the proliferation and migration of CRC cells. In contrast, adding KITLG alone in the medium of MC38 cell line did not produce any notable changes in CRC cell behavior (Fig. 5j and k). These results suggest that KITLG can inhibit CRC cell proliferation and migration through its action on MCs.

In summary, our findings indicate that the KITLG/KIT signaling pathway may be a key mechanism for activating MCs, which in turn inhibits the proliferation and migration of CRC cells.

Spatial colocalization of MCs with fibroblasts and endothelial cells

To explore the relationship between the spatial distribution of MCs and their activation, we analyzed the spatial distribution of MCs in CRC tissue using spatial transcriptomics data. Based on transcriptome expression and HE tissue information, 3,138 spots were divided into three regions: normal region, stromal region, and tumor region (Fig. 6a). The stromal region is enriched with fibroblasts (DCN, COL1A2) and endothelial cells (PECAM1), and the five MC signature genes and KIT, KITLG are also significantly enriched in the stromal region, indicating that MCs co-localize spatially with fibroblasts and endothelial cells (Fig. 6b).

Fig. 6

Spatial co-localization of MCs with fibroblasts and endothelial cells.

a. Spatial plots of H&E staining (column 1), tissue regions (column 2), fibroblasts (COL1A2) (column 3), endothelial cells (PECAM1) (column 4), and MCs (TPSB2) (column 5) in a CRC sample. b. Dot plots display the expression of the regional marker genes (left), MC signature genes (middle), and expression of KIT and KITLG (right) in different tissue regions. c. Immunofluorescence staining of human CRC tissue. KITLG (red), KIT (green), DAPI (blue). Dotted lines demarcate the boundary between the tumor and the stromal regions. d. Spatial visualization of the MC activation signature in the CRC sample. e. Immunofluorescence staining of human CRC tissue and paired NC tissue. TPSAB1 (pink), CMA1 (red), KIT (green), DAPI (blue), in individual and merged channels are shown

Further investigations revealed an additional spatial co-localization of MCs with mural cells (identified by RGS5 and ACTA2) within the stromal region (Fig. S4a-c). We validated these transcriptomic insights with immunofluorescence assays, specifically confirming the co-localization of KIT and KITLG proteins within the stromal compartment of CRC tissues (Fig. 6c). In addition, the MC activation signature also showed enrichment in the stromal region (Fig. 6d), suggesting that the stromal environment may serve as a critical niche for MC activation in CRC.

Lastly, our immunofluorescence results from the paired CRC and NC samples revealed that MCs were uniformly distributed within the stromal region interspersed among the endothelial glandular structures (Fig. 6e). Crucially, we observed a significant reduction in the expression of CMA1 in individual MCs within CRC tissues, along with an increase in TPSAB1 and KIT, compared to the NC paired samples (Fig. 6e). This confirmed the shift from a resting to an activated phenotype in MCs during tumor formation.