3. Mapping microenvironment cell composition and molecular changes along malignant transformation
Stem cells are recognized as pivotal contributors to cellular biological functions, possessing remarkable capacities for self-renewal and differentiation. However, dysregulation of self-renewal during aging and persistent microenvironmental and macroenvironmental stressors can lead to cancer 6. This prompts an exploration into whether stem-like cells within the epithelial compartment can exhibit a continuum across diverse tissues and how the states of other cells in the tissue microenvironment dynamically change along this malignant continuum. To address this, we established a transcriptional malignancy index for each tissue type by identifying the differentially expressed genes (DEGs) in precancerous and cancerous lesions. Specifically, we compared the gene expression profiles of stem-like cells in diseased tissues with those in healthy donors. The obtained DEGs between the stem-like cells of each specimen and their healthy controls allowed us to calculate the principal components of the log2FC of these DEGs. Subsequently, we ordered the samples based on their positions along a spline fit in this space (Fig. 3A). The position in this ordering can be interpreted as the location/pseudo-time along a continuum ranging from healthy tissues to malignancies. Our analysis revealed a stereotyped progression in gene expression variations between normal and diseased stem-like cells, culminating in invasive tumors. Based on the expression dynamics analysis, we observed that the expression of HNF4A gradually increases as precancerous lesions approach malignant transformation in the gastrointestinal tract (Fig. 3B). Notably, different patterns of HNF4A usage were observed in precancer stages (e.g., AD v.s. SER, liver cirrhosis v.s. NAFLD, and SIM v.s. CAG). In precancers, HNF4A decreases to drive WNT signaling, while in carcinomas, it is upregulated in CSCs to drive malignancy. This aligns with previous findings associating HNF4A knockdown with an upregulation of WNT pathway signaling molecules, and HNF4A loss drives metabolic reprogramming at an early stage of pancreatic cancer progression 10. Moreover, our SCENIC analysis 11 identified HNF4A as a transcription factor (TF) specifically enriched in hepatocytes of HCC and NAFLD but not in healthy tissues, with its target genes (e.g., IRS1 and FGFR3) exhibiting similar expression dynamics (Fig. 3C, Figure S2A, Table S3). Consistent with prior work, increased gene expression and chromatin accessibility of HNF4A motifs were reported only after the transformation to CRC 8. The forced re-expression of HNF4α, along with the application of all-trans retinoic acid, demonstrated effectiveness in reducing the number of liver CSCs and potentiating the chemotherapeutic effect 12. As HNF4A acts as a selective agonist of the peroxisome proliferator-activated receptor gamma (PPARγ), our analyses suggest its potential as a novel biomarker and target for cancer prevention. Furthermore, we observed upregulated regulon activities for NCOR1, a frequently mutated gene, in stem-like cells of both cirrhosis and HCC (Fig. 3C), showing increased expressions along malignant progression (Fig. 3B). Additionally, GPX2, a glutathione peroxidase acknowledged for its upregulation in CRC 8, exhibited elevated expression in HCC, ESCC, and cSCC, as well as their premalignant stages (Figure S2B). GPX2 functions to alleviate oxidative stress through the reduction of hydrogen peroxide, thereby facilitating both tumorigenesis and metastasis 13.
To gain a comprehensive understanding of the TME landscapes, we next examined the cellular compositions and molecular features of immune, stromal, epithelial and endothelial lineages across different stage groups. Within the immune compartment (Figure S2C), we found that plasma cells (PLA), NK, and inflammatory monocytes (INMON) were highly enriched in precancerous lesions, and regulatory T cells (TREG), neutrophils (NEUT) and germinal center B cells (GC) in both precancer and cancer stage. M2 macrophages (M2MAC), CD4+ follicular helper (TFH), T-17 helper (TH17), and CD8+ terminally exhausted cells (CD8TEREX) were observed to be specifically enriched in cancers when compared to benign and healthy tissues (Fig. 3D). These implicated the potential mechanisms in facilitating immune evasion and inducing immune tolerance in more malignant tissues. SCENIC analysis revealed the higher activity of GATA3, IKZF2, BATF, RORA and FOXP3 in Tregs of each stage across tissue types (Table S3). In addition, HLF and ARID5B were identified as TFs of Tregs in precancer stages (e.g., N_HPV, AAH, NEOLP, EOLP), while regulons such as FOXO1 and STAT5B were specifically detected in cancers (e.g., IDC, AIS, MIAC, IAC and PDAC). Similarly, for INMON, FOXO3 and IRF5 were common TFs in healthy and early-stage tissues (e.g., BRCA1-mut); specifically, PRDM4 was the TF in precancer stages (e.g., Goiters and HT) while ATF2 and RUNX1 were enriched in cancers (e.g., ESCC and IDC, Figure S2D, Table S3). Multiple types of fibroblasts and endothelial cells were identified for the stromal compartment (Fig. 3D, Figure S2C). Among distinct subpopulations of fibroblasts, myofibroblasts (MYOFIB) were found to be enriched in both PME and TME (Fig. 3D), while hepatic or pancreatic stellate cells (HSC/PSC), intermediate CAFs (INCAF) and pericytes (PERI) with high expression of RGS5 were significantly enriched in malignant stages (Fig. 3D). Additionally, upregulated TF activities in MYOFIB/HSC were noted for PATZ1 in preneoplastic stages (e.g., Cirrhosis, AK, and CAG), while CEBPZ and SOX6 demonstrated elevated activities in malignant (e.g., cSCC and HCC) stages (Figure S2D, Table S3).
We further calculated the fractional contributions of each cell type to each sample as a function of position in the malignancy index. Some cell subpopulations showed strong correlations with disease progression along the index. In the epithelia, stem-like cell fractions in the samples gradually increased throughout the malignant transformation process. Conversely, the number of enterocytes decreases as adenomas transformed into carcinomas. In the secretory compartment, we observe a decreasing trend in fractions of pit mucous cells (PMC) and chief cells (CHIEF) with neoplastic transformation of the pre-cancerous gastric mucosa. In carcinomas, there was a general lack of differentiation toward the secretory lineage, resulting in the elimination of goblet cells (GOB, Fig. 3E). This aligns with a previous study documenting a reduction in goblet cells in non-mucinous colon adenocarcinomas 8. Knocking down MUC2, a mucin protein associated with goblet cells, resulted in an increased formation of adenomas and carcinomas in mice 8, implying that the loss of immature and mature goblet cells could potentially contribute to tumorigenesis 14. In the immune cell compartment, the proportions of NK, dendritic cells (DC) and progenitor exhausted CD8 + T cells (CD8TEXP) in liver tissues decreased as the disease progressed. However, GC cells increased in the colorectum and esophagus, TREG in the breast and M2MAC in the pancreas increased (Fig. 3E, Figure S2E). To further evaluate the functional variations in these evolving microenvironment cells as disease escalated, we analyzed the curated gene signatures (See Methods) and mapped the related molecular changes along malignancy continuum. AUCell algorithm 11 was designed to score the activity of a gene set in each cell. Here we utilized AUCell to delineate the transcriptional patterns of five major cell types. For instance, as the malignant progression occurred in liver, lung, oral cavity and thyroid tissues, CD8 + T-cell activation & effector molecules (e.g., FGFBP2, CX3CR1, FCGR3A and KLRG1), high cytotoxicity gene signature and T cell receptor (TCR) signaling gradually increased from the initial stages (Fig. 3G, Figure S2F). Moreover, we noted a progressively accelerated expression of stress response gene signatures (e.g., stress-related heat shock genes such as HSPA1A and HSPA1B), and the expression of nuclear factor (NF)-κB signaling, a key regulator of cellular stress response, was also enhanced (Fig. 3G). Advanced tumor tissues displayed increased expression of adhesion, IFN response, anergy and exhaustion-related signatures (Fig. 3G). In line with prior findings 15, CD8 + TEX cells can express high levels of cytotoxicity markers, indicating likely antigen experience (Figure S2F, G). CD4 + T cells gradually expressed anti-apoptosis, effector function genes and TCR signaling signatures. They markedly expressed co-stimulatory molecules (e.g., TNFRSF4, TNFRSF9, TNFRSF18), adhesion, and IFN response signatures in the very late stage of cancer. The metabolic switch of T cells along disease progression was in accordance with their dynamic changes in activation and effector functions. Notably, CD4 + and CD8 + neoantigen-reactive tumor-infiltrating lymphocytes (TILs) exhibited distinct signatures across different stages (Figure S2G). In the examined tissues, advanced stages showed increased levels of neoantigen-specific T cell signatures along with higher mutational burden (Fig. 2F, Figure S1C, Figure S2F, G), indicating an active immune response against the tumor. T cell exhaustion, characterized by reduced cytokine production and increased expression of inhibitory receptors in response to chronic antigen stimulation (e.g., neoantigens, HBV, HPV and cancer testis antigens), has been recognized as a primary mechanism of immune escape by cancers 16,17. Intriguingly, this trend of T lymphocytes was reversed during disease progression in the cervix, endometrium, pancreas, and prostate tissues (Fig. 3G, Figure S2F), with higher levels of NeoTCR8/4 signatures and mutation load found in the precancerous stage compared to those in advanced tumors (Figure S2G). We inferred that heterogeneous distribution and functional patterns of T cells during multistage tumorigenesis across various solid tissues mirror inter- and intra-tumor heterogeneity, emphasizing the broad impact of tissue types on the states of T cell subpopulations. Additionally, NK cells exhibited the lowest cytotoxicity but the highest stress score in advanced tumors. Other TME cells also underwent drastic reprogramming during malignant transition, with macrophages shifted from the M1 to M2 state and fibroblasts transitioned from an inflammatory CAF (iCAF) to a dominant myofibroblastic CAF (myCAF) phenotype (Fig. 3F). These observations suggested that the stromal cell remodeling and immune suppression in the TME accompanied the acquisition of stemness by tumor cells and coincided with plasticity onset and progression.
In summary, our observations concerning the enrichment of innate immune cells, such as INMON and NK in the early precancerous phases reflect the activation of innate immunity and its rapid onset in response to harmful stimuli such as microbial infection. Various leukocytes (e.g., macrophages, NK, and DC) bridge innate and specific immunity and hold a dominant position in immune surveillance and clearance of abnormal epithelial cells. If pro-inflammatory stimuli persist during the epithelial wound-healing process, chronic inflammation may ensue, leading to the persistence of inflammatory factors and tissue damage. Various phenotypically distinct immune/stromal subclusters may dynamically shift and closely interact with other PME/TME cells, playing a critical role in tumorigenesis and cancer evolution.
4. Heterogeneity of transcriptional programs in stem-like cells.
Cancer cells exhibit heterogeneity in their degree of differentiation, ranging from stem- or progenitor-like to fully differentiated. Gene expression programs are often used to characterize different cell identities and cell activities. To better understand advanced and highly heterogeneous tumors, we investigated the intratumor heterogeneity (ITH) programs through transcriptomic profiling of stem-like cells from tumors and their originating lesions. We then characterized the expression programs for each diseased sample using non-negative matrix factorization (NMF) approach (see Method details) 9,18. Of note, the number of programs varies for each disease across the 13 tissues. To find the recurrent patterns of ITH in stem-like cells and decipher consensus modules consisting of coexpressed genes within all expression programs, we next clustered these coexpressed genes to gene modules using K-means clustering (Fig. 4A, Figure S3A). After filtering out gene modules annotated with ‘unknown’ (due to low quality data or doublets), we summarized a total of 30 ITH programs across 13 tissues (Fig. 4A, Figure S3A), named as meta-programs (MPs). For each sample, we used AUCell score to represent the activity of each MP and then calculated the correlations between 30 expression programs and stemness (measured by CytoTRACE). The stronger correlations (R > = 0.7) were observed for MYC, Hypoxia, Cellcycle-G2M, Translation-initiation, EMT-related, Epithelial-senescence, Metabolism and Cellcycle-G0Arrest (Fig. 4B). More specifically, Cellcycle-G2M in the oral cavity and skin tissues, Epithelial-senescence (EpiSen) and Hypoxia in oral cavity, MYC in colorectum, Translation-initiation as well as Cellcycle-G0Arrest in thyroid were found to concurrently and highly correlated with our defined malignancy index and stemness (Fig. 4B, Figure S3B). This indicated that these six MPs were crucial for acquisition of stemness and tumorigenic properties during malignant transformation. We found that MPs may be expressed variably or uniformly in a given precancer/cancer sample across different tissue types, but the extent to which the expression of these critical programs varies along the process of tumorigenesis remains to be uncovered.
To better understand the mechanisms governing the generation of precancer cells as well as their malignant transformation, we sought to unravel the distinctive features acquired or lost during the progression from benign to malignant stages. We first utilized Leiden algorithm 19 to cluster all diseased epithelial cells into multiple subclusters (e.g., 24 clusters identified in liver-diseased tissues, Figure S3C). Next, we employed RNA velocity analysis to delve into those subclusters of diseased epithelia and infer epithelial cell transition trajectories (Fig. 4C). However, when we combined data from multistage epithelial cells (e.g., HCCs, NAFLD, and cirrhosis), ITH increased significantly as liver epithelia transitioned from precancerous to cancerous states (as shown in Fig. 4C). It suggested that stem-like cells in HCC may arise from cirrhotic stem cells, subsequently differentiated into HCC hepatocytes. Malignant hepatocytes also appeared to acquire a degree of stemness and engaged in dedifferentiation into a stem cell-like state (Fig. 4C). We thus focused on the critical cell subpopulations (i.e., stem-like cells) involved in the malignant transformation of 13 diseased tissues. We turned to performing Leiden sub-clustering analysis using stem-cell expression data from precancer and cancer samples, thereby reducing the intrinsic complexity and heterogeneity of diseases by dissecting state transition in cancer cell-of-origin. The number of stem cell subpopulations identified in different tissues varied. For example, 20 stem-cell subclusters (ranging from C0 to C19) were identified in the liver (as shown in Fig. 4D). To gain further insights into the transition trajectories within these identified subclusters, we conducted a subsequent RNA velocity analysis to elucidate their directionality embedded on a diffusion map (Fig. 4D).
Similar analyses were carried out across 13 tissues, identifying one to three trajectories within each (Fig. 4E, Figure S3C). For instance, two pronounced directional flows were observed from cirrhotic stem-like cells (C1 and C11) towards HCC subclusters (C0 and C7), respectively (Fig. 4E). Along the C1 to C0 transition path (Trajectory 1, Fig. 4E), notable increases were observed in cellular senescence, damage response, MHC I & II processing and presentation, fetal signatures, and cell cycle G0 arrest signature scores (Fig. 4F). Interestingly, a notable uptick in mutational load within C0 compared with C1 was noted (Fig. 4G). This observation is consistent with that MHC-IIhigh stem cells represent a deeper quiescent state and are more resistant to stress-induced proliferation than MHC-IIlow stem cells 20,21. These findings underscore the pivotal role of cellular senescence in tumorigenesis and suggest that mutation-accumulated stem cells in C0 gained a clonal advantage during aging by upregulating their surface expression of MHC class II molecules. In addition, RNA velocity analysis showed cirrhotic stem cells in C11 exhibited another directional flow toward C15 and C7 (Fig. 4E). Along C11 to C7 (Trajectory 2), decreased cellular senescence, fetal and quiescence signature but increased S and G2M cell cycle scores were observed (Fig. 4F), indicating restoration of proliferative capacity. Escape from proliferation arrest had been reported in certain circumstances, such as reactivation of telomerase activity or deletion of CDKN2A 22. Genomic events regarding C7 HCC cells harboring CDKN2A mutations and displaying the highest mutation load further supported the resumption of proliferation of stem-like cells in a quiescence-(or G0)-like state (Fig. 4F, G).
Next, we shifted our attention from stem-like cell clusters to clinical implications within them. Using our integrated scRNA-seq profile as a reference, we profiled the cellular components of each sample in the corresponding TCGA cohorts using CIBERSORTx 23. The results showed that the cellular composition inferred by deconvolution correlated with clinical features of TCGA patients. Particularly, the percentage of several transitional populations showed significant difference in tumor and adjacent normal tissues and notable associations with clinical features, mutations, and TME components of TCGA tumors (Figure S3D, Table S4). These inferred stem-cell proportions could also predict the OS or DFS of TCGA patients, such as C9 and C14 in breast, C1 in colorectum, C0 and C15 in liver, etc. (Fig. 4H). Compared with healthy stem cells, several transitioning clusters from precancer tissues were enriched in signaling pathways that regulate inflammation initiation, such as NF-kB signaling, JAK-STAT, toll-like receptor (TLR) and mitogen-activated protein kinase (MAPK) pathways (Fig. 4I, Table S5). This implicated the crosstalk of precancer stem cells between tumorigenesis and inflammatory processes. Most intriguingly, EpiSen was inevitably enriched in critical transitional stem cell populations involved in malignant transformation across 13 tissues (Fig. 4J), which was reminiscent of senescence-associated stemness. In addition to EpiSen, many malignant stem-like cell clusters exhibited strong enrichment of EMT, Hypoxia, TGF-β pathway and Autophagy (e.g., C14 in the cervix, C0 in the liver, C14 in the lung, C4 in the pancreas and C5 in the thyroid, Table S5), indicating that hypoxia signaling appeared to involve in the activation of stem cell stress signaling by expressing some stemness regulation genes, such as TGF-β, to induce dormancy, stimulate autophagy, promote cell survival and maintenance of stem cell identity, as well as promote EMT 24. The inferred fractions of these dormant stem cells were significantly correlated with higher TGF-β response in TCGA bulk data, further validating the observations from our single-cell analyses (Figure S3E). Nevertheless, the mere detection of a senescence-like cellular state, e.g., in BRCA1+/mut breast tissues (C14), colon adenomas (C17) and other malignancies, does not indicate whether these cells have a tumor-suppressive or pro-tumorigenic functional role.
5. Alterations in cell-cell interactions over time and master pro-tumorigenic mediators.
To comprehensively compare and quantify the differences between cellular states in malignant transformation, we identified DEGs between each pair of core clusters from each trajectory, such as C4 v.s. C9, C9 v.s. C14 in breast, C0 v.s. C1 and C7 v.s C15 in the liver (avg_log2FC > = 0.25, pct.1 > = 0.25, and pct.ratio > = 1.5). In addition, we endeavored to detect the transition genes by calculating Spearman’s rank correlation between the inferred pseudo-time and gene expression of all stem-like cells along each trajectory (with FDR < 0.05). Subsequently, based on the consensus of related MPs, we generated gene signatures that specifically reflect characteristics of different stem cellular states. Spearman correlation analyses were performed on gene expression levels and AUCell scores of each MP. Genes significantly correlated with AUCell scores (FDR < 0.05) were considered meta-program genes. To obtain specific gene signatures involved in malignant transition, we identified (1) DEGs that distinguish the diseased stem-like cells from normal-like cells (Pseudo-Bulk – DEGs, FDR < 0.05 and log2FC > 0.25), (2) DEGs between stem-cell subcluster, (3) transition genes, and (4) MP genes, and then intersected them to generate malignant transformation-related (PCT) gene signatures (Fig. 5A). Subsequently, the geometric mean of the Spearman correlation coefficients for each gene from 13 tissues were calculated. Finally, 100 top-ranked genes sorted by geometric mean value were filtered into PCT genes for each MP (Fig. 5A).
To seek for gene(s) in 30 core sets of PCT signatures that may drive stem-cell inflammation and aging 6 during loss of tissue homeostasis and precancer development, we then analyzed the dynamics of stem cell crosstalk with the microenvironment along our identified pseudotime axis. We identified 989 significant ligand-receptor pairs which revealed intercellular communication between 62 transitioning stem cell subpopulations and corresponding immune and stromal cells using CellChat 25 (Table S6). Several predicted interactions related to cytokines/chemokines in inflammatory response were frequently observed in epithelial stem cells and other PME cells during the transition from tissue homeostasis to pathological lesions. For instance, we observed extensive crosstalk of TNF-TNFRSF1A, CCL2-ACKR1, CXCL2-ACKR1, and CXCL12-CXCR4 in Breast-C14, Cervix-C5, Liver-C15, and Prostate-C14, suggesting their role in proinflammation and leukocyte recruitment (Fig. 5B, Figure S4A). This contribution likely leads to macrophage M1 polarization while inhibiting the invasion and migration of damaged epithelial cells 26. It’s worth noting that these ligand-receptor interactions involved in wound repair after tissue damage (like CXCL12-CXCR4 axis), potentially exert a two-edged effect in regulating anti-tumor immune responses (e.g., Prostate-C14 and Oral cavity-C18) and promoting tumor cell growth (e.g., Prostate-C34 and Liver-C0, Fig. 5B, Figure S4A) 27. We also noticed unique cellular crosstalk at different stages. Strikingly, interactions between the TIGIT receptor in lymphocytes (e.g., TREG and CD8TEREX) and NECTIN2/3 ligands in dendritic, mesenchymal and endothelial cells were found at C4 of breast cancer, but not at C9 and C15 stages of preneoplastic BRCA1+/mut breast tissues, showing the multiple immunosuppressive signals in the Breast-C4 TME. Additionally, FIB, INMON, and DC expressed TGFB1 to communicate with ECM-receptor, indicating that TGF-β signaling of C4 tumors may regulate fibroblasts activity and modulate ECM production and components (Fig. 5C, Figure S4B). Furthermore, interaction between stem cells and ICAF, INCAF, VFIB, PERI, HSC, INFIB or END via FN1-ITGAV + ITGB1 and FN1-SDC1 were detected in Cervix-C8, Colon-C1, Liver-C7, Oral calvity-C3 and Prostate-C34 (Fig. 5C). This indicates that extensive stromal remodeling occurs during tumorigenesis in infected or injured tissues, acquiring EMT-like features (e.g., FN1) in transitioning stem-like cells and showing high expression levels of SDC1 in fibroblasts, which is linked to aggressive phenotypes and poor survival (Fig. 4I, Fig. 5D, Figure S4C). FN1 from the Stress and Platelet-activation program was highly expressed in those ‘hybrid’ cells undergoing plasticity changes, but somewhat less in other transitional stem clusters (Fig. 5A, Figure S4C). This finding suggests that during physiochemical stress situations, pro-tumorigenic fibroblasts in close proximity can provide the fertile ‘soil’ to the cancer ‘seed’, further influence platelet activation, which produced pro-angiogenic and growth factors that facilitate tumor growth and survival, as well as promoting the metastatic potential of tumor cells 28. In addition, dynamic shift of HSC/VFIB to SMC/MYOFIB that highly expressed ACTA2 and MYH11 was observed in the liver, colon, and prostate diseased tissues (Fig. 5E), consistent with the aforementioned accumulation of MYOFIB and HSC in precancerous stage (Fig. 3E). These observations were also in line with that type I collagen, enriched in activated myofibroblastic HSC, promoted proliferation, tumor development and related to elevated HCC risk in patients 29.
Notably, several ligand-receptor pairs such as LGALS9-CD44/CD45, laminins-integrins and MIF-CD74 + CXCR4 were significantly enriched in the communication between transitional stem cells and myeloid subsets (e.g., NEUT, INMON, M2MAC, etc., Fig. 5F). The MIF-CD74 + CXCR4 interaction can exert an immunosuppressive role by affecting downstream MAPK signaling pathway effectors 30. Most interestingly, ANXA1-FPR1/2 of the ANNEXIN pathway was activated in the malignant transition of 11 tissues (Fig. 5F). ANXA1, one of the representative genes in four stem-cell MPs (i.e., EpiSen, Stress, Unfolded-protein-response and Interferon/MHC-II), is highly expressed in most healthy human tissues, but silenced in nonalcoholic steatohepatitis 31 and early onset of esophageal and prostate carcinoma 32. The expression levels of ANXA1 in stem cells decreased from the healthy to the precancerous stage but gradually increased along our identified malignant transformation paths (Fig. 5G, Figure S4D). The expression of ANXA1 receptor genes, including FPR1, FPR2, and FPR3, was enriched in myeloid cells and elevated in most malignant stages (e.g., M2MAC and INMON, Fig. 5H, Figure S4D). Additionally, we observed a pronounced elevation in the inflammatory response scores of transitioning stem subpopulations at the starting point of trajectories when ANXA1 expression was depressed during tissue damage (Fig. 5I). Conversely, as malignancy progression, positive correlations were noted between higher ANXA1 expression levels and the GSVA score of the TGF-β signaling, implying the involvement of ANXA1 in promoting TGF-β activation from epithelial stem cells (Fig. 5J, Figure S4E).
6. Spatiotemporal ANXA1 expression patterns in heterogeneous tumor ecosystems facilitate immunosuppression of the microenvironment.
To further investigate the spatial architecture of transitioning stem subpopulations and distinct TME clusters during disease progression, we analyzed the spatial transcriptome (ST) data from breast, cervix and liver-diseased tissue samples. The cell components in each spot were determined by SpaCET using our integrated scRNA-seq data as the reference. We also deconvoluted the cell-type composition of each region in one HPV-negative normal cervix sample based on scRNA-seq data from the healthy cervix. The dynamic expression levels of ANXA1 were confirmed by ST data analysis, demonstrating a substantial increase in the bulk levels of ANXA1 RNA (located in both healthy epithelium and tumor spots) along malignant transition routes (Fig. 6A). Since ST data has not yet reached single-cell resolution, ANXA1 also exhibited strong expression in myeloid and mesenchymal spots of cancers and elevated in advanced stages (Figure S5A). It is still challenging to validate several transitional stem clusters with ST data (e.g., Liver-C0). Under these circumstances, tumor epithelial regions may have both stem-like and CAF features. As expected, activation of the ANNEXIN pathway and Retinoic Acid (RA) were observed among stem-like, myeloid, and mesenchymal subclusters, as inferred through spatially-proximal cell-cell communication analysis using CellChat V2 (Fig. 6B, Figure S5B). These findings are consistent with previous studies highlighting the role of CAF-secreted ANXA1 in imparting stem cell-like properties to cancer cells 33. Additionally, research has shown that the N-terminal peptide of AnxA1 initiates a signaling cascade through FPR2, leading to the polarization of macrophages towards the M2 phenotype 34. Subsequently, we aimed to delve into the spatial heterogeneity of tumors in more depth. Our results showed that the spots of FIB-INMON interactions with co-expression of ANXA1-FPR3 and TGFB1-TGFBR1/2 in the breast tumor sample are close to the C9 but distant to the C4 population (Fig. 6C, D, Figure S5C). Two different cancer cell states in ST data displayed distinct enrichment pathways, consistent with GSVA analysis of scRNA-seq data. Proximal C9 cells clustered separately from other tumor cells and showed senescent phenotype, while distal C4 spots exhibited hypermetabolic features and higher activity of oncogenic pathways, such as KRAS and EGFR signaling as well as DAP12 interactions (with myeloid cells) (Fig. 6E). We speculate that the senescent state of C9 cells in breast tumors is triggered by hypoxia and dysregulation of microenvironmental growth factors like TGF-β and NOTCH signaling. This, in turn, enhances ANXA1-FPR1/2/3 crosstalk, facilitating the reprogramming of fibroblasts and promoting the shift into an immunosuppressive environment from the PME. Last, we examined the bulk levels of ANXA1 RNA in TCGA data and observed a significant correlation of ANXA1 with CAF activation marker genes (FAP, ACTA2, COL6A3, MMP1 and TGF-β), wound healing-myCAF, ecm-myCAF, and TGFβ-myCAF gene signatures, as well as M2-like signatures of tumor-associated macrophages (TAMs) (Fig. 6F). This further confirmed that ANXA1 plays an important role in inducing immunosuppressive TME by mediating the transformation of normal fibroblasts into CAFs and macrophages M2 polarization. Additionally, our defined EpiSen signature was also significantly and positively associated with stromal cell fractions and TGF-β response in TCGA tumors (Figure S5D).