Integrating spatial and single-cell transcriptomics reveals tumor heterogeneity and intercellular networks in colorectal cancer

Integrating spatial and single-cell transcriptomics reveals tumor heterogeneity and intercellular networks in colorectal cancer | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Article Integrating spatial and single-cell transcriptomics reveals tumor heterogeneity and intercellular networks in colorectal cancer Ligong Lu, Jing Xiao, Xinyang Yu, Fanlin Meng, Yuncong Zhang, and 7 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-3820516/v1 This work is licensed under a CC BY 4.0 License Status: Under Review Version 1 posted 9 You are reading this latest preprint version Abstract Single cell RNA sequencing (scRNA-seq), a powerful tool for studying the tumor microenvironment (TME), does not preserve/provide spatial information on tissue morphology and cellular interactions. To understand the crosstalk between diverse cellular components in proximity in the TME, we performed scRNA-seq coupled with spatial transcriptomic (ST) assay to profile 41,700 cells from three colorectal cancer（CRC） tumor-normal-blood pairs. Standalone scRNA-seq analyses revealed eight major cell populations, including B cells, T cells, Monocytes, NK cells, Epithelial cells, Fibroblasts, Mast cells, Endothelial cells. After the identification of malignant cells from epithelial cells, we observed seven subtypes of malignant cells that reflect heterogeneous status in tumor, including tumor_CAV1, tumor_ATF3_JUN|FOS, tumor_ZEB2, tumor_VIM, tumor_WSB1, tumor_LXN, and tumor_PGM1. By transferring the cellular annotations obtained by scRNA-seq to ST spots, we annotated four regions in a cryosection from CRC patients, including tumor, stroma, immune infiltration, and colon epithelium regions. Furthermore, we observed intensive intercellular interactions between stroma and tumor regions which were extremely proximal in the cryosection. Specially, one pair of ligands and receptors (C5AR1 and RPS19) was inferred to play key roles in the crosstalk of stroma and tumor regions. For the tumor region, a typical feature of TMSB4X -high expression was identified, which could be a potential marker of CRC. The stroma region was found to be characterized by VIM -high expression, suggesting it fostered a stromal niche in the TME. Collectively, single cell and spatial analysis in our study reveal the tumor heterogeneity and molecular interactions in CRC TME, which provides insights into the mechanisms underlying CRC progression and may contribute to the development of anticancer therapies targeting on non-tumor components, such as the extracellular matrix (ECM) in CRC. The typical genes we identified may facilitate to new molecular subtypes of CRC. Biological sciences/Cancer/Cancer microenvironment Biological sciences/Cell biology/Cell adhesion Figures Figure 1 Figure 2 Figure 3 Figure 4 Introduction Colorectal cancer (CRC) is a molecularly heterogeneous disease [ 1 , 2 ]. The heterogeneity of cell types involved in CRC carcinogenesis makes it difficult to elucidate cell lineages using traditional developmental biology techniques such as bulk transcriptomics methods [ 3 ]. Through single-cell transcriptomics technology, it is now possible to deconstruct a tumor into its diverse cell subpopulations and thus gain a better understanding of the underlying biology like subtyping [ 4 – 6 ]. However, spatial or anatomical information inherent in the tissue architecture is lost using single-cell transcriptomic technology only. Spatial transcriptomics (ST) is an emerging technology that adds spatial dimensionality and tissue morphology information to the single-cell transcriptomics data of cells in an undissociated tissue, thus helping to preserve precise spatial or anatomical information. Overcoming the throughput limitation of in situ hybridization (ISH) methods, ST allows for unbiased mapping of transcripts in individual tissue sections with spatial resolution by using spatially barcoded oligo-deoxythymidine microarrays [ 7 ]. As a high-throughput spatially resolved transcriptomics tool, ST has been used to study architecturally complex tissues and diseases including melanoma [ 8 ], prostate cancer [ 9 ], gingival tissue [ 10 ], adult human heart tissue [ 11 ], mouse and human spinal cord tissue [ 12 ], model plant species [ 13 ], as well as their spatiotemporal characterizations during development [ 14 ] [ 15 ]. Extensive multimodal studies have unraveled molecular landscape of diverse diseases [ 16 ]. Combining these two complementary and powerful technologies has been confirmed to be scalable to study architecturally complex tissues and to provide meaningful biological insight across a range of pathologies, such as bone marrow [ 17 ], pancreatic ductal adenocarcinomas [ 18 ], lung fbroblasts [ 19 ], and human dorsolateral prefrontal cortex [ 20 ]. The tumor microenvironment (TME) comprises various cell types (immune cells, fibroblasts, endothelial cells, etc.) and extracellular components (growth factors, cytokines, extracellular matrix, hormones, etc.) that surround cancerous/tumor cells [ 21 ]. Since many currently used anticancer therapies target non-tumor components, such as the extracellular matrix (ECM), immune system and vascular system [ 22 ], understanding cellular components and how their dynamic interactions to shape the tumor landscape are particularly important. In this study, we aim to provide a comprehensive global view of tumor heterogeneity and intercellular interaction networks of CRC using single-cell transcriptional profiles coupled with spatial transcriptional profiles. By analyzing the single-cell and spatial transcriptional profiles of 41,700 cells from 3 treatment-naïve patients with CRC, we generated a molecular map of all major CRC populations based on single-cell RNA sequencing (scRNA-seq). The malignant cells in epithelial cells were identified and categorized into seven subclasses (tumor_CAV1, tumor_ATF3_JUN|FOS, tumor_ZEB2, tumor_VIM, tumor_WSB1, tumor_LXN, tumor_PGM1), which may help to the molecular subtyping of colorectal cancers. In addition, we used spatially resolved transcriptomics in combination with computational tools to attribute cell types to different CRC niches. Annotated tumor regions based on the cryosection sections represented high TMSB4X expression, and suggested a typical marker of tumorgenesis. The stroma region was characterized by VIM gene, which was also used as a typical feature of one subtype of malignant cells in CRC scRNA-seq. Furthermore, we inferred the important interaction between tumor and stroma regions mediated by gene pair of C5AR1 and RPS19, which played roles of ligand and receptor, respectively. Results Landscape view of cell composition in tumors, adjacent tissues and peripheral blood in patients with CRC To shed light on the complexity of the TME in CRC, we performed scRNA-seq along with spatial trancriptome sequencing on viable cells derived from matched tumor and adjacent tissues, as well as peripheral blood mononuclear cells (PBMCs) of 3 patients with CRC (Fig. 1 a, Supplementary Table S1 ). On average, we obtained more than 150 G sequencing reads for each sample, with a median sequencing saturation of 91.40% (87.0%-95.5%). After quality control filters (few detected features in cells and few expressed cells associated with detected features), a total of 35,666 cells were identified in 9 samples derived from 3 patients (including 10347, 13241 and 18112 in tumor tissues, adjacent tissues and peripheral blood, respectively; Supplementary Table S2 ). We obtained approximately 1000 genes and 2500 unique molecular identifiers (UMIs) for each cell, indicating sufficient coverage and transcript representations. To define cell clusters with similar expression profiles, we performed dimensionality reduction of t-distributed stochastic neighbor embedding (tSNE) implemented in the Seurat package. Each cluster was further identified as a specific cell subpopulation on the basis of the expression of the most variable genes and the canonical markers, including those in epithelial cells (with gene markers: EPCAM, KRT5, PHGR1, LGALS4, and TFF3), T cells (CD4 + T cells: PTPRC, CD3D, and CD4; CD8 + T cells: PTPRC, CD3D, and CD8A), B cells (CD19 and MS4A1), monocytes (CD14, ITGAX for CD11C), natural killer (NK) cells (FCGR3A and NCAM1), endothelial cells (CDH5, PLVAP, CLDN5, VWF), fibroblasts (LUM, DCN, COL1A1), and mast cells (KIT, CPA3, MS4A2, and TPSAB1) (Fig. 1 b). In addition to these well-known markers, we also analyzed cluster-specific genes via differential gene expression analysis (Supplementary Table S3 ). These cluster-specific marker genes included FBLN1 for fibroblasts, as well as MT1A and PLN for smooth muscle cells (Fig. 1 c, Fig S1 e ). In total, eight cell types in CRC were identified based on canonical markers and cluster-specific genes: epithelial cells, fibroblasts, endothelial cells, monocytes, T cells, NK cells, B cells, and mast cells. The heterogeneous compositions of the TME in CRC across tumor tissues, normal tissues and peripheral blood are consistent with a recent single-cell transcriptome study of CRC [ 23 ]. To characterize different cell compositions in tumor tissues, normal tissues and peripheral blood in CRC, the proportions of each cell type were investigated. An overall increase in myeloid cell populations and decrease in B cell populations were observed in tumor tissues compared to normal tissues (Fig. 1 b, bottom; Fig S1 d), suggesting a redirected immune response in CRC patients. In details, we observed that the proportion of monocytes was increased with approximately 2.5-fold, whereas that of NK cells and B cells was decreased (about 0.3–0.4 times) in tumors compared to normal tissues, suggesting a myeloid immunosuppression in the CRC TME (Fig S1 a,b, Supplementary Table S4 ). To further explore the distinct cell compositions in the TME across individuals, more detailed proportions were assessed (Fig. 1 d). These results showed, for example, that in patient T0602, the proportion of epithelial cells decreased in contrast with patient T0529 and increased compared to that in patient T0609 (Fig. 1 d, left; Fig S1 c1; Supplementary Table S5 ). Since the transition from normal epithelium to intraepithelial neoplasia were found to be associated with CRC patient survival [ 24 ], the difference in epithelial cells across individuals may be important for survival and worthy of further investigation. Considering that the cellular proportion determined by scRNA-seq may be biased toward an underrepresentation of malignant cells derived from epithelial cells [ 25 ], we also explored the proportions of immune and stromal cells to account for all cells except epithelial cells which includes tumor cells. The results showed that myeloid cell-driven immune response in patient T0529 was stronger than that in the other two patients (Fig. 1 d, right; Fig S1 c2 ). Epithelial cells represents multilineages including a lineage of malignant cells It has been suggested that human colon cancer cells recapitulate the multilineage differentiation processes of normal colon epithelia. To investigate each lineage contributing to the CRC heterogeneity at single cell resolution, we subclustered cell populations for each cell type to identify subpopulations. To annotate these subpopulations, we combined another published CRC cohort consisting of 6 CRC patients in tumor regions as well as matched normal mucosa [ 23 ], and transferred the annotations of subtypes to our datasets in this study with the Seurat R package (Fig. 2 a). Since the transition from normal epithelium to intraepithelial neoplasia were found to be associated with CRC patient survival [ 24 ], we focused on epithelial cells and found 9 subpopulations, namely CD19 + CD20 + B cells, crypt cells, enterocytes, goblet cells, intermediate, mature colonotypes, proinflammatory, stem-like, and tumor cells (Fig. 2 b). To distinguish malignant cells and nonmalignant cells in epithelial cells, we performed scRNA-seq-based copy number variation (CNV) and subclustering analysis (Fig S2 a, b). The proportions of malignant cells in each subcluster of epithelial cells were shown in Fig S2 c. The trajectory revealed a transcriptional hierarchy, defining seven molecular states (Fig. 2 c, top). The cells from tumor tissues dominated the divergent differentiation states 2 and 5, suggesting the tissue arrangement along pseudotimes (Fig. 2 c, bottom). To illustrate the differentiation paths across the multilineages among the epithelial cell populations, the semisupervised trajectory inferred by monocle2 [ 26 ] revealed a transcriptional hierarchy defining three branches. The hierarchy was dominated by malignant epithelial cells, as well as normal epithelial cells (including goblet cells and (stem-like/transit amplifying cells) and immune-related cell types (including proinflammatory and mature colonotypes), which originated from normal epithelial cells with branching toward malignant epithelial cells (grey) and immune-related cell types (light green), respectively (Fig. 2 d). Projection of malignant epithelial cells along the epithelial cell differentiation trajectory revealed segregation of tumor cells from normal epithelial cell types and stem-like populations. The greater stemness of malignant epithelial cells suggested the regenerative/proliferative potential of these tumor cells (Fig. 2 e). The hypoxia and epithelial mesenchymal transition (EMT) were also investigated in the malignant epithelial cell populations (Fig S2 e). Transcriptional and functional features of malignant cells reveal heterogeneity in CRC patients To characterize the malignant cell populations, we scrutinized the transcriptional features between malignant and nonmalignant cells. The known malignant epithelial cell populations characterized by upregulated expression of S100A4, VEGFA, MYC, and ICAM1 (intercellular adhesion molecule-1), according to their significant differential expression (log e |fold change| > 0.25, p value < 0.05) (Fig. 2 f, left). The most differentially expressed gene EMP3 (Epithelial membrane protein 3), which has been identified as an oncogene and therapeutic target in breast cancer [ 27 , 28 ], remains to be elucidated in colon cancer. The characteristic genes in malignant epithelial cell populations were found to be involved in biological processes such as matrix remodelling, cell proliferation and apoptosis (Fig. 2 f, right, Supplementary Table S6), hinting the occurrence of EMT. Moreover, the terms ‘positive regulation of angiogenesis’, ‘cellular response to decreased oxygen levels’ and ‘extracellular matrix organization’ were enriched with differentially expressed genes, also suggestive of the malignant tendency of the cell populations. Using genes that characterizing malignant epithelial cell populations, the biological pathways that those genes implicated in were shown in Fig S2 d. In addition, we predicted the regulons for malignant and non-malignant cells, respectively. The list of top 5 regulons for the two cell populations were shown in the heatmap (Fig. 2 g). In malignant cells, a famous oncogene MYC was shown to be one of the key regulons. The transcriptional regulation role of ATF1 in CRC cell lines has been characterized by a study of combing RNA-seq and ChIP-seq assays, in which found rs7017386 allele-specifcally enhanced the binding affnity of ATF1 and promotion of two oncogenic lncRNAs via forming a long-range chromatin loop [ 29 ]. Therefore, both the characteristics of the Stemness, hypoxia, and EMT and the prediction of oncogenic regulons consistently reflected the malignancy and tumorgenesis roles of the identified malignant epithelial cells aformentioned. To focus on transcriptional programs for subcategorizing tumor cells, we performed subclustering and trajectory analysis for malignant cell populations. We assigned seven tumor cell subclusters to all the malignant cell populations, namely tumor_CAV1, tumor_ATF3_JUN|FOS, tumor_ZEB2, tumor_VIM, tumor_WSB1, tumor_LXN, and tumor_PGM1 (Fig. 2 h, Fig S2 f). The trajectory revealed a transcriptional hierarchy, converging into three discrete tumor subclusters. One subcluster was highly enriched with gene response to histone deacetylase (HDAC) inhibitors (like ATF3 and CAV1 ), and another subcluster was highly enriched with inflammatory gene ( LXN and PGM1 ). The subcluster enriched with tumor metastasis-related hypoxia ( WSB1 ) originated from the aforementioned differentiation paths. Sample arrangement along the differentiation trajectory of tumor cells revealed malignancy in patient T0602. Transcriptional features in a spatial resolution on cryosections from patients with CRC To examine spatially transcriptional differences within colorectal tissues, we mounted cryosections of unfixed CRC tissues originating from the same CRC patients onto spatially barcoded ST microarray slides to generate unbiased transcriptome maps. After haematoxylin and eosin (H&E) staining and brightfield imaging, we annotated the slides according to the distinct histological features (Fig. 3 a). The samples were then processed for ST analysis. We demultiplexed the sequenced reads and identified their spatial location within tissues using the ST location-specific barcodes of the array. For patient T0602, we detected approximately/average 3,413 median UMIs and approximately/average 1,660 median gene numbers per ST spot for both the tumor tissue section (named CRC5_1) and the normal tissue section (named CRCN5_1). First, the spatial transcriptomics data were integrated with the scRNA-seq data using Seurat-v3 anchor-based integration to annotate each region in the corresponding section [ 30 , 31 ]. Every spot in the spatial data was considered a weighted mix of cell-types identified by scRNA-seq. For each spot, the cell type with the maximum prediction score among all possible cell types and thus transferred from the scRNA-seq dataset is illustrated (Fig. 3 b, f). After further adjustment on the basis of annotated histological features, we annotated four and two anatomical regions in the CRC5_1 section (derived from a tumor tissue, Fig. 3 c), and in the CRCN5_1 section (derived from an adjacent tissue), separately (Fig. 3 g). We observed many obviously characteristic genes, which represented higher expression in annotated regions especially in the tumor tissues (Fig. 3 i) compared to normal tissues (Fig. 3 j). It is noted that five DEGs in comparision of malignant and non-malignant cells in the CRC scRNA-seq dataset were included IFITM1, CXCL1, CXCL8, S100A4, and TGFBI. The higher expression in tumor or stroma regions were shown in Fig S3 . IFITM1 was highly expressed and spatially restricted relative to the annotated tumor regions. IFITM1 is an interferon-induced transmembrane protein family member. Fang et al. investigated that over-expression of IFITM1 promoted the aggressiveness of CRC cells, whereas knockdown of IFITM1 expression inhibited cell migration, invasion or tumorigenicity in vitro [ 32 ]. Pauline et al. found a highly significant increase in IFITM mRNA levels in 154 patients with colon and rectal carcinomas, compared to corresponding normal tissues [ 33 ]. CXCL1 and CXCL8 , members of the angiogenic CXC chemokine family [ 34 ] are highly expressed only on the bjojobsorder of annotated tumor regions, suggesting that an immune exclusion and blockage by the tumor. S100A4 (S100 calcium-binding protein A4) was a typical feature in the annotated stromal region and it has been confirmed to be associated with metastatic and recurrent CRCs [ 35 , 36 ]. TGFB1 is highly expressed and spatially restricted relative to the annotated stromal region. TGFβ-activated stroma was one of the key hallmarks of microsatellite-stable colorectal cancers displayed by a quadruple-mutant mouse model [ 36 ]. The TGF-β activity in stromal cells functionally increased the efficiency of organ colonization by CRC cells. In the oral squamous cell carcinoma, THBS1 expression was induced by TGFB1 in the stroma, and the dependency on the TGF-β stromal programme for metastasis initiation may be exploited to improve the diagnosis and treatment of patients of CRC [ 37 ]. In addition, the genes characterizing classical phenotypes of cancer stem-like ( ASCL2 ), hypoxia ( WSB1 ) and apoptosis ( ATF3 ) were all highly expressed in the tumor or stromal region, suggesting cell differentiation programs. To facilitate comparisons, the expression of the aforementioned genes in the adjacent normal tissue samples is presented (Fig. 3 j). Different anatomical regions in a spatial resolution represents cell compositions on tumor, stroma, immune and epithelium Second, standalone analysis of the spatial transcriptome along with the annotated anatomical regions was performed to identify spatially differential transcriptional programs. Interestingly, we found region-specific transcriptional differences. Unsupervised clustering of spatial transcriptome of CRC5_1 along with annotated anatomical regions revealed that the tumor region was characterized by heterogeneous cell cluster enriched high TMSB4X expressing cell populations (Fig. 3 d, Fig S3 ), which was not observed in CRCN5_1 (Fig. 3 f). TMSB4X encodes thymosin beta-4, a well-known secreted small peptide, which plays key roles in the cell motility, cytoskeletal organization, cell survival, and cell differentiation [ 38 ]. To the best of our knowledge, little is known about the spatiotemporal expression profile of TMSB4X in CRC. We observed a tumor region-restriction feature within high TMSB4X expressing cell populations (Fig. 3 d). TMSB4X has been validated as a therapeutic target in colon cancer stem cells in a previous study [ 39 ]. Additionally, in other architecture-dependent tissues whose spatial locations are deeply intertwined with their functions, such as brain and heart, TMSB4X has been reported to be involved in tumor progression via neovascularization, cell adhesion and the epithelial-mesenchymal transition [ 39 – 41 ]. When TMSB4X was silenced both in vitro and in vivo , differentiation and tumorigenicity were diminished [ 40 ]. Targeting highly expressed TMSB4X or TMSB4X-high cell populations identified in spatial trascriptome in CRC has the potential to become a new therapeutic strategy in CRC. Moreover, the annotated stromal region overexpressed classical marker genes for fibroblasts or endothelial cells, such as LUM , VIM , COL1A1 , and COL1A2 . It is noted that the differential expression of VIM gene in comparison of stroma regions and other regions in CRC ST-seq dataset was also identified in CRC scRNA-seq dataset,in which VIM gene was one of the DEGs in malignant cells compared to non-malignant cells. As proximity is a necessity for physical interactions among cells, anatomical regions or cell-type proximity/interactive maps can be used to guide the discovery of interactions between anatomical regions or cell types in the same or different lineages. We first carried out a pseudotemporal trajectory analysis of the four anatomical regions in CRC5_1 section. Our data revealed that the lineage stemmed from the colon epithelium, and then went through an infiltration to divide into two major lineages, corresponding to stroma and tumor spots, respectively (Fig. 4 a, top). The trajectory reconstruction confirmed two major terminal/branching cell fates at tumor and stoma spots, respectively (Fig S4 a). Moreover, seven continuous states were identified during along pseudotime trajectory (Fig S4 b). Tumor spots were observed mainly in state 5, and stroma spots were found primarily in state 6, which indicated different states during tumorgenesis, even though these two anatomical regions were proximal to each other. The colon epithelium was observed mostly in state 1, and immune infiltration was distributed in states 2,3,4 and 7, sporadically in state 6. Consistent with the trajectory, the tumor and stroma spots exhibited greatest pseudotime meaning the most extent of differentiation and mature/terminal programme (Fig S4 b). To explore the characteristics of the branching cell fates, cell cycle phase analysis was performed. Approximately 70% of stroma spots arrested in the G1 phase of cell cycle, but this outcome was not observed in tumor anatomical regions (Fig. 4 b). G1 phase cell cycle arrest may be responsible for the inhibition of colorectal cell proliferation in the stroma region. Compounds that promote G1 cell cycle arrest were used and confirmed to be a treatment of colorectal cancer [ 42 ]. Spatially resolved interactions of tumor and stroma regions To investigate the regulators critical for the branching cell fates of tumor and stroma regions, SCENIC-based regulon analysis was performed. The results showed a spatial-resolved, specific regulon set (Fig S4 c). Four of the top 5 regulons in stroma or tumor region were overlapped, including CDX1 , IRF8 , HNF4A ,and CREB3lL1 (Fig S4 d). The stroma region was specifically regulated by IRF3 . The previous research suggested that overexpression of IRF3 causes cell-cycle arrest in the G1/S phase thereby resulting in inhibition of DNA synthesis [ 43 , 44 ]. However, the tumor region was found to be specifically regulated by POU2F2 . A recent research suggested that POU2F2 played tumorigenic roles in glioblastoma by leading to a metabolic shift towards aerobic glycolysis [ 45 ], but the roles in colorectal cancer remain poorly understood. In conclusion, the transcription factor IRF3 plays a key role in the generation of stromal cells in the stroma region, but not tumor cells in the tumor region, and the underlying programme involves in inhibition of cell proliferations by arresting cells in the G1/S phase. To further investigate the interaction between the neighboring anatomical region, especially between the tumor region and other regions, we secondly performed CellphoneDB-based cell interaction analysis to investigate the underlying ligand-receptor pairs in different anatomical regions derived from CRC tumors. The aforementioned trajectory suggested the crosstalk between stroma region and tumor anatomical region (Fig. 4 a, top). We observed intensive cellular interactions between the stroma region and tumor regions. (Fig. 4 c). For example, the stroma region was predicted to interact with the tumor region in tumor tissues through C5AR1-RPS19, which is known to promote tumor growth by facilitating recruitment of these cells to tumors [ 46 ]. C5AR1 is known to activate and recruit myeloid-derived suppressor cells to tumors [ 47 ]. Accordingly, research studies have shown that the activated complement system has a tumor-promoting effect, including angiogenesis, trophoblastic invasion and tissue remodeling, which includes processes favorable for tumor establishment and progression [ 48 , 49 ]. RPS19 , which is released from apoptotic tumor cells and upregulated in colon carcinoma compared to normal epithelial colon cells, has been reported by its novel immunosuppressive properties [ 50 ]. To evaluate the refined characterization of the tumor region, we reclustered the spots in the tumor region to discern any spatial differentiation process. Reclustering the spots in the tumor region led to the identification of five subregions corresponding to five gene modules (Fig. 4 d): C0 spots expressed high levels of collagen, which is the major component of the TME and participates in cancer fibrosis [ 51 ], including collagen type III ( COL3A1 ), COL4A2 , COL4A1 , COL6A2 , COL1A2 , and COL1A1 . In contrast to C1 spots, C0 spots in the outermost layer in the CRC5_1 section, constituted the primary structural element of the ECM. C1 spots expressed high levels of ECM transcripts implicated in cell migration, including fibronectin (FN1), tumor protein, transnationally-controlled 1 (TPT1), transforming growth factor beta induced (TGFBI), and whey acidic protein (WAP) 4-disulfide core domain protein 2 (WFDC2).In addition, C1 spots maintain an intermediate phase with both undifferentiated and differentiated phenotypes. Transferring cell type labels to spatial data suggested that only the epithelial cells in scRNA-seq were spatially restricted to the outermost layer of the section (tumor subcluster C0), and C1 was spatially restricted to the centre of the section. C2 spots expressed extremely high levels of proline rich and Gla domain 3 (PRRG3), which is a member of a family of vitamin K-dependent transmembrane proteins that contain a glutamate-rich extracellular domain. C2 spots characterized by a senescence-like proliferation, are represented by MAP3K12. MAP3K12 is an upstream activator of the mitogen-activated protein kinase (MAPK) pathway [ 52 ]. C3 spots expressed extremely high levels of fatty acid-binding protein 1 (FABP1), which is essential for proper lipid metabolism in differentiated enterocytes [ 53 ]. The creatine kinase B (CKB) promotes metastatic survival by modulating intra- and extracellular energetics [ 54 ]. C3 spots are characterized by metabolism. C4 spots expressed high levels of mitochondrial phosphoenolpyruvate carboxykinase 1 (PCK1) gene, which increases colon cancer cell growth in part by promoting the consumption of both glucose and glutamine in the tricarboxylic acid (TCA) cycle [ 18 ]. The roles of hypoxia-reprogrammed TCA cycles in promoting human breast cancer cell growth via a HIF-1α-mediated PCK2 pathway have been reported [ 55 ]. C4 spots are characterized by hypoxia-response. Based on characteristic genes of five subregions of tumor region, we defined five gene modules, including focal adhesion dynamics (C0), intermediate (C1), ECM (C2), metabolic (C3), hypoxia-response (C4) modules. The results suggested roles for these gene modules in tumor progression, implying a need for a combination of an anticancer therapy with corresponding modulators. Discussion Colorectal cancer is a complex and heterogeneous malignant tumor of the colon and rectum. According to the degree of tumor differentiation and invasion, CRC can be classified into different subtypes. However, the published CRC studies usually were performed with whole tissues, such as bulk RNA-seq, blurring the heterogeneous characteristics of different cell types and limiting the ability to capture tumor heterogeneity. Cancerous tissue is composed of a mixture of various components, such as tumor cells, stromal cells, immune cells and ECM, leading to a complex TME [ 56 ]. The TME components exhibit interactive crosstalk with tumor cells and there surrounding factors, which in turn shapes tumor structure, metabolism and secretion, thus affects tumor development and/or metastasis. Immune cells within the TME play crucial roles during tumorgenesis. Immunotherapy aims to fight against cancer, infection and other diseases stimulating or suppressing the immune system. Immunotherapy displays promising therapeutic outcomes and limited side effects [ 57 ]. An increasing number of clinical trials have proven the effects of immunotherapy in certain types of solid tumors, such as melanoma, non-small-cell lung cancer, renal cancer and prostate cancer [ 58 ]. Recently, immunotherapy drugs such as CAR-T drugs have been approved by the FDA for clinical application. However, not all patients respond favorably to immunotherapy. Researchers have begun to examine the complexity and diversity of the TME and are realizing its importance in immunotherapy. Focusing on the TME facilitates to better understanding of the occurrence, development and metastasis of tumor and may lead to better diagnosis and treatment [ 59 ]. The expression pattern and function of tumor cell-associated immune molecules from in the TME provide useful information to determine whether a patient will might benefit from immunotherapy. There is an urgent need for improved techniques to better understand the TME and analyze the composition of immune cells and various other various cell types in tumor tissues. Single-cell sequencing and spatial transcriptomics have satisfied this requirement, and these two new-newly emerging techniques can be used to analyze samples at the single -cell level and monitor the in situ spatial information of tumor tissues. Single-cell RNA sequencing enables to investigate into the transcriptional regulation of highly heterogeneous cell populations or subpopulations and facilitates the discovery of genes that indicate cell subtypes, or that mark intermediate states during a biological process, as well as bifurcation between two alternative cellular fates. Spatial RNA sequencing enables the anatomical, in situ locations to be preserved, but cannot achieve rigorous single cell resolution. Actually, 10–20 cells are typically identified in each spot in the ST assay. Since in situ locations representing tissue sections are lacking in the scRNA-seq and single-cell resolution is lacking in the spatial RNA-seq, it is necessary to combine single-cell RNA-seq and spatial RNA-seq to reflect both cellular locations at true single-cell resolution. This allows inferences on functional relationships between scRNA-seq-defined populations based on their colocalization in space, and ultimately provides a more comprehensive characterization of cell types in their native environment than can be gained from either modality alone. In addition to transferring annotation at single cell level to spatial locations, we think that the current combinatory analysis between single-cell and spatial transcriptome remain to be more closely and complex. In this study, we combined single-cell and spatial transcriptomics to create a hierarchical map of cellular lineages in CRC. There CRC patients are involved, as well as single-cell transcriptomes of ~ 47,000 cells, and spatial gene expression maps. We firstly constructed CRC single-cell maps consist of epithelial cells, mast cells, monocytes, T cells, B cells, Endothelial cells, and NK cells. Then we focused on the epithelial cells, from which we identified malignant cells. Also, we performed analyses of differential expression, functional enrichments, transcription factors to characterize the features of malignant cells compared to non-malignant cells. By sub-classification of malignant cells, we found seven subtypes according to transcriptional features, which could be helpful to molecular subtyping of CRC. Combining spatial and scRNA-seq datasets using bioinformatics approaches, we transferred cell type annotations at single cell level to those at spot level. The results showed four regions including tumor, stroma, immune-infiltration, and colon epithelium. Then we mainly focused on proximate cellular interactions (tumor and stroma) within colorectal tissue, to quantify anatomically restricted gene expression and explore crosstalks between tumor and stroma regions. Tumor and stroma regions were characterized by TMSB4X and VIM high expressions, respectively. The cellular interactions or crosstalks were inferred to be mediated by C5AR1 and RPS19, which remained to be further validations. In summary, we comprehensively explored the cellular landscape and reconstructed the putative interaction network consisting of tumor cells and their microenvironment. This collective view allowed us to elucidate how diverse cellular components jointly determine CRC molecular subtypes in individual patients. Materials and Methods Subjects and clinical characteristics Patient samples were collected from 3 CRC patients (CRC0529, CRC0602, CRC0609), see detailed information in the Supplementary Table S7. Matched adjacent normal tissues and primary tumors and peripheral blood were obtained from all 3 patients. All patients were treatment-naive before tumour resection. All sampling and experimental steps in this study were approved by the Ethics Committee of Zhuhai People's Hospital Affiliated with Jinan University (Research projects IRB Review Approval Notice: LW-[2022]#1). Preparation of single-cell suspensions All tissue samples were washed twice with cold PBS. Tissue samples were cut into 1 mm 3 in size and placed in petri dish with cold PBS, then transferred into centrifuge tube, adding appropriate amount of enzyme and shaking at a certain temperature for a period of time. After 2–3 minutes' standing, supernatant were collected and then use a filter membrane to remove large clumps. After centrifuge the cells were collected, and then we resuspended the cells with red blood cell lysis buffer and incubate 2–3 min at room temperature and then centrifuge at 120×g under 4°C for 3min. Samples were resuspended again with cold PBS. Droplet-based single-cell sequencing Using the Single Cell 5’ Library and Gel Bead Kit (10X Genomics, 120237) and Chromium Single Cell A Chip Kit (10X Genomics, 120236), the cell suspension was loaded onto the Chromium single-cell controller (10X Genomics) to generate single-cell gel beads in the emulsion (GEMs) according to the manufacturer’s protocol. Briefly, single cells were suspended in PBS containing 0.04% bovine serum albumin. Approximately 10,000 cells were added to each channel, and about 6000 cells were recovered. The captured cells were lysed, and the released RNA was barcoded via reverse transcription in individual GEMs. Reverse transcription was performed at 53°C for 45 min, followed by 85°C for 5 min, and then the temperature was held at 4°C in a C1000 Touch Thermal Cycler (Bio Rad). After reverse transcription, single-cell droplets were broken and the single-strand cDNA was isolated and cleaned with Cleanup Mix containing DynaBeads (Thermo Fisher Scientific). cDNA was generated and amplified, and quality was assessed using the Agilent 4200. Single-cell RNA-seq libraries were prepared using Single Cell 5’ Library Gel Bead Kit V2 following the manufacture’s introduction. Next generation sequencing was performed on an Illumina Novaseq6000 with a sequencing depth of at least 100,000 reads per cell and pair end 150bp (performed by CapitalBio Technology, Beijing). Single cell RNA-seq (scRNA-seq) data processing Sequencing data were aligned to the human reference genome (GRCh38) and processed using the CellRanger (version 4.0.0). The gene expression matrix from the CellRanger pipeline was filtered, normalized using the Seurat R package (v3.2) [ 30 ]. Cells were selected if they met the following criteria: (i) top 99% of cells in unique molecular identifier counts; (ii) > 200 genes; and (iii) < 25% of mitochondrial gene expression in UMI counts. After removal of low-quality cells, the gene expression matrices were normalized to the total UMI counts per cell and transformed to the natural log scale. Then all the datasets of individual sample were integrated using the “FindIntegrationAnchors” and “IntegrateData” function in Seurat. Louvain algorithm was applied to iteratively group proximal cells together by “FindClusters” function with resolution of 0.6. Visualization was achieved by both the t-Distributed Stochastic Neighbor Embedding (tSNE) projection and Uniform Manifold Approximation and Projection (UMAP). Cell type annotations were performed on Blueprint and Encode reference dataset via SingleR [ 60 ], along with the marker-based correction. We classified all cells into eight major cell types, including T cells, B cells, NK cells, Monocytes, Epithelial cells, Fibroblasts, Endothelial cells, and Mast cells. 10x Visium Spatial transcriptomics (ST) Cryosections were cut at 10-µm thickness, mounted onto the GEX arrays. Sections were placed on Thermocycler Adaptor with the active surface facing up and incubated for 1 min at 37°C, and fixed for 30 min with methyl alcohol under − 20°C, and then stained with H&E (Eosin, Dako CS701, Hematoxylin Dako S3309, bluing buffer CS702). The brightfield Images were taken on a Leica DMI8 whole-slide scanner at 10X resolution. Visium spatial gene expression was processed using Visium spatial gene expression slide and Reagent Kit (10× Genomics, PN-1000184). For each well, Slide Cassette was used to create leakproof wells for adding reagents. 70 µL Permeabilization enzyme was added and incubated at 37°C for 20 min. Each well was washed with 100 µL SSC, and 75 µL reverse transcription Master Mix was added for cDNA Synthesis. cDNA library were prepared for sequencing. After the first-strand synthesis finished, reverse transcription Master Mix was removed from the wells, and then 75 µL 0.08 M KOH was added and incubated for 5 min at room temperature, then we removed the KOH from wells and washed with 100 µL EB buffer. 75 µL Second Strand Mix was added into each well for second-strand synthesis. cDNA amplification was performed on a S1000TM Touch Thermal Cycler (Bio Rad). According to the manufacture’s introduction, Visium spatial libraries were constructed using Visium spatial Library construction kit (10× Genomics, PN-1000184). The libraries were sequenced using an Illumina Novaseq6000 sequencer with a sequencing depth of at least 100,000 reads per spot with pair-end 150 bp (PE150) reading strategy (performed by CapitalBio Technology, Beijing). Spatial transcriptome sequencing (ST-seq) data processing The sequencing reads were mapped to the GRCh38 human genome and expression was quantified with the spaceranger-1.0.0. Further analysis was performed with Seurat (version 3.0.2). To annotate spots, we applied the integration workflow introduced in Seurat v3, which enabled the probabilistic transfer of cell types from the scRNA-seq data to the ST data. Specifically, we first identified pairwise correspondences between single cells and single spots to quantify the batch effect. Each spot was then annotated based on the transcriptomic similarity between spots and cell types in the scRNA-seq dataset. This probabilistic transfer procedure was implemented using the FindTransferAnchors (dims = 1:30) and TransferData (dims = 1:30) functions in Seurat with the combination of top 100 DEGs of each cell type. Differential expression and functional enrichment analysis After dimensional reduction and projection of all cells into two-dimensional space by tSNE and UMAP, cells were clustered together according to common features. The “FindAllMarkers” function in Seurat was used to find markers for each of the identified clusters. Using differentially expressed genes (DEGs) of each cluster, we performed functional enrichment analysis which were implemented by clusterprofiler (v3.10.1) with |log2Foldchange|>0 and p.adj < 0.05 as thresholds. The enrichment analysis of comprehensive functions including Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways, Reactome and Disease. Gene sets enrichment analysis was performed by GSEA application version of JAVA (v2.2.2.4), which used predefined gene sets from the Molecular Signatures Database (MSigDB, v6.2). Regulon analyses Regulon scores for individual cells were computed using the SCENIC (single-cell regulatory network inference and clustering) pipeline [ 61 ]. A log-normalized expression matrix of neuronal cells was used as an input into the pySCENIC workflow ( https://pyscenic.readthedocs.io/en/latest/ index.html) with default settings to infer regulons (master TFs and their target genes). CNV estimation and identification of malignant cells The chromosomal CNA profile of single cells was inferred by the R package inferCNV (version 1.0.4) [ 62 ]. Average signal were used as reference to define a baseline of normal karyotype such that their average copy number value was subtracted from all cells. The following parameters were applied: cutoff = 0.1, cluster_by_groups = TRUE, HMM = TRUE, and denoise = TRUE. Cell-cell communication analysis In order to explore cell-cell communications via ligand–receptor interactions, we employed the strategy proposed by Vento-Tormo et. al. [ 63 ] based on a public repository of ligands, receptors and interactions database CellPhoneDB (v2.0) [ 64 ]. The interaction score between two different cell types was mediated by a specific ligand-receptor pair based on the mean gene expression of ligand from one cell type and the corresponding receptor from another cell type. To identify the significant cell-cell interaction, we permuted the change of cell type label for each cell at 1,000 times to calculate the significance of each pair (p-value < 0.01). This procedure was performed between all pairs of cell types. The interactions between distinct cell subpopulations via putative ligand-receptor pairs were visualized using the ggplot2 package. Single-cell trajectory analysis We used Monocle v.2 [ 26 , 65 , 66 ] to illustrate the cell state transition in total epithelial cells, tumor cells in the CRC scRNA-seq dataset and in CRC5_1 tumor cryosection in the ST-seq dataset. This R package applied a reversed graph embedding technique to reconstruct single-cell trajectories. UMI count matrices and the negbinomial.size parameter were used to create a CellDataSet object in the default setting. We filtered variable genes with the following cutoff criteria: (1) genes expressed in more than 10 cells; (2) average expression value > 0.1; and (3) Qval < 0.01. These variable genes were used for semisupervised trajectory reconstruction. Dimensional reduction and cell ordering were performed using the DDRTree method and the orderCells function. Declarations Data availability All data used in this study have been deposited in the GSA (Genome Sequence Archive in BIG Data Center, Beijing Institute of Genomics, Chinese Academy of Sciences), under the accession number xxx (we have been uploading for three days, still not finished yet). The additional data that support the findings of this study are available from the corresponding author on request. Competing interests The authors declare no competing interests. Acknowledgments This work was supported by project grants from the National Natural Science Foundation of China (81901857, 32100447), Guangdong Provincial Key Laboratory of Tumor Interventional Diagnosis and Treatment (2021B1212040004), the Science and Technology Development Fund, Macau SAR (File No. 0085/2023/AMJ and 0002/2023/RIB1). References Binefa, G., et al., Colorectal cancer: from prevention to personalized medicine . World J Gastroenterol, 2014. 20(22): p. 6786–808. Subramanian, S., et al., Personalized medicine for prevention: can risk stratified screening decrease colorectal cancer mortality at an acceptable cost? Cancer Causes Control, 2017. 28(4): p. 299–308. Fanelli, G.N., et al., The heterogeneous clinical and pathological landscapes of metastatic Braf-mutated colorectal cancer . Cancer Cell Int, 2020. 20: p. 30. Ochocka, N., et al., Single-cell RNA sequencing reveals functional heterogeneity of glioma-associated brain macrophages . Nat Commun, 2021. 12(1): p. 1151. Zhang, Y., et al., Single-cell transcriptome analysis reveals tumor immune microenvironment heterogenicity and granulocytes enrichment in colorectal cancer liver metastases . Cancer Lett, 2020. 470: p. 84–94. Li, H., et al., Reference component analysis of single-cell transcriptomes elucidates cellular heterogeneity in human colorectal tumors . Nat Genet, 2017. 49(5): p. 708–718. Stahl, P.L., et al., Visualization and analysis of gene expression in tissue sections by spatial transcriptomics . Science, 2016. 353(6294): p. 78–82. Thrane, K., et al., Spatially Resolved Transcriptomics Enables Dissection of Genetic Heterogeneity in Stage III Cutaneous Malignant Melanoma . Cancer Res, 2018. 78(20): p. 5970–5979. Berglund, E., et al., Spatial maps of prostate cancer transcriptomes reveal an unexplored landscape of heterogeneity . Nat Commun, 2018. 9(1): p. 2419. Lundmark, A., et al., Gene expression profiling of periodontitis-affected gingival tissue by spatial transcriptomics . Sci Rep, 2018. 8(1): p. 9370. Asp, M., et al., Spatial detection of fetal marker genes expressed at low level in adult human heart tissue . Sci Rep, 2017. 7(1): p. 12941. Maniatis, S., et al., Spatiotemporal dynamics of molecular pathology in amyotrophic lateral sclerosis . Science, 2019. 364(6435): p. 89–93. Giacomello, S., et al., Spatially resolved transcriptome profiling in model plant species . Nat Plants, 2017. 3: p. 17061. Mantri, M., et al., Spatiotemporal single-cell RNA sequencing of developing chicken hearts identifies interplay between cellular differentiation and morphogenesis . Nat Commun, 2021. 12(1): p. 1771. Fawkner-Corbett, D., et al., Spatiotemporal analysis of human intestinal development at single-cell resolution . Cell, 2021. 184(3): p. 810–826 e23. Boehm, K.M., et al., Harnessing multimodal data integration to advance precision oncology . Nat Rev Cancer, 2021. Baccin, C., et al., Combined single-cell and spatial transcriptomics reveal the molecular, cellular and spatial bone marrow niche organization . Nat Cell Biol, 2020. 22(1): p. 38–48. Moncada, R., et al., Integrating microarray-based spatial transcriptomics and single-cell RNA-seq reveals tissue architecture in pancreatic ductal adenocarcinomas . Nat Biotechnol, 2020. 38(3): p. 333–342. Boyd, D.F., et al., Exuberant fibroblast activity compromises lung function via ADAMTS4 . Nature, 2020. 587(7834): p. 466–471. Maynard, K.R., et al., Transcriptome-scale spatial gene expression in the human dorsolateral prefrontal cortex . Nat Neurosci, 2021. 24(3): p. 425–436. Wu, T. and Y. Dai, Tumor microenvironment and therapeutic response . Cancer Lett, 2017. 387: p. 61–68. Binnewies, M., et al., Understanding the tumor immune microenvironment (TIME) for effective therapy . Nat Med, 2018. 24(5): p. 541–550. Lee, H.O., et al., Lineage-dependent gene expression programs influence the immune landscape of colorectal cancer . Nat Genet, 2020. 52(6): p. 594–603. Shi, X., et al., Genes involved in the transition from normal epithelium to intraepithelial neoplasia are associated with colorectal cancer patient survival . Biochem Biophys Res Commun, 2013. 435(2): p. 282–8. Lambrechts, D., et al., Phenotype molding of stromal cells in the lung tumor microenvironment . Nat Med, 2018. 24(8): p. 1277–1289. Trapnell, C., et al., The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells . Nat Biotechnol, 2014. 32(4): p. 381–386. Hong, X.C., et al., Epithelial membrane protein 3 functions as an oncogene and is regulated by microRNA-765 in primary breast carcinoma . Mol Med Rep, 2015. 12(5): p. 6445–50. Wang, Y.W., et al., EMP1, EMP 2, and EMP3 as novel therapeutic targets in human cancer . Biochim Biophys Acta Rev Cancer, 2017. 1868(1): p. 199–211. Wang, H., et al., Colorectal cancer risk variant rs7017386 modulates two oncogenic lncRNAs expression via ATF1-mediated long-range chromatin loop . Cancer Lett, 2021. 518: p. 140–151. Butler, A., et al., Integrating single-cell transcriptomic data across different conditions, technologies, and species . Nat Biotechnol, 2018. 36(5): p. 411–420. Stuart, T., et al., Comprehensive Integration of Single-Cell Data . Cell, 2019. 177(7): p. 1888–1902 e21. Yu, F., et al., IFITM1 promotes the metastasis of human colorectal cancer via CAV-1 . Cancer Lett, 2015. 368(1): p. 135–143. Andreu, P., et al., Identification of the IFITM family as a new molecular marker in human colorectal tumors . Cancer Res, 2006. 66(4): p. 1949–55. Russo, R.C., et al., The CXCL8/IL-8 chemokine family and its receptors in inflammatory diseases . Expert Rev Clin Immunol, 2014. 10(5): p. 593–619. Sack, U., et al., Novel effect of antihelminthic Niclosamide on S100A4-mediated metastatic progression in colon cancer . J Natl Cancer Inst, 2011. 103(13): p. 1018–36. Kwak, J.M., et al., Expression of protein S100A4 is a predictor of recurrence in colorectal cancer . World J Gastroenterol, 2010. 16(31): p. 3897–904. Calon, A., et al., Dependency of colorectal cancer on a TGF-beta-driven program in stromal cells for metastasis initiation . Cancer Cell, 2012. 22(5): p. 571–84. Crockford, D., et al., Thymosin beta4: structure, function, and biological properties supporting current and future clinical applications . Ann N Y Acad Sci, 2010. 1194: p. 179–89. Ricci-Vitiani, L., et al., Thymosin beta4 targeting impairs tumorigenic activity of colon cancer stem cells . FASEB J, 2010. 24(11): p. 4291–301. Wirsching, H.G., et al., Thymosin beta 4 gene silencing decreases stemness and invasiveness in glioblastoma . Brain, 2014. 137(Pt 2): p. 433–48. Makowiecka, A., et al., Changes in Biomechanical Properties of A375 Cells Due to the Silencing of TMSB4X Expression Are Not Directly Correlated with Alterations in Their Stemness Features . Cells, 2021. 10(4). Ripani, P., et al., Thiazolides promote G1 cell cycle arrest in colorectal cancer cells by targeting the mitochondrial respiratory chain . Oncogene, 2020. 39(11): p. 2345–2357. Kim, T.Y., et al., Oncogenic potential of a dominant negative mutant of interferon regulatory factor 3 . J Biol Chem, 2003. 278(17): p. 15272–8. Basit, A., et al., The cGAS/STING/TBK1/IRF3 innate immunity pathway maintains chromosomal stability through regulation of p21 levels . Exp Mol Med, 2020. 52(4): p. 643–657. Yang, R., et al., POU2F2 regulates glycolytic reprogramming and glioblastoma progression via PDPK1-dependent activation of PI3K/AKT/mTOR pathway . Cell Death Dis, 2021. 12(5): p. 433. Markiewski, M.M., et al., The Ribosomal Protein S19 Suppresses Antitumor Immune Responses via the Complement C5a Receptor 1 . J Immunol, 2017. 198(7): p. 2989–2999. Markiewski, M.M., et al., Modulation of the antitumor immune response by complement . Nat Immunol, 2008. 9(11): p. 1225–35. Bossi, F., et al., C1q as a unique player in angiogenesis with therapeutic implication in wound healing . Proc Natl Acad Sci U S A, 2014. 111(11): p. 4209–14. Bareke, H. and J. Akbuga, Complement system's role in cancer and its therapeutic potential in ovarian cancer . Scand J Immunol, 2018. 88(1): p. e12672. Kondoh, N., et al., Differential expression of S19 ribosomal protein, laminin-binding protein, and human lymphocyte antigen class I messenger RNAs associated with colon carcinoma progression and differentiation . Cancer Res, 1992. 52(4): p. 791–6. Xu, S., et al., The role of collagen in cancer: from bench to bedside . J Transl Med, 2019. 17(1): p. 309. Daviau, A., J.P. Couture, and R. Blouin, Loss of DLK expression in WI-38 human diploid fibroblasts induces a senescent-like proliferation arrest . Biochem Biophys Res Commun, 2011. 413(2): p. 282–7. Rodriguez Sawicki, L., et al., FABP1 knockdown in human enterocytes impairs proliferation and alters lipid metabolism . Biochim Biophys Acta Mol Cell Biol Lipids, 2017. 1862(12): p. 1587–1594. Loo, J.M., et al., Extracellular metabolic energetics can promote cancer progression . Cell, 2015. 160(3): p. 393–406. Tang, K., et al., Hypoxia-reprogrammed tricarboxylic acid cycle promotes the growth of human breast tumorigenic cells . Oncogene, 2019. 38(44): p. 6970–6984. Baghban, R., et al., Tumor microenvironment complexity and therapeutic implications at a glance . Cell Commun Signal, 2020. 18(1): p. 59. Bejarano, L., M.J.C. Jordao, and J.A. Joyce, Therapeutic Targeting of the Tumor Microenvironment. Cancer Discov, 2021. 11(4): p. 933–959. Cao, W., et al., Immunotherapy in endometrial cancer: rationale, practice and perspectives . Biomark Res, 2021. 9(1): p. 49. Li, H., et al., Underlying mechanisms and drug intervention strategies for the tumour microenvironment . J Exp Clin Cancer Res, 2021. 40(1): p. 97. Aran, D., et al., Reference-based analysis of lung single-cell sequencing reveals a transitional profibrotic macrophage . Nat Immunol, 2019. 20(2): p. 163–172. Aibar, S., et al., SCENIC: single-cell regulatory network inference and clustering . Nat Methods, 2017. 14(11): p. 1083–1086. Tirosh, I., et al., Dissecting the multicellular ecosystem of metastatic melanoma by single-cell RNA-seq. Science, 2016. 352(6282): p. 189–96. Vento-Tormo, R., et al., Single-cell reconstruction of the early maternal-fetal interface in humans . Nature, 2018. 563(7731): p. 347–353. Efremova, M., et al., CellPhoneDB: inferring cell-cell communication from combined expression of multi-subunit ligand-receptor complexes . Nat Protoc, 2020. 15(4): p. 1484–1506. Qiu, X., et al., Single-cell mRNA quantification and differential analysis with Census . Nat Methods, 2017. 14(3): p. 309–315. Qiu, X., et al., Reversed graph embedding resolves complex single-cell trajectories . Nat Methods, 2017. 14(10): p. 979–982. Supplementary Table 7 Supplementary Table 7 is not available with this version. Additional Declarations (Not answered) Supplementary Files FigureS1.pdf Figure S1. The annotations and cell proportions of eight cell type in CRC scRNA-seq. FigureS2.pdf Figure S2. The transcriptomic and functional features of malignant cells identified from epithelial cells. FigureS3.pdf Figure S3. The spatial distribution of typical genes on the cyrosections from paitents with CRC. FigureS4.pdf Figure S4. The trajectory and regulons of four spatial regions on the cyrosections from paitents with CRC. SupplementaryTableS16.xlsx Table S1. Gene and UMI of each patient and cryosections. Table S2. Cell number information of each major cluster. Table S3. All upregulated-gene of each major cluster. Table S4. Major clusters distribution in CRC and normal tissues. Table S5. Major clusters distribution in CRC and normal tissues of each patient. Table S6. The list of enriched GO terms using DEGs of malignant cells compared to non-malignant cells. Table S7. Patient samples and clinical characteristics. Cite Share Download PDF Status: Under Review Version 1 posted Editorial decision: revise 15 Jan, 2024 Review # 1 received at journal 11 Jan, 2024 Review # 2 received at journal 10 Jan, 2024 Reviewer # 2 agreed at journal 08 Jan, 2024 Reviewer # 1 agreed at journal 08 Jan, 2024 Reviewers invited by journal 08 Jan, 2024 Submission checks completed at journal 29 Dec, 2023 First submitted to journal 29 Dec, 2023 Editor assigned by journal 29 Dec, 2023 You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-3820516","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":265891993,"identity":"0b2ecf76-f379-4726-8d10-c369b4865a2d","order_by":0,"name":"Ligong Lu","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA/0lEQVRIiWNgGAWjYBACPgbmhgMgBht788EHH3jAggZ4tbAxMEK08PEcSzacQawWMENOIkdNGqKDkBaJxMaDXyru2LVJ5LBJ28jcSWxgb94mwVBzB5+WhsMyZ54lt/G8PWydw/MssYHnWJkEw7Fn+LVIth1OZmPPS7ydw3M4sUEix0yCseEwEVoYcgykLUBa5N8Q1nLwY9thOzaOHCNpBrAtPAS08DwEyp45nMAGCuQensPGbTxpxRYJx3Br4WdPPvzxR8Vhe/l2YFT+7Dks289+eOONDzW4tYAAMzA6EhtALMYeUEwBQQJeDUCFPxgY7CHMHwSUjoJRMApGwYgEAKaEWBGz5Ou6AAAAAElFTkSuQmCC","orcid":"https://orcid.org/0000-0003-1405-0052","institution":"Zhuhai Interventional Medical Center, Guangdong Provincial Key Laboratory of Tumour Interventional Diagnosis and Treatment, Zhuhai People's Hospital, Zhuhai Hospital affiliated with Jinan University","correspondingAuthor":true,"prefix":"","firstName":"Ligong","middleName":"","lastName":"Lu","suffix":""},{"id":265891994,"identity":"eb89509e-12bc-4ff8-9532-9d0214c58d89","order_by":1,"name":"Jing Xiao","email":"","orcid":"","institution":"","correspondingAuthor":false,"prefix":"","firstName":"Jing","middleName":"","lastName":"Xiao","suffix":""},{"id":265891995,"identity":"42903333-e34f-491b-a93e-0b9bfc833074","order_by":2,"name":"Xinyang Yu","email":"","orcid":"https://orcid.org/0000-0001-7431-3913","institution":"Zhuhai People’s Hospital","correspondingAuthor":false,"prefix":"","firstName":"Xinyang","middleName":"","lastName":"Yu","suffix":""},{"id":265891996,"identity":"5328f4ab-7955-4e29-8677-28231b9f85ad","order_by":3,"name":"Fanlin Meng","email":"","orcid":"","institution":"CapitalBio Technology Corporation","correspondingAuthor":false,"prefix":"","firstName":"Fanlin","middleName":"","lastName":"Meng","suffix":""},{"id":265891997,"identity":"c225405f-564c-4a2a-b72a-39a4ec38becd","order_by":4,"name":"Yuncong Zhang","email":"","orcid":"","institution":"Zhuhai People's Hospital","correspondingAuthor":false,"prefix":"","firstName":"Yuncong","middleName":"","lastName":"Zhang","suffix":""},{"id":265891998,"identity":"ba21847c-bea4-46f2-a9d2-d05180435e3b","order_by":5,"name":"Wenbin Zhou","email":"","orcid":"","institution":"Harbin Medical University","correspondingAuthor":false,"prefix":"","firstName":"Wenbin","middleName":"","lastName":"Zhou","suffix":""},{"id":265891999,"identity":"075a14ea-3c27-452f-a6f8-4a969d44e436","order_by":6,"name":"Yonghong Ren","email":"","orcid":"","institution":"Marketing and Management Department","correspondingAuthor":false,"prefix":"","firstName":"Yonghong","middleName":"","lastName":"Ren","suffix":""},{"id":265892000,"identity":"d9c80c03-65e8-4be7-b29d-64063939f78b","order_by":7,"name":"Jingxia Li","email":"","orcid":"","institution":"Zhuhai People's Hospital","correspondingAuthor":false,"prefix":"","firstName":"Jingxia","middleName":"","lastName":"Li","suffix":""},{"id":265892001,"identity":"40293b74-5149-46a5-9cf5-54839ae11653","order_by":8,"name":"Yimin Sun","email":"","orcid":"","institution":"CapitalBio Technology","correspondingAuthor":false,"prefix":"","firstName":"Yimin","middleName":"","lastName":"Sun","suffix":""},{"id":265892002,"identity":"34cda3eb-d178-4de4-9fe2-1e2f789c4836","order_by":9,"name":"Hongwei Sun","email":"","orcid":"","institution":"Zhuhai People’s Hospital","correspondingAuthor":false,"prefix":"","firstName":"Hongwei","middleName":"","lastName":"Sun","suffix":""},{"id":265892003,"identity":"fcf630b7-6a92-4367-825f-c40aa3fa4533","order_by":10,"name":"Guokai Chen","email":"","orcid":"https://orcid.org/0000-0003-0921-8244","institution":"University of Macau","correspondingAuthor":false,"prefix":"","firstName":"Guokai","middleName":"","lastName":"Chen","suffix":""},{"id":265892004,"identity":"1214e02d-9be6-4e85-8acb-2da329cb1ece","order_by":11,"name":"Ke He","email":"","orcid":"","institution":"Guangdong Second Provincial General Hospital","correspondingAuthor":false,"prefix":"","firstName":"Ke","middleName":"","lastName":"He","suffix":""}],"badges":[],"createdAt":"2023-12-29 09:20:35","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-3820516/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-3820516/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":49439966,"identity":"81f7f0d4-b9be-4686-aaf9-a6ec8f21b976","added_by":"auto","created_at":"2024-01-10 21:51:39","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":1034712,"visible":true,"origin":"","legend":"\u003cp\u003eCell type identification in human CRC by 10X Genomics scRNA-seq (a). Workflow of sample collection and single-cell transcriptome analysis from Chinese patients with CRC. (b). t-distributed stochastic neighbor embedding (t-SNE) plot of 35,666 cells from CRC patients (the CRC scRNA-seq dataset), grouped into eight major cell types (left panel, top). Proportions of the global cell types in tumor tissues, adjacent tissues and blood on average (left panel, bottom). The normalized expression of marker genes for each cell type (right panel) (c). Gene expression heatmap analyzed by 10X Genomics scRNA-seq. (d). Proportions of the global cell types in individual samples with CRC.\u003c/p\u003e","description":"","filename":"Figure1.png","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/00934aa99b012b186a3c486d.png"},{"id":49440840,"identity":"22477b0a-83d7-4a73-bcc2-dc8d9a90bfae","added_by":"auto","created_at":"2024-01-10 21:59:39","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":137513,"visible":true,"origin":"","legend":"\u003cp\u003eTranscriptome signatures and heterogeneity in normal and tumor epithelial cells. (a). t-SNE plot of the CRC scRNA-seq dataset color-coded by colorectal subtypes. (b). t-SNE plot of all 5,887 epithelial cells (tumor/malignant cells are included) of the CRC scRNA-seq dataset color-coded by subtypes. (c-e). The semisupervised trajectory of all epithelial cells inferred by Monocle v2, color-coded by state (c) or subtypes (d) or stemness. Stemness levels were calculated as the mean expression of stem-like signature (e) (f). Volcano plot showing differentially expressed genes between tumor cells and other normal epithelial cells (non-malignant cells). (\u003cem\u003eP\u003c/em\u003e-value \u0026lt; 0.05, Wilcoxon rank sum test, log\u003csub\u003ee\u003c/sub\u003e(fold change) \u0026gt; 0.25. (g). Significant biological processes (GO terms) enriched in tumor/malignant cells by clusterProfiler. (h). t-SNE plot of 3,150 tumor cells derived from the CRC scRNA-seq dataset, color-coded by cell subtypes (h). (i-j). The trajectory of tumor cells inferred by Monocle v2, color-coded by cell subtypes (i) and sample origins (j).\u003c/p\u003e","description":"","filename":"Figure2.png","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/cf915d99cd9eeaffc562b648.png"},{"id":49439968,"identity":"b71ea027-24c5-4a65-8bf1-ebc68835d03d","added_by":"auto","created_at":"2024-01-10 21:51:39","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":1292372,"visible":true,"origin":"","legend":"\u003cp\u003eSpatial transcriptome (ST) of CRC and mapping of cell types at spatial resolution. (a). A pathologic section from tumor tissues of one CRC patient (T0602). (b) Annotations obtained by integration analysis of the CRC scRNA-seq dataset and CRC5_1 in ST-seq dataset using seruat labeltransfer. (c). Clustering of the CRC5_1 ST spots and annotating CRC5_1 tumor cryosection on the ST slide. CRC5_1 cryosection was obtained from tumor tissues of patient T0602. (d). Expression levels for genes with subtype-specific patterns in CRC5_1 ST spots. (i). Standardized expression levels of five genes in the CRC5_1 in ST-seq datasets. (e). A pathologic section from normal, adjacent tissues of the CRC patient (T0602). (f) Annotations obtained by integration analysis of the CRC scRNA-seq dataset and CRC5N_1 in ST-seq dataset using seruat labeltransfer. (g). Clustering of the CRC5N_1 ST spots and annotating CRCN5_1 normal cryosection on the ST slide. CRC5N_1 cryosection was obtained from adjacent tissues of patient T0602. (h). Expression levels for genes with subtype-specific patterns in CRC5N_1 ST spots. (j). Standardized expression levels of five genes in the CRC5N_1 in ST-seq datasets.\u003c/p\u003e","description":"","filename":"Figure3.png","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/b4745d35635e5c5f1b2c9d59.png"},{"id":49439969,"identity":"b080b406-2b53-410d-b795-362f21759267","added_by":"auto","created_at":"2024-01-10 21:51:39","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":275239,"visible":true,"origin":"","legend":"\u003cp\u003eThe trajectory and interactions of cell lineages at spatial resolution. (a). The trajectory of all ST spots of CRC5_1 cryosection in the ST-seq dataset, color-coded by four annotated, spatial regions (tumor, stroma, immune infiltration, colon epithelium, as shown in Fig. 3c. (b). The distribution of ST spots of CRC5_1 cryosection on cell cycle phases. (c). Receptor-ligand pair expression in each pair of spatial regions in CRC5_1 cryosection using CellPhoneDB. (d). Hierarchical clustering of ST spots from the tumor region in CRC5_1 cryosection and indicative of five tumor subtypes based on five transcriptomic signatures. (e). The trajectory of five tumor subregions.\u003c/p\u003e","description":"","filename":"Figure4.png","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/a0b3b2743a0c962147c6d0c8.png"},{"id":49441054,"identity":"f3db4cb8-6d47-4c20-acbb-0de5ea322e81","added_by":"auto","created_at":"2024-01-10 22:07:40","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":2584156,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/83a1bd2c-be04-4242-9eec-a14e6caba229.pdf"},{"id":49439971,"identity":"893a562c-3457-41f5-83e9-27ee241b5788","added_by":"auto","created_at":"2024-01-10 21:51:40","extension":"pdf","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":6599475,"visible":true,"origin":"","legend":"\u003cp\u003eFigure S1. The annotations and cell proportions of eight cell type in CRC scRNA-seq.\u003c/p\u003e","description":"","filename":"FigureS1.pdf","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/bc4162d57cb0b39dc6a6b17f.pdf"},{"id":49439974,"identity":"35195e32-562f-4ba2-8166-1332e0e12964","added_by":"auto","created_at":"2024-01-10 21:51:40","extension":"pdf","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":4211313,"visible":true,"origin":"","legend":"\u003cp\u003eFigure S2. The transcriptomic and functional features of malignant cells identified from epithelial cells.\u003c/p\u003e","description":"","filename":"FigureS2.pdf","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/8db2abd8f2a8faf13877a1fc.pdf"},{"id":49439975,"identity":"ceaab1e6-5c8d-47e4-b907-68b6fb615b34","added_by":"auto","created_at":"2024-01-10 21:51:40","extension":"pdf","order_by":3,"title":"","display":"","copyAsset":false,"role":"supplement","size":7773714,"visible":true,"origin":"","legend":"\u003cp\u003eFigure S3. The spatial distribution of typical genes on the cyrosections from paitents with CRC.\u003c/p\u003e","description":"","filename":"FigureS3.pdf","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/dc49e2e9985d287103c7c8c5.pdf"},{"id":49439973,"identity":"c9da1ec3-31f5-4f21-b2f9-1aa0a1e75b6d","added_by":"auto","created_at":"2024-01-10 21:51:40","extension":"pdf","order_by":4,"title":"","display":"","copyAsset":false,"role":"supplement","size":1538830,"visible":true,"origin":"","legend":"\u003cp\u003eFigure S4. The trajectory and regulons of four spatial regions on the cyrosections from paitents with CRC.\u003c/p\u003e","description":"","filename":"FigureS4.pdf","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/822462673ec1f088bde19977.pdf"},{"id":49439970,"identity":"dd40909a-3a04-41a6-9b38-0751ee738880","added_by":"auto","created_at":"2024-01-10 21:51:40","extension":"xlsx","order_by":5,"title":"","display":"","copyAsset":false,"role":"supplement","size":755112,"visible":true,"origin":"","legend":"\u003cp\u003eTable S1. Gene and UMI of each patient and cryosections.\u003c/p\u003e\n\u003cp\u003eTable S2. Cell number information of each major cluster.\u003c/p\u003e\n\u003cp\u003eTable S3. All upregulated-gene of each major cluster.\u003c/p\u003e\n\u003cp\u003eTable S4. Major clusters distribution in CRC and normal tissues.\u003c/p\u003e\n\u003cp\u003eTable S5. Major clusters distribution in CRC and normal tissues of each patient.\u003c/p\u003e\n\u003cp\u003eTable S6. The list of enriched GO terms using DEGs of malignant cells compared to non-malignant cells.\u003c/p\u003e\n\u003cp\u003eTable S7. Patient samples and clinical characteristics.\u003c/p\u003e","description":"","filename":"SupplementaryTableS16.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-3820516/v1/5623eaf540d0d3f099d43744.xlsx"}],"financialInterests":"(Not answered)","formattedTitle":"Integrating spatial and single-cell transcriptomics reveals tumor heterogeneity and intercellular networks in colorectal cancer","fulltext":[{"header":"Introduction","content":"\u003cp\u003eColorectal cancer (CRC) is a molecularly heterogeneous disease [\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e]. The heterogeneity of cell types involved in CRC carcinogenesis makes it difficult to elucidate cell lineages using traditional developmental biology techniques such as bulk transcriptomics methods [\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e]. Through single-cell transcriptomics technology, it is now possible to deconstruct a tumor into its diverse cell subpopulations and thus gain a better understanding of the underlying biology like subtyping [\u003cspan additionalcitationids=\"CR5\" citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e]. However, spatial or anatomical information inherent in the tissue architecture is lost using single-cell transcriptomic technology only.\u003c/p\u003e \u003cp\u003eSpatial transcriptomics (ST) is an emerging technology that adds spatial dimensionality and tissue morphology information to the single-cell transcriptomics data of cells in an undissociated tissue, thus helping to preserve precise spatial or anatomical information. Overcoming the throughput limitation of in situ hybridization (ISH) methods, ST allows for unbiased mapping of transcripts in individual tissue sections with spatial resolution by using spatially barcoded oligo-deoxythymidine microarrays [\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e]. As a high-throughput spatially resolved transcriptomics tool, ST has been used to study architecturally complex tissues and diseases including melanoma [\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e], prostate cancer [\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e], gingival tissue [\u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e], adult human heart tissue [\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e], mouse and human spinal cord tissue [\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e], model plant species [\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e], as well as their spatiotemporal characterizations during development [\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e] [\u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eExtensive multimodal studies have unraveled molecular landscape of diverse diseases [\u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e]. Combining these two complementary and powerful technologies has been confirmed to be scalable to study architecturally complex tissues and to provide meaningful biological insight across a range of pathologies, such as bone marrow [\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e], pancreatic ductal adenocarcinomas [\u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e], lung fbroblasts [\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e], and human dorsolateral prefrontal cortex [\u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eThe tumor microenvironment (TME) comprises various cell types (immune cells, fibroblasts, endothelial cells, etc.) and extracellular components (growth factors, cytokines, extracellular matrix, hormones, etc.) that surround cancerous/tumor cells [\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e]. Since many currently used anticancer therapies target non-tumor components, such as the extracellular matrix (ECM), immune system and vascular system [\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e], understanding cellular components and how their dynamic interactions to shape the tumor landscape are particularly important.\u003c/p\u003e \u003cp\u003eIn this study, we aim to provide a comprehensive global view of tumor heterogeneity and intercellular interaction networks of CRC using single-cell transcriptional profiles coupled with spatial transcriptional profiles. By analyzing the single-cell and spatial transcriptional profiles of 41,700 cells from 3 treatment-na\u0026iuml;ve patients with CRC, we generated a molecular map of all major CRC populations based on single-cell RNA sequencing (scRNA-seq). The malignant cells in epithelial cells were identified and categorized into seven subclasses (tumor_CAV1, tumor_ATF3_JUN|FOS, tumor_ZEB2, tumor_VIM, tumor_WSB1, tumor_LXN, tumor_PGM1), which may help to the molecular subtyping of colorectal cancers. In addition, we used spatially resolved transcriptomics in combination with computational tools to attribute cell types to different CRC niches. Annotated tumor regions based on the cryosection sections represented high \u003cem\u003eTMSB4X\u003c/em\u003e expression, and suggested a typical marker of tumorgenesis. The stroma region was characterized by \u003cem\u003eVIM\u003c/em\u003e gene, which was also used as a typical feature of one subtype of malignant cells in CRC scRNA-seq.\u0026nbsp;Furthermore, we inferred the important interaction between tumor and stroma regions mediated by gene pair of C5AR1 and RPS19, which played roles of ligand and receptor, respectively.\u003c/p\u003e"},{"header":"Results","content":"\u003cp\u003e \u003cb\u003eLandscape view of cell composition in tumors, adjacent tissues and peripheral blood in patients with CRC\u003c/b\u003e \u003c/p\u003e \u003cp\u003eTo shed light on the complexity of the TME in CRC, we performed scRNA-seq along with spatial trancriptome sequencing on viable cells derived from matched tumor and adjacent tissues, as well as peripheral blood mononuclear cells (PBMCs) of 3 patients with CRC (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea, Supplementary Table \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e). On average, we obtained more than 150 G sequencing reads for each sample, with a median sequencing saturation of 91.40% (87.0%-95.5%). After quality control filters (few detected features in cells and few expressed cells associated with detected features), a total of 35,666 cells were identified in 9 samples derived from 3 patients (including 10347, 13241 and 18112 in tumor tissues, adjacent tissues and peripheral blood, respectively; Supplementary Table \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003e). We obtained approximately 1000 genes and 2500 unique molecular identifiers (UMIs) for each cell, indicating sufficient coverage and transcript representations.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo define cell clusters with similar expression profiles, we performed dimensionality reduction of t-distributed stochastic neighbor embedding (tSNE) implemented in the Seurat package. Each cluster was further identified as a specific cell subpopulation on the basis of the expression of the most variable genes and the canonical markers, including those in epithelial cells (with gene markers: EPCAM, KRT5, PHGR1, LGALS4, and TFF3), T cells (CD4\u0026thinsp;+\u0026thinsp;T cells: PTPRC, CD3D, and CD4; CD8\u0026thinsp;+\u0026thinsp;T cells: PTPRC, CD3D, and CD8A), B cells (CD19 and MS4A1), monocytes (CD14, ITGAX for CD11C), natural killer (NK) cells (FCGR3A and NCAM1), endothelial cells (CDH5, PLVAP, CLDN5, VWF), fibroblasts (LUM, DCN, COL1A1), and mast cells (KIT, CPA3, MS4A2, and TPSAB1) (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eb). In addition to these well-known markers, we also analyzed cluster-specific genes via differential gene expression analysis (Supplementary Table \u003cspan refid=\"MOESM3\" class=\"InternalRef\"\u003eS3\u003c/span\u003e). These cluster-specific marker genes included FBLN1 for fibroblasts, as well as MT1A and PLN for smooth muscle cells (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ec, \u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003eFig \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003ee\u003c/span\u003e). In total, eight cell types in CRC were identified based on canonical markers and cluster-specific genes: epithelial cells, fibroblasts, endothelial cells, monocytes, T cells, NK cells, B cells, and mast cells. The heterogeneous compositions of the TME in CRC across tumor tissues, normal tissues and peripheral blood are consistent with a recent single-cell transcriptome study of CRC [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e].\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo characterize different cell compositions in tumor tissues, normal tissues and peripheral blood in CRC, the proportions of each cell type were investigated. An overall increase in myeloid cell populations and decrease in B cell populations were observed in tumor tissues compared to normal tissues (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eb, bottom; Fig \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e d), suggesting a redirected immune response in CRC patients. In details, we observed that the proportion of monocytes was increased with approximately 2.5-fold, whereas that of NK cells and B cells was decreased (about 0.3\u0026ndash;0.4 times) in tumors compared to normal tissues, suggesting a myeloid immunosuppression in the CRC TME (Fig \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e a,b, Supplementary Table \u003cspan refid=\"MOESM4\" class=\"InternalRef\"\u003eS4\u003c/span\u003e). To further explore the distinct cell compositions in the TME across individuals, more detailed proportions were assessed (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ed). These results showed, for example, that in patient T0602, the proportion of epithelial cells decreased in contrast with patient T0529 and increased compared to that in patient T0609 (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ed, left; Fig \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e c1; Supplementary Table \u003cspan refid=\"MOESM5\" class=\"InternalRef\"\u003eS5\u003c/span\u003e). Since the transition from normal epithelium to intraepithelial neoplasia were found to be associated with CRC patient survival [\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e], the difference in epithelial cells across individuals may be important for survival and worthy of further investigation. Considering that the cellular proportion determined by scRNA-seq may be biased toward an underrepresentation of malignant cells derived from epithelial cells [\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e], we also explored the proportions of immune and stromal cells to account for all cells except epithelial cells which includes tumor cells. The results showed that myeloid cell-driven immune response in patient T0529 was stronger than that in the other two patients (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ed, right; \u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003eFig \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e c2\u003c/span\u003e).\u003c/p\u003e \u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003eEpithelial cells represents multilineages including a lineage of malignant cells\u003c/h2\u003e \u003cp\u003eIt has been suggested that human colon cancer cells recapitulate the multilineage differentiation processes of normal colon epithelia. To investigate each lineage contributing to the CRC heterogeneity at single cell resolution, we subclustered cell populations for each cell type to identify subpopulations. To annotate these subpopulations, we combined another published CRC cohort consisting of 6 CRC patients in tumor regions as well as matched normal mucosa [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e], and transferred the annotations of subtypes to our datasets in this study with the Seurat R package (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ea). Since the transition from normal epithelium to intraepithelial neoplasia were found to be associated with CRC patient survival [\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e], we focused on epithelial cells and found 9 subpopulations, namely CD19\u0026thinsp;+\u0026thinsp;CD20\u0026thinsp;+\u0026thinsp;B cells, crypt cells, enterocytes, goblet cells, intermediate, mature colonotypes, proinflammatory, stem-like, and tumor cells (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo distinguish malignant cells and nonmalignant cells in epithelial cells, we performed scRNA-seq-based copy number variation (CNV) and subclustering analysis (Fig \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003ea, b). The proportions of malignant cells in each subcluster of epithelial cells were shown in Fig \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003ec. The trajectory revealed a transcriptional hierarchy, defining seven molecular states (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ec, top). The cells from tumor tissues dominated the divergent differentiation states 2 and 5, suggesting the tissue arrangement along pseudotimes (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ec, bottom).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eTo illustrate the differentiation paths across the multilineages among the epithelial cell populations, the semisupervised trajectory inferred by monocle2 [\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e] revealed a transcriptional hierarchy defining three branches. The hierarchy was dominated by malignant epithelial cells, as well as normal epithelial cells (including goblet cells and (stem-like/transit amplifying cells) and immune-related cell types (including proinflammatory and mature colonotypes), which originated from normal epithelial cells with branching toward malignant epithelial cells (grey) and immune-related cell types (light green), respectively (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ed). Projection of malignant epithelial cells along the epithelial cell differentiation trajectory revealed segregation of tumor cells from normal epithelial cell types and stem-like populations. The greater stemness of malignant epithelial cells suggested the regenerative/proliferative potential of these tumor cells (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ee). The hypoxia and epithelial mesenchymal transition (EMT) were also investigated in the malignant epithelial cell populations (Fig \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003ee).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec4\" class=\"Section2\"\u003e \u003ch2\u003eTranscriptional and functional features of malignant cells reveal heterogeneity in CRC patients\u003c/h2\u003e \u003cp\u003eTo characterize the malignant cell populations, we scrutinized the transcriptional features between malignant and nonmalignant cells. The known malignant epithelial cell populations characterized by upregulated expression of S100A4, VEGFA, MYC, and ICAM1 (intercellular adhesion molecule-1), according to their significant differential expression (log\u003csub\u003ee\u003c/sub\u003e|fold change| \u0026gt; 0.25, \u003cem\u003ep\u003c/em\u003e value\u0026thinsp;\u0026lt;\u0026thinsp;0.05) (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ef, left). The most differentially expressed gene EMP3 (Epithelial membrane protein 3), which has been identified as an oncogene and therapeutic target in breast cancer [\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e, \u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e], remains to be elucidated in colon cancer. The characteristic genes in malignant epithelial cell populations were found to be involved in biological processes such as matrix remodelling, cell proliferation and apoptosis (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003ef, right, Supplementary Table S6), hinting the occurrence of EMT. Moreover, the terms \u0026lsquo;positive regulation of angiogenesis\u0026rsquo;, \u0026lsquo;cellular response to decreased oxygen levels\u0026rsquo; and \u0026lsquo;extracellular matrix organization\u0026rsquo; were enriched with differentially expressed genes, also suggestive of the malignant tendency of the cell populations. Using genes that characterizing malignant epithelial cell populations, the biological pathways that those genes implicated in were shown in Fig \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003ed. In addition, we predicted the regulons for malignant and non-malignant cells, respectively. The list of top 5 regulons for the two cell populations were shown in the heatmap (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003eg). In malignant cells, a famous oncogene \u003cem\u003eMYC\u003c/em\u003e was shown to be one of the key regulons. The transcriptional regulation role of \u003cem\u003eATF1\u003c/em\u003e in CRC cell lines has been characterized by a study of combing RNA-seq and ChIP-seq assays, in which found rs7017386 allele-specifcally enhanced the binding affnity of \u003cem\u003eATF1\u003c/em\u003e and promotion of two oncogenic lncRNAs via forming a long-range chromatin loop [\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e]. Therefore, both the characteristics of the Stemness, hypoxia, and EMT and the prediction of oncogenic regulons consistently reflected the malignancy and tumorgenesis roles of the identified malignant epithelial cells aformentioned.\u003c/p\u003e \u003cp\u003eTo focus on transcriptional programs for subcategorizing tumor cells, we performed subclustering and trajectory analysis for malignant cell populations. We assigned seven tumor cell subclusters to all the malignant cell populations, namely tumor_CAV1, tumor_ATF3_JUN|FOS, tumor_ZEB2, tumor_VIM, tumor_WSB1, tumor_LXN, and tumor_PGM1 (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e2\u003c/span\u003eh, Fig \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003ef). The trajectory revealed a transcriptional hierarchy, converging into three discrete tumor subclusters. One subcluster was highly enriched with gene response to histone deacetylase (HDAC) inhibitors (like \u003cem\u003eATF3\u003c/em\u003e and \u003cem\u003eCAV1\u003c/em\u003e), and another subcluster was highly enriched with inflammatory gene (\u003cem\u003eLXN\u003c/em\u003e and \u003cem\u003ePGM1\u003c/em\u003e). The subcluster enriched with tumor metastasis-related hypoxia (\u003cem\u003eWSB1\u003c/em\u003e) originated from the aforementioned differentiation paths. Sample arrangement along the differentiation trajectory of tumor cells revealed malignancy in patient T0602.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec5\" class=\"Section2\"\u003e \u003ch2\u003eTranscriptional features in a spatial resolution on cryosections from patients with CRC\u003c/h2\u003e \u003cp\u003eTo examine spatially transcriptional differences within colorectal tissues, we mounted cryosections of unfixed CRC tissues originating from the same CRC patients onto spatially barcoded ST microarray slides to generate unbiased transcriptome maps. After haematoxylin and eosin (H\u0026amp;E) staining and brightfield imaging, we annotated the slides according to the distinct histological features (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ea). The samples were then processed for ST analysis. We demultiplexed the sequenced reads and identified their spatial location within tissues using the ST location-specific barcodes of the array. For patient T0602, we detected approximately/average 3,413 median UMIs and approximately/average 1,660 median gene numbers per ST spot for both the tumor tissue section (named CRC5_1) and the normal tissue section (named CRCN5_1).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eFirst, the spatial transcriptomics data were integrated with the scRNA-seq data using Seurat-v3 anchor-based integration to annotate each region in the corresponding section [\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e, \u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e]. Every spot in the spatial data was considered a weighted mix of cell-types identified by scRNA-seq.\u0026nbsp;For each spot, the cell type with the maximum prediction score among all possible cell types and thus transferred from the scRNA-seq dataset is illustrated (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003eb, f). After further adjustment on the basis of annotated histological features, we annotated four and two anatomical regions in the CRC5_1 section (derived from a tumor tissue, Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ec), and in the CRCN5_1 section (derived from an adjacent tissue), separately (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003eg). We observed many obviously characteristic genes, which represented higher expression in annotated regions especially in the tumor tissues (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ei) compared to normal tissues (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ej). It is noted that five DEGs in comparision of malignant and non-malignant cells in the CRC scRNA-seq dataset were included IFITM1, CXCL1, CXCL8, S100A4, and TGFBI. The higher expression in tumor or stroma regions were shown in Fig \u003cspan refid=\"MOESM3\" class=\"InternalRef\"\u003eS3\u003c/span\u003e. \u003cem\u003eIFITM1\u003c/em\u003e was highly expressed and spatially restricted relative to the annotated tumor regions. \u003cem\u003eIFITM1\u003c/em\u003e is an interferon-induced transmembrane protein family member. Fang \u003cem\u003eet al.\u003c/em\u003e investigated that over-expression of \u003cem\u003eIFITM1\u003c/em\u003e promoted the aggressiveness of CRC cells, whereas knockdown of \u003cem\u003eIFITM1\u003c/em\u003e expression inhibited cell migration, invasion or tumorigenicity \u003cem\u003ein vitro\u003c/em\u003e [\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e]. Pauline \u003cem\u003eet al.\u003c/em\u003e found a highly significant increase in \u003cem\u003eIFITM\u003c/em\u003e mRNA levels in 154 patients with colon and rectal carcinomas, compared to corresponding normal tissues [\u003cspan citationid=\"CR33\" class=\"CitationRef\"\u003e33\u003c/span\u003e]. \u003cem\u003eCXCL1\u003c/em\u003e and \u003cem\u003eCXCL8\u003c/em\u003e, members of the angiogenic CXC chemokine family [\u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e] are highly expressed only on the bjojobsorder of annotated tumor regions, suggesting that an immune exclusion and blockage by the tumor. \u003cem\u003eS100A4\u003c/em\u003e (S100 calcium-binding protein A4) was a typical feature in the annotated stromal region and it has been confirmed to be associated with metastatic and recurrent CRCs [\u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e, \u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e]. \u003cem\u003eTGFB1\u003c/em\u003e is highly expressed and spatially restricted relative to the annotated stromal region. TGFβ-activated stroma was one of the key hallmarks of microsatellite-stable colorectal cancers displayed by a quadruple-mutant mouse model [\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e]. The TGF-β activity in stromal cells functionally increased the efficiency of organ colonization by CRC cells. In the oral squamous cell carcinoma, \u003cem\u003eTHBS1\u003c/em\u003e expression was induced by \u003cem\u003eTGFB1\u003c/em\u003e in the stroma, and the dependency on the TGF-β stromal programme for metastasis initiation may be exploited to improve the diagnosis and treatment of patients of CRC [\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e]. In addition, the genes characterizing classical phenotypes of cancer stem-like (\u003cem\u003eASCL2\u003c/em\u003e), hypoxia (\u003cem\u003eWSB1\u003c/em\u003e) and apoptosis (\u003cem\u003eATF3\u003c/em\u003e) were all highly expressed in the tumor or stromal region, suggesting cell differentiation programs. To facilitate comparisons, the expression of the aforementioned genes in the adjacent normal tissue samples is presented (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ej).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003cb\u003eDifferent anatomical regions in a spatial resolution represents cell compositions on tumor, stroma, immune and epithelium\u003c/b\u003e \u003c/p\u003e \u003cp\u003eSecond, standalone analysis of the spatial transcriptome along with the annotated anatomical regions was performed to identify spatially differential transcriptional programs. Interestingly, we found region-specific transcriptional differences. Unsupervised clustering of spatial transcriptome of CRC5_1 along with annotated anatomical regions revealed that the tumor region was characterized by heterogeneous cell cluster enriched high TMSB4X expressing cell populations (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ed, Fig \u003cspan refid=\"MOESM3\" class=\"InternalRef\"\u003eS3\u003c/span\u003e), which was not observed in CRCN5_1 (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ef). \u003cem\u003eTMSB4X\u003c/em\u003e encodes thymosin beta-4, a well-known secreted small peptide, which plays key roles in the cell motility, cytoskeletal organization, cell survival, and cell differentiation [\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e]. To the best of our knowledge, little is known about the spatiotemporal expression profile of \u003cem\u003eTMSB4X\u003c/em\u003e in CRC. We observed a tumor region-restriction feature within high TMSB4X expressing cell populations (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e3\u003c/span\u003ed). \u003cem\u003eTMSB4X\u003c/em\u003e has been validated as a therapeutic target in colon cancer stem cells in a previous study [\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e]. Additionally, in other architecture-dependent tissues whose spatial locations are deeply intertwined with their functions, such as brain and heart, TMSB4X has been reported to be involved in tumor progression via neovascularization, cell adhesion and the epithelial-mesenchymal transition [\u003cspan additionalcitationids=\"CR40\" citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e]. When TMSB4X was silenced both \u003cem\u003ein vitro\u003c/em\u003e and \u003cem\u003ein vivo\u003c/em\u003e, differentiation and tumorigenicity were diminished [\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e]. Targeting highly expressed TMSB4X or TMSB4X-high cell populations identified in spatial trascriptome in CRC has the potential to become a new therapeutic strategy in CRC. Moreover, the annotated stromal region overexpressed classical marker genes for fibroblasts or endothelial cells, such as \u003cem\u003eLUM\u003c/em\u003e, \u003cem\u003eVIM\u003c/em\u003e, \u003cem\u003eCOL1A1\u003c/em\u003e, and \u003cem\u003eCOL1A2\u003c/em\u003e. It is noted that the differential expression of \u003cem\u003eVIM\u003c/em\u003e gene in comparison of stroma regions and other regions in CRC ST-seq dataset was also identified in CRC scRNA-seq dataset,in which \u003cem\u003eVIM\u003c/em\u003e gene was one of the DEGs in malignant cells compared to non-malignant cells.\u003c/p\u003e \u003cp\u003eAs proximity is a necessity for physical interactions among cells, anatomical regions or cell-type proximity/interactive maps can be used to guide the discovery of interactions between anatomical regions or cell types in the same or different lineages. We \u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003efirst\u003c/span\u003e carried out a pseudotemporal trajectory analysis of the four anatomical regions in CRC5_1 section. Our data revealed that the lineage stemmed from the colon epithelium, and then went through an infiltration to divide into two major lineages, corresponding to stroma and tumor spots, respectively (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e4\u003c/span\u003ea, top). The trajectory reconstruction confirmed two major terminal/branching cell fates at tumor and stoma spots, respectively (Fig \u003cspan refid=\"MOESM4\" class=\"InternalRef\"\u003eS4\u003c/span\u003ea). Moreover, seven continuous states were identified during along pseudotime trajectory (Fig \u003cspan refid=\"MOESM4\" class=\"InternalRef\"\u003eS4\u003c/span\u003eb). Tumor spots were observed mainly in state 5, and stroma spots were found primarily in state 6, which indicated different states during tumorgenesis, even though these two anatomical regions were proximal to each other. The colon epithelium was observed mostly in state 1, and immune infiltration was distributed in states 2,3,4 and 7, sporadically in state 6.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eConsistent with the trajectory, the tumor and stroma spots exhibited greatest pseudotime meaning the most extent of differentiation and mature/terminal programme (Fig \u003cspan refid=\"MOESM4\" class=\"InternalRef\"\u003eS4\u003c/span\u003eb). To explore the characteristics of the branching cell fates, cell cycle phase analysis was performed. Approximately 70% of stroma spots arrested in the G1 phase of cell cycle, but this outcome was not observed in tumor anatomical regions (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e4\u003c/span\u003eb). G1 phase cell cycle arrest may be responsible for the inhibition of colorectal cell proliferation in the stroma region. Compounds that promote G1 cell cycle arrest were used and confirmed to be a treatment of colorectal cancer [\u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e].\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec6\" class=\"Section2\"\u003e \u003ch2\u003eSpatially resolved interactions of tumor and stroma regions\u003c/h2\u003e \u003cp\u003eTo investigate the regulators critical for the branching cell fates of tumor and stroma regions, SCENIC-based regulon analysis was performed. The results showed a spatial-resolved, specific regulon set (Fig \u003cspan refid=\"MOESM4\" class=\"InternalRef\"\u003eS4\u003c/span\u003ec). Four of the top 5 regulons in stroma or tumor region were overlapped, including \u003cem\u003eCDX1\u003c/em\u003e, \u003cem\u003eIRF8\u003c/em\u003e, \u003cem\u003eHNF4A\u003c/em\u003e,and \u003cem\u003eCREB3lL1\u003c/em\u003e (Fig \u003cspan refid=\"MOESM4\" class=\"InternalRef\"\u003eS4\u003c/span\u003ed). The stroma region was specifically regulated by \u003cem\u003eIRF3\u003c/em\u003e. The previous research suggested that overexpression of \u003cem\u003eIRF3\u003c/em\u003e causes cell-cycle arrest in the G1/S phase thereby resulting in inhibition of DNA synthesis [\u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e, \u003cspan citationid=\"CR44\" class=\"CitationRef\"\u003e44\u003c/span\u003e]. However, the tumor region was found to be specifically regulated by \u003cem\u003ePOU2F2\u003c/em\u003e. A recent research suggested that \u003cem\u003ePOU2F2\u003c/em\u003e played tumorigenic roles in glioblastoma by leading to a metabolic shift towards aerobic glycolysis [\u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e], but the roles in colorectal cancer remain poorly understood. In conclusion, the transcription factor \u003cem\u003eIRF3\u003c/em\u003e plays a key role in the generation of stromal cells in the stroma region, but not tumor cells in the tumor region, and the underlying programme involves in inhibition of cell proliferations by arresting cells in the G1/S phase.\u003c/p\u003e \u003cp\u003eTo further investigate the interaction between the neighboring anatomical region, especially between the tumor region and other regions, we \u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003esecondly\u003c/span\u003e performed CellphoneDB-based cell interaction analysis to investigate the underlying ligand-receptor pairs in different anatomical regions derived from CRC tumors. The aforementioned trajectory suggested the crosstalk between stroma region and tumor anatomical region (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e4\u003c/span\u003ea, top). We observed intensive cellular interactions between the stroma region and tumor regions. (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e4\u003c/span\u003ec). For example, the stroma region was predicted to interact with the tumor region in tumor tissues through C5AR1-RPS19, which is known to promote tumor growth by facilitating recruitment of these cells to tumors [\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e]. \u003cem\u003eC5AR1\u003c/em\u003e is known to activate and recruit myeloid-derived suppressor cells to tumors [\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e]. Accordingly, research studies have shown that the activated complement system has a tumor-promoting effect, including angiogenesis, trophoblastic invasion and tissue remodeling, which includes processes favorable for tumor establishment and progression [\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e, \u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e]. \u003cem\u003eRPS19\u003c/em\u003e, which is released from apoptotic tumor cells and upregulated in colon carcinoma compared to normal epithelial colon cells, has been reported by its novel immunosuppressive properties [\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eTo evaluate the refined characterization of the tumor region, we reclustered the spots in the tumor region to discern any spatial differentiation process. Reclustering the spots in the tumor region led to the identification of five subregions corresponding to five gene modules (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e4\u003c/span\u003ed): C0 spots expressed high levels of collagen, which is the major component of the TME and participates in cancer fibrosis [\u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e], including collagen type III (\u003cem\u003eCOL3A1\u003c/em\u003e), \u003cem\u003eCOL4A2\u003c/em\u003e, \u003cem\u003eCOL4A1\u003c/em\u003e, \u003cem\u003eCOL6A2\u003c/em\u003e, \u003cem\u003eCOL1A2\u003c/em\u003e, and \u003cem\u003eCOL1A1\u003c/em\u003e. In contrast to C1 spots, C0 spots in the outermost layer in the CRC5_1 section, constituted the primary structural element of the ECM. C1 spots expressed high levels of ECM transcripts implicated in cell migration, including fibronectin (FN1), tumor protein, transnationally-controlled 1 (TPT1), transforming growth factor beta induced (TGFBI), and whey acidic protein (WAP) 4-disulfide core domain protein 2 (WFDC2).In addition, C1 spots maintain an intermediate phase with both undifferentiated and differentiated phenotypes. Transferring cell type labels to spatial data suggested that only the epithelial cells in scRNA-seq were spatially restricted to the outermost layer of the section (tumor subcluster C0), and C1 was spatially restricted to the centre of the section. C2 spots expressed extremely high levels of proline rich and Gla domain 3 (PRRG3), which is a member of a family of vitamin K-dependent transmembrane proteins that contain a glutamate-rich extracellular domain. C2 spots characterized by a senescence-like proliferation, are represented by MAP3K12. MAP3K12 is an upstream activator of the mitogen-activated protein kinase (MAPK) pathway [\u003cspan citationid=\"CR52\" class=\"CitationRef\"\u003e52\u003c/span\u003e]. C3 spots expressed extremely high levels of fatty acid-binding protein 1 (FABP1), which is essential for proper lipid metabolism in differentiated enterocytes [\u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e53\u003c/span\u003e]. The creatine kinase B (CKB) promotes metastatic survival by modulating intra- and extracellular energetics [\u003cspan citationid=\"CR54\" class=\"CitationRef\"\u003e54\u003c/span\u003e]. C3 spots are characterized by metabolism. C4 spots expressed high levels of mitochondrial phosphoenolpyruvate carboxykinase 1 (PCK1) gene, which increases colon cancer cell growth in part by promoting the consumption of both glucose and glutamine in the tricarboxylic acid (TCA) cycle [\u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e]. The roles of hypoxia-reprogrammed TCA cycles in promoting human breast cancer cell growth via a HIF-1α-mediated PCK2 pathway have been reported [\u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e55\u003c/span\u003e]. C4 spots are characterized by hypoxia-response.\u003c/p\u003e \u003cp\u003eBased on characteristic genes of five subregions of tumor region, we defined five gene modules, including focal adhesion dynamics (C0), intermediate (C1), ECM (C2), metabolic (C3), hypoxia-response (C4) modules. The results suggested roles for these gene modules in tumor progression, implying a need for a combination of an anticancer therapy with corresponding modulators.\u003c/p\u003e \u003c/div\u003e"},{"header":"Discussion","content":"\u003cp\u003eColorectal cancer is a complex and heterogeneous malignant tumor of the colon and rectum. According to the degree of tumor differentiation and invasion, CRC can be classified into different subtypes. However, the published CRC studies usually were performed with whole tissues, such as bulk RNA-seq, blurring the heterogeneous characteristics of different cell types and limiting the ability to capture tumor heterogeneity. Cancerous tissue is composed of a mixture of various components, such as tumor cells, stromal cells, immune cells and ECM, leading to a complex TME [\u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e56\u003c/span\u003e]. The TME components exhibit interactive crosstalk with tumor cells and there surrounding factors, which in turn shapes tumor structure, metabolism and secretion, thus affects tumor development and/or metastasis. Immune cells within the TME play crucial roles during tumorgenesis. Immunotherapy aims to fight against cancer, infection and other diseases stimulating or suppressing the immune system. Immunotherapy displays promising therapeutic outcomes and limited side effects [\u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e57\u003c/span\u003e]. An increasing number of clinical trials have proven the effects of immunotherapy in certain types of solid tumors, such as melanoma, non-small-cell lung cancer, renal cancer and prostate cancer [\u003cspan citationid=\"CR58\" class=\"CitationRef\"\u003e58\u003c/span\u003e]. Recently, immunotherapy drugs such as CAR-T drugs have been approved by the FDA for clinical application. However, not all patients respond favorably to immunotherapy. Researchers have begun to examine the complexity and diversity of the TME and are realizing its importance in immunotherapy. Focusing on the TME facilitates to better understanding of the occurrence, development and metastasis of tumor and may lead to better diagnosis and treatment [\u003cspan citationid=\"CR59\" class=\"CitationRef\"\u003e59\u003c/span\u003e]. The expression pattern and function of tumor cell-associated immune molecules from in the TME provide useful information to determine whether a patient will might benefit from immunotherapy. There is an urgent need for improved techniques to better understand the TME and analyze the composition of immune cells and various other various cell types in tumor tissues. Single-cell sequencing and spatial transcriptomics have satisfied this requirement, and these two new-newly emerging techniques can be used to analyze samples at the single -cell level and monitor the \u003cem\u003ein situ\u003c/em\u003e spatial information of tumor tissues.\u003c/p\u003e \u003cp\u003eSingle-cell RNA sequencing enables to investigate into the transcriptional regulation of highly heterogeneous cell populations or subpopulations and facilitates the discovery of genes that indicate cell subtypes, or that mark intermediate states during a biological process, as well as bifurcation between two alternative cellular fates. Spatial RNA sequencing enables the anatomical, \u003cem\u003ein situ\u003c/em\u003e locations to be preserved, but cannot achieve rigorous single cell resolution. Actually, 10\u0026ndash;20 cells are typically identified in each spot in the ST assay. Since \u003cem\u003ein situ\u003c/em\u003e locations representing tissue sections are lacking in the scRNA-seq and single-cell resolution is lacking in the spatial RNA-seq, it is necessary to combine single-cell RNA-seq and spatial RNA-seq to reflect both cellular locations at true single-cell resolution. This allows inferences on functional relationships between scRNA-seq-defined populations based on their colocalization in space, and ultimately provides a more comprehensive characterization of cell types in their native environment than can be gained from either modality alone. In addition to transferring annotation at single cell level to spatial locations, we think that the current combinatory analysis between single-cell and spatial transcriptome remain to be more closely and complex.\u003c/p\u003e \u003cp\u003eIn this study, we combined single-cell and spatial transcriptomics to create a hierarchical map of cellular lineages in CRC. There CRC patients are involved, as well as single-cell transcriptomes of ~\u0026thinsp;47,000 cells, and spatial gene expression maps. We firstly constructed CRC single-cell maps consist of epithelial cells, mast cells, monocytes, T cells, B cells, Endothelial cells, and NK cells. Then we focused on the epithelial cells, from which we identified malignant cells. Also, we performed analyses of differential expression, functional enrichments, transcription factors to characterize the features of malignant cells compared to non-malignant cells. By sub-classification of malignant cells, we found seven subtypes according to transcriptional features, which could be helpful to molecular subtyping of CRC. Combining spatial and scRNA-seq datasets using bioinformatics approaches, we transferred cell type annotations at single cell level to those at spot level. The results showed four regions including tumor, stroma, immune-infiltration, and colon epithelium. Then we mainly focused on proximate cellular interactions (tumor and stroma) within colorectal tissue, to quantify anatomically restricted gene expression and explore crosstalks between tumor and stroma regions. Tumor and stroma regions were characterized by TMSB4X and VIM high expressions, respectively. The cellular interactions or crosstalks were inferred to be mediated by C5AR1 and RPS19, which remained to be further validations.\u003c/p\u003e \u003cp\u003eIn summary, we comprehensively explored the cellular landscape and reconstructed the putative interaction network consisting of tumor cells and their microenvironment. This collective view allowed us to elucidate how diverse cellular components jointly determine CRC molecular subtypes in individual patients.\u003c/p\u003e"},{"header":"Materials and Methods","content":"\u003cdiv id=\"Sec9\" class=\"Section2\"\u003e \u003ch2\u003eSubjects and clinical characteristics\u003c/h2\u003e \u003cp\u003ePatient samples were collected from 3 CRC patients (CRC0529, CRC0602, CRC0609), see detailed information in the Supplementary Table S7. Matched adjacent normal tissues and primary tumors and peripheral blood were obtained from all 3 patients. All patients were treatment-naive before tumour resection. All sampling and experimental steps in this study were approved by the Ethics Committee of Zhuhai People's Hospital Affiliated with Jinan University (Research projects IRB Review Approval Notice: LW-[2022]#1).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec10\" class=\"Section2\"\u003e \u003ch2\u003ePreparation of single-cell suspensions\u003c/h2\u003e \u003cp\u003eAll tissue samples were washed twice with cold PBS. Tissue samples were cut into 1 mm\u003csup\u003e3\u003c/sup\u003e in size and placed in petri dish with cold PBS, then transferred into centrifuge tube, adding appropriate amount of enzyme and shaking at a certain temperature for a period of time. After 2\u0026ndash;3 minutes' standing, supernatant were collected and then use a filter membrane to remove large clumps. After centrifuge the cells were collected, and then we resuspended the cells with red blood cell lysis buffer and incubate 2\u0026ndash;3 min at room temperature and then centrifuge at 120\u0026times;g under 4\u0026deg;C for 3min. Samples were resuspended again with cold PBS.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec11\" class=\"Section2\"\u003e \u003ch2\u003eDroplet-based single-cell sequencing\u003c/h2\u003e \u003cp\u003eUsing the Single Cell 5\u0026rsquo; Library and Gel Bead Kit (10X Genomics, 120237) and Chromium Single Cell A Chip Kit (10X Genomics, 120236), the cell suspension was loaded onto the Chromium single-cell controller (10X Genomics) to generate single-cell gel beads in the emulsion (GEMs) according to the manufacturer\u0026rsquo;s protocol. Briefly, single cells were suspended in PBS containing 0.04% bovine serum albumin. Approximately 10,000 cells were added to each channel, and about 6000 cells were recovered. The captured cells were lysed, and the released RNA was barcoded via reverse transcription in individual GEMs. Reverse transcription was performed at 53\u0026deg;C for 45 min, followed by 85\u0026deg;C for 5 min, and then the temperature was held at 4\u0026deg;C in a C1000 Touch Thermal Cycler (Bio Rad). After reverse transcription, single-cell droplets were broken and the single-strand cDNA was isolated and cleaned with Cleanup Mix containing DynaBeads (Thermo Fisher Scientific). cDNA was generated and amplified, and quality was assessed using the Agilent 4200. Single-cell RNA-seq libraries were prepared using Single Cell 5\u0026rsquo; Library Gel Bead Kit V2 following the manufacture\u0026rsquo;s introduction. Next generation sequencing was performed on an Illumina Novaseq6000 with a sequencing depth of at least 100,000 reads per cell and pair end 150bp (performed by CapitalBio Technology, Beijing).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003eSingle cell RNA-seq (scRNA-seq) data processing\u003c/h2\u003e \u003cp\u003eSequencing data were aligned to the human reference genome (GRCh38) and processed using the CellRanger (version 4.0.0). The gene expression matrix from the CellRanger pipeline was filtered, normalized using the Seurat R package (v3.2) [\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e]. Cells were selected if they met the following criteria: (i) top 99% of cells in unique molecular identifier counts; (ii)\u0026thinsp;\u0026gt;\u0026thinsp;200 genes; and (iii)\u0026thinsp;\u0026lt;\u0026thinsp;25% of mitochondrial gene expression in UMI counts. After removal of low-quality cells, the gene expression matrices were normalized to the total UMI counts per cell and transformed to the natural log scale. Then all the datasets of individual sample were integrated using the \u0026ldquo;FindIntegrationAnchors\u0026rdquo; and \u0026ldquo;IntegrateData\u0026rdquo; function in Seurat. Louvain algorithm was applied to iteratively group proximal cells together by \u0026ldquo;FindClusters\u0026rdquo; function with resolution of 0.6. Visualization was achieved by both the t-Distributed Stochastic Neighbor Embedding (tSNE) projection and Uniform Manifold Approximation and Projection (UMAP).\u003c/p\u003e \u003cp\u003eCell type annotations were performed on Blueprint and Encode reference dataset via SingleR [\u003cspan citationid=\"CR60\" class=\"CitationRef\"\u003e60\u003c/span\u003e], along with the marker-based correction. We classified all cells into eight major cell types, including T cells, B cells, NK cells, Monocytes, Epithelial cells, Fibroblasts, Endothelial cells, and Mast cells.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec13\" class=\"Section2\"\u003e \u003ch2\u003e10x Visium Spatial transcriptomics (ST)\u003c/h2\u003e \u003cp\u003eCryosections were cut at 10-\u0026micro;m thickness, mounted onto the GEX arrays. Sections were placed on Thermocycler Adaptor with the active surface facing up and incubated for 1 min at 37\u0026deg;C, and fixed for 30 min with methyl alcohol under \u0026minus;\u0026thinsp;20\u0026deg;C, and then stained with H\u0026amp;E (Eosin, Dako CS701, Hematoxylin Dako S3309, bluing buffer CS702). The brightfield Images were taken on a Leica DMI8 whole-slide scanner at 10X resolution.\u003c/p\u003e \u003cp\u003eVisium spatial gene expression was processed using Visium spatial gene expression slide and Reagent Kit (10\u0026times; Genomics, PN-1000184). For each well, Slide Cassette was used to create leakproof wells for adding reagents. 70 \u0026micro;L Permeabilization enzyme was added and incubated at 37\u0026deg;C for 20 min. Each well was washed with 100 \u0026micro;L SSC, and 75 \u0026micro;L reverse transcription Master Mix was added for cDNA Synthesis.\u003c/p\u003e \u003cp\u003ecDNA library were prepared for sequencing. After the first-strand synthesis finished, reverse transcription Master Mix was removed from the wells, and then 75 \u0026micro;L 0.08 M KOH was added and incubated for 5 min at room temperature, then we removed the KOH from wells and washed with 100 \u0026micro;L EB buffer. 75 \u0026micro;L Second Strand Mix was added into each well for second-strand synthesis. cDNA amplification was performed on a S1000TM Touch Thermal Cycler (Bio Rad). According to the manufacture\u0026rsquo;s introduction, Visium spatial libraries were constructed using Visium spatial Library construction kit (10\u0026times; Genomics, PN-1000184). The libraries were sequenced using an Illumina Novaseq6000 sequencer with a sequencing depth of at least 100,000 reads per spot with pair-end 150 bp (PE150) reading strategy (performed by CapitalBio Technology, Beijing).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003eSpatial transcriptome sequencing (ST-seq) data processing\u003c/h2\u003e \u003cp\u003eThe sequencing reads were mapped to the GRCh38 human genome and expression was quantified with the spaceranger-1.0.0. Further analysis was performed with Seurat (version 3.0.2). To annotate spots, we applied the integration workflow introduced in Seurat v3, which enabled the probabilistic transfer of cell types from the scRNA-seq data to the ST data. Specifically, we first identified pairwise correspondences between single cells and single spots to quantify the batch effect. Each spot was then annotated based on the transcriptomic similarity between spots and cell types in the scRNA-seq dataset. This probabilistic transfer procedure was implemented using the FindTransferAnchors (dims\u0026thinsp;=\u0026thinsp;1:30) and TransferData (dims\u0026thinsp;=\u0026thinsp;1:30) functions in Seurat with the combination of top 100 DEGs of each cell type.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec15\" class=\"Section2\"\u003e \u003ch2\u003eDifferential expression and functional enrichment analysis\u003c/h2\u003e \u003cp\u003eAfter dimensional reduction and projection of all cells into two-dimensional space by tSNE and UMAP, cells were clustered together according to common features. The \u0026ldquo;FindAllMarkers\u0026rdquo; function in Seurat was used to find markers for each of the identified clusters. Using differentially expressed genes (DEGs) of each cluster, we performed functional enrichment analysis which were implemented by clusterprofiler (v3.10.1) with |log2Foldchange|\u0026gt;0 and p.adj\u0026thinsp;\u0026lt;\u0026thinsp;0.05 as thresholds. The enrichment analysis of comprehensive functions including Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways, Reactome and Disease. Gene sets enrichment analysis was performed by GSEA application version of JAVA (v2.2.2.4), which used predefined gene sets from the Molecular Signatures Database (MSigDB, v6.2).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec16\" class=\"Section2\"\u003e \u003ch2\u003eRegulon analyses\u003c/h2\u003e \u003cp\u003eRegulon scores for individual cells were computed using the SCENIC (single-cell regulatory network inference and clustering) pipeline [\u003cspan citationid=\"CR61\" class=\"CitationRef\"\u003e61\u003c/span\u003e]. A log-normalized expression matrix of neuronal cells was used as an input into the pySCENIC workflow (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pyscenic.readthedocs.io/en/latest/\u003c/span\u003e\u003cspan address=\"https://pyscenic.readthedocs.io/en/latest/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e index.html) with default settings to infer regulons (master TFs and their target genes).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec17\" class=\"Section2\"\u003e \u003ch2\u003eCNV estimation and identification of malignant cells\u003c/h2\u003e \u003cp\u003eThe chromosomal CNA profile of single cells was inferred by the R package inferCNV (version 1.0.4) [\u003cspan citationid=\"CR62\" class=\"CitationRef\"\u003e62\u003c/span\u003e]. Average signal were used as reference to define a baseline of normal karyotype such that their average copy number value was subtracted from all cells. The following parameters were applied: cutoff\u0026thinsp;=\u0026thinsp;0.1, cluster_by_groups\u0026thinsp;=\u0026thinsp;TRUE, HMM\u0026thinsp;=\u0026thinsp;TRUE, and denoise\u0026thinsp;=\u0026thinsp;TRUE.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec18\" class=\"Section2\"\u003e \u003ch2\u003eCell-cell communication analysis\u003c/h2\u003e \u003cp\u003eIn order to explore cell-cell communications via ligand\u0026ndash;receptor interactions, we employed the strategy proposed by Vento-Tormo et. al. [\u003cspan citationid=\"CR63\" class=\"CitationRef\"\u003e63\u003c/span\u003e] based on a public repository of ligands, receptors and interactions database CellPhoneDB (v2.0) [\u003cspan citationid=\"CR64\" class=\"CitationRef\"\u003e64\u003c/span\u003e]. The interaction score between two different cell types was mediated by a specific ligand-receptor pair based on the mean gene expression of ligand from one cell type and the corresponding receptor from another cell type. To identify the significant cell-cell interaction, we permuted the change of cell type label for each cell at 1,000 times to calculate the significance of each pair (p-value\u0026thinsp;\u0026lt;\u0026thinsp;0.01). This procedure was performed between all pairs of cell types. The interactions between distinct cell subpopulations via putative ligand-receptor pairs were visualized using the ggplot2 package.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec19\" class=\"Section2\"\u003e \u003ch2\u003eSingle-cell trajectory analysis\u003c/h2\u003e \u003cp\u003eWe used Monocle v.2 [\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e, \u003cspan citationid=\"CR65\" class=\"CitationRef\"\u003e65\u003c/span\u003e, \u003cspan citationid=\"CR66\" class=\"CitationRef\"\u003e66\u003c/span\u003e] to illustrate the cell state transition in total epithelial cells, tumor cells in the CRC scRNA-seq dataset and in CRC5_1 tumor cryosection in the ST-seq dataset. This R package applied a reversed graph embedding technique to reconstruct single-cell trajectories. UMI count matrices and the negbinomial.size parameter were used to create a CellDataSet object in the default setting. We filtered variable genes with the following cutoff criteria: (1) genes expressed in more than 10 cells; (2) average expression value\u0026thinsp;\u0026gt;\u0026thinsp;0.1; and (3) Qval\u0026thinsp;\u0026lt;\u0026thinsp;0.01. These variable genes were used for semisupervised trajectory reconstruction. Dimensional reduction and cell ordering were performed using the DDRTree method and the orderCells function.\u003c/p\u003e \u003c/div\u003e "},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eData availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eAll data used in this study have been deposited in the GSA (Genome Sequence Archive in BIG Data Center, Beijing Institute of Genomics, Chinese Academy of Sciences), under the accession number xxx (we have been uploading for three days, still not finished yet). The additional data that support the findings of this study are available from the corresponding author on request.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCompeting interests\u0026nbsp;\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe authors declare no competing interests.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAcknowledgments\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThis work was supported by project grants from the National Natural Science Foundation of China (81901857, 32100447), Guangdong Provincial Key Laboratory of Tumor Interventional Diagnosis and Treatment (2021B1212040004), the Science and Technology Development Fund, Macau SAR (File No. 0085/2023/AMJ and\u0026nbsp;0002/2023/RIB1).\u003cbr\u003e\u0026nbsp;\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eBinefa, G., et al., \u003cem\u003eColorectal cancer: from prevention to personalized medicine\u003c/em\u003e. World J Gastroenterol, 2014. 20(22): p.\u0026nbsp;6786\u0026ndash;808.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSubramanian, S., et al., \u003cem\u003ePersonalized medicine for prevention: can risk stratified screening decrease colorectal cancer mortality at an acceptable cost?\u003c/em\u003e Cancer Causes Control, 2017. 28(4): p.\u0026nbsp;299\u0026ndash;308.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFanelli, G.N., et al., \u003cem\u003eThe heterogeneous clinical and pathological landscapes of metastatic Braf-mutated colorectal cancer\u003c/em\u003e. Cancer Cell Int, 2020. 20: p.\u0026nbsp;30.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOchocka, N., et al., \u003cem\u003eSingle-cell RNA sequencing reveals functional heterogeneity of glioma-associated brain macrophages\u003c/em\u003e. Nat Commun, 2021. 12(1): p.\u0026nbsp;1151.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang, Y., et al., \u003cem\u003eSingle-cell transcriptome analysis reveals tumor immune microenvironment heterogenicity and granulocytes enrichment in colorectal cancer liver metastases\u003c/em\u003e. Cancer Lett, 2020. 470: p.\u0026nbsp;84\u0026ndash;94.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi, H., et al., \u003cem\u003eReference component analysis of single-cell transcriptomes elucidates cellular heterogeneity in human colorectal tumors\u003c/em\u003e. Nat Genet, 2017. 49(5): p.\u0026nbsp;708\u0026ndash;718.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStahl, P.L., et al., \u003cem\u003eVisualization and analysis of gene expression in tissue sections by spatial transcriptomics\u003c/em\u003e. Science, 2016. 353(6294): p.\u0026nbsp;78\u0026ndash;82.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eThrane, K., et al., \u003cem\u003eSpatially Resolved Transcriptomics Enables Dissection of Genetic Heterogeneity in Stage III Cutaneous Malignant Melanoma\u003c/em\u003e. Cancer Res, 2018. 78(20): p.\u0026nbsp;5970\u0026ndash;5979.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBerglund, E., et al., \u003cem\u003eSpatial maps of prostate cancer transcriptomes reveal an unexplored landscape of heterogeneity\u003c/em\u003e. Nat Commun, 2018. 9(1): p.\u0026nbsp;2419.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLundmark, A., et al., \u003cem\u003eGene expression profiling of periodontitis-affected gingival tissue by spatial transcriptomics\u003c/em\u003e. Sci Rep, 2018. 8(1): p.\u0026nbsp;9370.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAsp, M., et al., \u003cem\u003eSpatial detection of fetal marker genes expressed at low level in adult human heart tissue\u003c/em\u003e. Sci Rep, 2017. 7(1): p.\u0026nbsp;12941.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eManiatis, S., et al., \u003cem\u003eSpatiotemporal dynamics of molecular pathology in amyotrophic lateral sclerosis\u003c/em\u003e. Science, 2019. 364(6435): p.\u0026nbsp;89\u0026ndash;93.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGiacomello, S., et al., \u003cem\u003eSpatially resolved transcriptome profiling in model plant species\u003c/em\u003e. Nat Plants, 2017. 3: p.\u0026nbsp;17061.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMantri, M., et al., \u003cem\u003eSpatiotemporal single-cell RNA sequencing of developing chicken hearts identifies interplay between cellular differentiation and morphogenesis\u003c/em\u003e. Nat Commun, 2021. 12(1): p.\u0026nbsp;1771.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFawkner-Corbett, D., et al., \u003cem\u003eSpatiotemporal analysis of human intestinal development at single-cell resolution\u003c/em\u003e. Cell, 2021. 184(3): p.\u0026nbsp;810\u0026ndash;826 e23.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBoehm, K.M., et al., \u003cem\u003eHarnessing multimodal data integration to advance precision oncology\u003c/em\u003e. Nat Rev Cancer, 2021.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBaccin, C., et al., \u003cem\u003eCombined single-cell and spatial transcriptomics reveal the molecular, cellular and spatial bone marrow niche organization\u003c/em\u003e. Nat Cell Biol, 2020. 22(1): p.\u0026nbsp;38\u0026ndash;48.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMoncada, R., et al., \u003cem\u003eIntegrating microarray-based spatial transcriptomics and single-cell RNA-seq reveals tissue architecture in pancreatic ductal adenocarcinomas\u003c/em\u003e. Nat Biotechnol, 2020. 38(3): p.\u0026nbsp;333\u0026ndash;342.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBoyd, D.F., et al., \u003cem\u003eExuberant fibroblast activity compromises lung function via ADAMTS4\u003c/em\u003e. Nature, 2020. 587(7834): p.\u0026nbsp;466\u0026ndash;471.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMaynard, K.R., et al., \u003cem\u003eTranscriptome-scale spatial gene expression in the human dorsolateral prefrontal cortex\u003c/em\u003e. Nat Neurosci, 2021. 24(3): p.\u0026nbsp;425\u0026ndash;436.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWu, T. and Y. Dai, \u003cem\u003eTumor microenvironment and therapeutic response\u003c/em\u003e. Cancer Lett, 2017. 387: p.\u0026nbsp;61\u0026ndash;68.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBinnewies, M., et al., \u003cem\u003eUnderstanding the tumor immune microenvironment (TIME) for effective therapy\u003c/em\u003e. Nat Med, 2018. 24(5): p.\u0026nbsp;541\u0026ndash;550.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLee, H.O., et al., \u003cem\u003eLineage-dependent gene expression programs influence the immune landscape of colorectal cancer\u003c/em\u003e. Nat Genet, 2020. 52(6): p.\u0026nbsp;594\u0026ndash;603.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShi, X., et al., \u003cem\u003eGenes involved in the transition from normal epithelium to intraepithelial neoplasia are associated with colorectal cancer patient survival\u003c/em\u003e. Biochem Biophys Res Commun, 2013. 435(2): p.\u0026nbsp;282\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLambrechts, D., et al., \u003cem\u003ePhenotype molding of stromal cells in the lung tumor microenvironment\u003c/em\u003e. Nat Med, 2018. 24(8): p.\u0026nbsp;1277\u0026ndash;1289.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTrapnell, C., et al., \u003cem\u003eThe dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells\u003c/em\u003e. Nat Biotechnol, 2014. 32(4): p.\u0026nbsp;381\u0026ndash;386.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHong, X.C., et al., \u003cem\u003eEpithelial membrane protein 3 functions as an oncogene and is regulated by microRNA-765 in primary breast carcinoma\u003c/em\u003e. Mol Med Rep, 2015. 12(5): p.\u0026nbsp;6445\u0026ndash;50.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang, Y.W., et al., \u003cem\u003eEMP1, EMP 2, and EMP3 as novel therapeutic targets in human cancer\u003c/em\u003e. Biochim Biophys Acta Rev Cancer, 2017. 1868(1): p.\u0026nbsp;199\u0026ndash;211.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang, H., et al., \u003cem\u003eColorectal cancer risk variant rs7017386 modulates two oncogenic lncRNAs expression via ATF1-mediated long-range chromatin loop\u003c/em\u003e. Cancer Lett, 2021. 518: p.\u0026nbsp;140\u0026ndash;151.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eButler, A., et al., \u003cem\u003eIntegrating single-cell transcriptomic data across different conditions, technologies, and species\u003c/em\u003e. Nat Biotechnol, 2018. 36(5): p.\u0026nbsp;411\u0026ndash;420.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStuart, T., et al., \u003cem\u003eComprehensive Integration of Single-Cell Data\u003c/em\u003e. Cell, 2019. 177(7): p.\u0026nbsp;1888\u0026ndash;1902 e21.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYu, F., et al., \u003cem\u003eIFITM1 promotes the metastasis of human colorectal cancer via CAV-1\u003c/em\u003e. Cancer Lett, 2015. 368(1): p.\u0026nbsp;135\u0026ndash;143.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAndreu, P., et al., \u003cem\u003eIdentification of the IFITM family as a new molecular marker in human colorectal tumors\u003c/em\u003e. Cancer Res, 2006. 66(4): p.\u0026nbsp;1949\u0026ndash;55.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRusso, R.C., et al., \u003cem\u003eThe CXCL8/IL-8 chemokine family and its receptors in inflammatory diseases\u003c/em\u003e. Expert Rev Clin Immunol, 2014. 10(5): p.\u0026nbsp;593\u0026ndash;619.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSack, U., et al., \u003cem\u003eNovel effect of antihelminthic Niclosamide on S100A4-mediated metastatic progression in colon cancer\u003c/em\u003e. J Natl Cancer Inst, 2011. 103(13): p.\u0026nbsp;1018\u0026ndash;36.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKwak, J.M., et al., \u003cem\u003eExpression of protein S100A4 is a predictor of recurrence in colorectal cancer\u003c/em\u003e. World J Gastroenterol, 2010. 16(31): p.\u0026nbsp;3897\u0026ndash;904.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCalon, A., et al., \u003cem\u003eDependency of colorectal cancer on a TGF-beta-driven program in stromal cells for metastasis initiation\u003c/em\u003e. Cancer Cell, 2012. 22(5): p.\u0026nbsp;571\u0026ndash;84.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCrockford, D., et al., \u003cem\u003eThymosin beta4: structure, function, and biological properties supporting current and future clinical applications\u003c/em\u003e. Ann N Y Acad Sci, 2010. 1194: p.\u0026nbsp;179\u0026ndash;89.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRicci-Vitiani, L., et al., \u003cem\u003eThymosin beta4 targeting impairs tumorigenic activity of colon cancer stem cells\u003c/em\u003e. FASEB J, 2010. 24(11): p.\u0026nbsp;4291\u0026ndash;301.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWirsching, H.G., et al., \u003cem\u003eThymosin beta 4 gene silencing decreases stemness and invasiveness in glioblastoma\u003c/em\u003e. Brain, 2014. 137(Pt 2): p.\u0026nbsp;433\u0026ndash;48.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMakowiecka, A., et al., \u003cem\u003eChanges in Biomechanical Properties of A375 Cells Due to the Silencing of TMSB4X Expression Are Not Directly Correlated with Alterations in Their Stemness Features\u003c/em\u003e. Cells, 2021. 10(4).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRipani, P., et al., \u003cem\u003eThiazolides promote G1 cell cycle arrest in colorectal cancer cells by targeting the mitochondrial respiratory chain\u003c/em\u003e. Oncogene, 2020. 39(11): p.\u0026nbsp;2345\u0026ndash;2357.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKim, T.Y., et al., \u003cem\u003eOncogenic potential of a dominant negative mutant of interferon regulatory factor 3\u003c/em\u003e. J Biol Chem, 2003. 278(17): p.\u0026nbsp;15272\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBasit, A., et al., \u003cem\u003eThe cGAS/STING/TBK1/IRF3 innate immunity pathway maintains chromosomal stability through regulation of p21 levels\u003c/em\u003e. Exp Mol Med, 2020. 52(4): p.\u0026nbsp;643\u0026ndash;657.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYang, R., et al., \u003cem\u003ePOU2F2 regulates glycolytic reprogramming and glioblastoma progression via PDPK1-dependent activation of PI3K/AKT/mTOR pathway\u003c/em\u003e. Cell Death Dis, 2021. 12(5): p.\u0026nbsp;433.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMarkiewski, M.M., et al., \u003cem\u003eThe Ribosomal Protein S19 Suppresses Antitumor Immune Responses via the Complement C5a Receptor 1\u003c/em\u003e. J Immunol, 2017. 198(7): p.\u0026nbsp;2989\u0026ndash;2999.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMarkiewski, M.M., et al., \u003cem\u003eModulation of the antitumor immune response by complement\u003c/em\u003e. Nat Immunol, 2008. 9(11): p.\u0026nbsp;1225\u0026ndash;35.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBossi, F., et al., \u003cem\u003eC1q as a unique player in angiogenesis with therapeutic implication in wound healing\u003c/em\u003e. Proc Natl Acad Sci U S A, 2014. 111(11): p.\u0026nbsp;4209\u0026ndash;14.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBareke, H. and J. Akbuga, \u003cem\u003eComplement system's role in cancer and its therapeutic potential in ovarian cancer\u003c/em\u003e. Scand J Immunol, 2018. 88(1): p.\u0026nbsp;e12672.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKondoh, N., et al., \u003cem\u003eDifferential expression of S19 ribosomal protein, laminin-binding protein, and human lymphocyte antigen class I messenger RNAs associated with colon carcinoma progression and differentiation\u003c/em\u003e. Cancer Res, 1992. 52(4): p.\u0026nbsp;791\u0026ndash;6.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXu, S., et al., \u003cem\u003eThe role of collagen in cancer: from bench to bedside\u003c/em\u003e. J Transl Med, 2019. 17(1): p.\u0026nbsp;309.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDaviau, A., J.P. Couture, and R. Blouin, \u003cem\u003eLoss of DLK expression in WI-38 human diploid fibroblasts induces a senescent-like proliferation arrest\u003c/em\u003e. Biochem Biophys Res Commun, 2011. 413(2): p.\u0026nbsp;282\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRodriguez Sawicki, L., et al., \u003cem\u003eFABP1 knockdown in human enterocytes impairs proliferation and alters lipid metabolism\u003c/em\u003e. Biochim Biophys Acta Mol Cell Biol Lipids, 2017. 1862(12): p.\u0026nbsp;1587\u0026ndash;1594.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLoo, J.M., et al., \u003cem\u003eExtracellular metabolic energetics can promote cancer progression\u003c/em\u003e. Cell, 2015. 160(3): p.\u0026nbsp;393\u0026ndash;406.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTang, K., et al., \u003cem\u003eHypoxia-reprogrammed tricarboxylic acid cycle promotes the growth of human breast tumorigenic cells\u003c/em\u003e. Oncogene, 2019. 38(44): p.\u0026nbsp;6970\u0026ndash;6984.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBaghban, R., et al., \u003cem\u003eTumor microenvironment complexity and therapeutic implications at a glance\u003c/em\u003e. Cell Commun Signal, 2020. 18(1): p.\u0026nbsp;59.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBejarano, L., M.J.C. Jordao, and J.A. Joyce, Therapeutic Targeting of the Tumor Microenvironment. Cancer Discov, 2021. 11(4): p.\u0026nbsp;933\u0026ndash;959.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCao, W., et al., \u003cem\u003eImmunotherapy in endometrial cancer: rationale, practice and perspectives\u003c/em\u003e. Biomark Res, 2021. 9(1): p.\u0026nbsp;49.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi, H., et al., \u003cem\u003eUnderlying mechanisms and drug intervention strategies for the tumour microenvironment\u003c/em\u003e. J Exp Clin Cancer Res, 2021. 40(1): p.\u0026nbsp;97.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAran, D., et al., \u003cem\u003eReference-based analysis of lung single-cell sequencing reveals a transitional profibrotic macrophage\u003c/em\u003e. Nat Immunol, 2019. 20(2): p.\u0026nbsp;163\u0026ndash;172.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAibar, S., et al., \u003cem\u003eSCENIC: single-cell regulatory network inference and clustering\u003c/em\u003e. Nat Methods, 2017. 14(11): p.\u0026nbsp;1083\u0026ndash;1086.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTirosh, I., et al., \u003cem\u003eDissecting the multicellular ecosystem of metastatic melanoma by single-cell RNA-seq.\u003c/em\u003e Science, 2016. 352(6282): p.\u0026nbsp;189\u0026ndash;96.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eVento-Tormo, R., et al., \u003cem\u003eSingle-cell reconstruction of the early maternal-fetal interface in humans\u003c/em\u003e. Nature, 2018. 563(7731): p.\u0026nbsp;347\u0026ndash;353.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eEfremova, M., et al., \u003cem\u003eCellPhoneDB: inferring cell-cell communication from combined expression of multi-subunit ligand-receptor complexes\u003c/em\u003e. Nat Protoc, 2020. 15(4): p.\u0026nbsp;1484\u0026ndash;1506.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eQiu, X., et al., \u003cem\u003eSingle-cell mRNA quantification and differential analysis with Census\u003c/em\u003e. Nat Methods, 2017. 14(3): p.\u0026nbsp;309\u0026ndash;315.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eQiu, X., et al., \u003cem\u003eReversed graph embedding resolves complex single-cell trajectories\u003c/em\u003e. Nat Methods, 2017. 14(10): p.\u0026nbsp;979\u0026ndash;982.\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"},{"header":"Supplementary Table 7","content":"\u003cp\u003eSupplementary Table 7 is not available with this version.\u003c/p\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":false,"highlight":"","institution":"","isAcceptedByJournal":true,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"cell-death-and-disease","isNatureJournal":false,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"cddis","sideBox":"Learn more about [Cell Death \u0026 Disease](http://www.nature.com/cddis/)","snPcode":"41419","submissionUrl":"https://mts-cddis.nature.com/cgi-bin/main.plex","title":"Cell Death \u0026 Disease","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"Nature AJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-3820516/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-3820516/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eSingle cell RNA sequencing (scRNA-seq), a powerful tool for studying the tumor microenvironment (TME), does not preserve/provide spatial information on tissue morphology and cellular interactions. To understand the crosstalk between diverse cellular components in proximity in the TME, we performed scRNA-seq coupled with spatial transcriptomic (ST) assay to profile 41,700 cells from three colorectal cancer（CRC） tumor-normal-blood pairs. Standalone scRNA-seq analyses revealed eight major cell populations, including B cells, T cells, Monocytes, NK cells, Epithelial cells, Fibroblasts, Mast cells, Endothelial cells. After the identification of malignant cells from epithelial cells, we observed seven subtypes of malignant cells that reflect heterogeneous status in tumor, including tumor_CAV1, tumor_ATF3_JUN|FOS, tumor_ZEB2, tumor_VIM, tumor_WSB1, tumor_LXN, and tumor_PGM1. By transferring the cellular annotations obtained by scRNA-seq to ST spots, we annotated four regions in a cryosection from CRC patients, including tumor, stroma, immune infiltration, and colon epithelium regions. Furthermore, we observed intensive intercellular interactions between stroma and tumor regions which were extremely proximal in the cryosection. Specially, one pair of ligands and receptors (C5AR1 and RPS19) was inferred to play key roles in the crosstalk of stroma and tumor regions. For the tumor region, a typical feature of \u003cem\u003eTMSB4X\u003c/em\u003e-high expression was identified, which could be a potential marker of CRC. The stroma region was found to be characterized by \u003cem\u003eVIM\u003c/em\u003e-high expression, suggesting it fostered a stromal niche in the TME.\u003c/p\u003e\n\u003cp\u003eCollectively, single cell and spatial analysis in our study reveal the tumor heterogeneity and molecular interactions in CRC TME, which provides insights into the mechanisms underlying CRC progression and may contribute to the development of anticancer therapies targeting on non-tumor components, such as the extracellular matrix (ECM) in CRC. The typical genes we identified may facilitate to new molecular subtypes of CRC.\u003c/p\u003e","manuscriptTitle":"Integrating spatial and single-cell transcriptomics reveals tumor heterogeneity and intercellular networks in colorectal cancer","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2024-01-10 21:51:35","doi":"10.21203/rs.3.rs-3820516/v1","editorialEvents":[{"type":"communityComments","content":0},{"type":"decision","content":"revise","date":"2024-01-15T13:01:54+00:00","index":"","fulltext":""},{"type":"editorInvitedReview","content":"This content is not available.","date":"2024-01-11T05:43:26+00:00","index":1,"fulltext":"This content is not available."},{"type":"editorInvitedReview","content":"This content is not available.","date":"2024-01-10T08:31:51+00:00","index":2,"fulltext":"This content is not available."},{"type":"reviewerAgreed","content":"This content is not available.","date":"2024-01-08T16:47:48+00:00","index":2,"fulltext":"This content is not available."},{"type":"reviewerAgreed","content":"This content is not available.","date":"2024-01-08T14:47:01+00:00","index":1,"fulltext":"This content is not available."},{"type":"reviewersInvited","content":"","date":"2024-01-08T13:43:00+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2023-12-29T11:13:40+00:00","index":"","fulltext":""},{"type":"submitted","content":"Cell Death \u0026 Disease","date":"2023-12-29T09:17:13+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2023-12-29T09:17:13+00:00","index":"","fulltext":""}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"cell-death-and-disease","isNatureJournal":false,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"cddis","sideBox":"Learn more about [Cell Death \u0026 Disease](http://www.nature.com/cddis/)","snPcode":"41419","submissionUrl":"https://mts-cddis.nature.com/cgi-bin/main.plex","title":"Cell Death \u0026 Disease","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"Nature AJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"6b98e039-1976-410e-a76f-a63a15a0ddb0","owner":[],"postedDate":"January 10th, 2024","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"under-review","subjectAreas":[{"id":28010577,"name":"Biological sciences/Cancer/Cancer microenvironment"},{"id":28010578,"name":"Biological sciences/Cell biology/Cell adhesion"}],"tags":[],"updatedAt":"2024-03-07T14:26:35+00:00","versionOfRecord":[],"versionCreatedAt":"2024-01-10 21:51:35","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-3820516","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-3820516","identity":"rs-3820516","version":["v1"]},"buildId":"qtupq5eGEP_6zYnWcrvyt","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}