Spatial Gene Expression Analysis Reveals Drivers of Extremely Early Lymph Node Metastasis in Breast Cancer

preprint OA: closed CC-BY-4.0
📄 Open PDF Full text JSON View at publisher

Abstract

Abstract Lymph node metastasis critically influences breast cancer prognosis; however, the cellular mechanisms underlying the earliest metastatic events remain unclear. In spatial transcriptomic analysis of a patient with breast cancer at single cell resolution, we identified 30 tumor cells representing the initial metastatic seeding in a lymph node. These cells originated from multiple epithelial–mesenchymal (EM) transition status and included six distinct subpopulations with biological significance. Only cells exhibiting metabolic shift toward fatty acid metabolism successfully established lymph node colonies, implicating this shift in metastatic fitness. The tumor microenvironment surrounding these cells showed immunosuppressive and tumor-promoting features, supporting metastasis establishment. Cross-referencing these expression profiles with public datasets revealed that poor prognosis correlated not with fully mesenchymal or metastatic populations, but with hybrid EM cells exhibiting epithelial and mesenchymal traits. These findings highlight the metabolic and phenotypic plasticity of metastatic cells and serve as translational bridges between the spatial evolution of tumor cells in the extremely early stages of lymph nodes metastasis and clinical prognosis in breast cancer.
Full text 144,922 characters · extracted from preprint-html · click to expand
Spatial Gene Expression Analysis Reveals Drivers of Extremely Early Lymph Node Metastasis in Breast 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 Spatial Gene Expression Analysis Reveals Drivers of Extremely Early Lymph Node Metastasis in Breast Cancer Satoi Nagasawa, Keiko Kajiya, Erina Ishikawa, Akinori Kanai, Ayako Suzuki, and 13 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-7371816/v1 This work is licensed under a CC BY 4.0 License Status: Published Journal Publication published 23 Jan, 2026 Read the published version in npj Breast Cancer → Version 1 posted 9 You are reading this latest preprint version Abstract Lymph node metastasis critically influences breast cancer prognosis; however, the cellular mechanisms underlying the earliest metastatic events remain unclear. In spatial transcriptomic analysis of a patient with breast cancer at single cell resolution, we identified 30 tumor cells representing the initial metastatic seeding in a lymph node. These cells originated from multiple epithelial–mesenchymal (EM) transition status and included six distinct subpopulations with biological significance. Only cells exhibiting metabolic shift toward fatty acid metabolism successfully established lymph node colonies, implicating this shift in metastatic fitness. The tumor microenvironment surrounding these cells showed immunosuppressive and tumor-promoting features, supporting metastasis establishment. Cross-referencing these expression profiles with public datasets revealed that poor prognosis correlated not with fully mesenchymal or metastatic populations, but with hybrid EM cells exhibiting epithelial and mesenchymal traits. These findings highlight the metabolic and phenotypic plasticity of metastatic cells and serve as translational bridges between the spatial evolution of tumor cells in the extremely early stages of lymph nodes metastasis and clinical prognosis in breast cancer. Biological sciences/Cancer Biological sciences/Computational biology and bioinformatics Health sciences/Oncology Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Introduction Axillary lymph nodes are considered the “first stop” in breast cancer metastasis 1-2 . The presence of metastatic tumor cells in these nodes is a well-recognized marker of poor prognosis 3-7 . Tumor cell profiles in metastases have been studied extensively. For example, Lee et al. used a mouse model and found that a shift toward fatty acid oxidation (FAO), driven by the selective activation of the transcription coactivator YAP in lymph node metastases, is essential for their establishment 8 . Liu et al. used single-cell RNA sequencing (RNA-seq) and spatial transcriptomics to compare the microenvironments of primary and lymph node metastatic tumors, reporting that in lymph node metastases, compared with primary tumors, T cell activation, cytotoxicity, and proliferation are markedly suppressed while dendritic cells show reduced ability to prime and activate T cells 9 . Although the requirements for lymph node metastasis and the niche environment are becoming clearer, how these molecular characteristics relate to cellular plasticity, possibly involving EMT and MET, remains uncertain. Notably, almost no studies have reported the first metastatic event in the lymph nodes. This is because it is difficult to track the trajectory of metastasis in humans, and it is extremely difficult to observe transient and plastic EMT states in vivo. Substantial evidence supports the role of epithelial-to-mesenchymal transition (EMT) and mesenchymal-to-epithelial transition (MET) in metastasis 10-11 . EMT is a cellular process in which epithelial cells acquire mesenchymal cell characteristics, causing dramatic changes in tissue structure and function. Cells take on a fibroblast-like morphology and increase their motility (ability to move) and invasiveness (ability to invade surrounding tissues). MET is a process that reverses the changes in cellular phenotype induced by EMT. In cancer cells, suppression of EMT-Transcription factors like TWIST1 is necessary to promote MET, which is required for proliferation after metastasis 12 .Despite extensive debate, it remains unclear how EMT or MET contributes to metastatic potential and whether it is a necessary condition for metastasis 10,13-16 . A major challenge in addressing this question lies in EMT’s broad and evolving definition. EMT was once considered a binary switch between epithelial and mesenchymal states. However, several in vitro studies have shown that EMT progresses gradually, forming hybrid subpopulations that occupy intermediate states between epithelial and mesenchymal identities 17-18 . Interestingly, tumor cell subpopulations nearest to the mesenchymal state do not necessarily exhibit the highest metastatic potential 19 . Indeed, high metastatic potential often lies in the hybrid state. Conversely, evidence suggests that metastatic potential and the reverse process, MET, are not correlated, indicating that EMT and MET may need to be considered separately 17 . To understand the molecular events driving metastasis, it is essential to examine how EMT progresses in the primary tumor by dissecting its intermediate states and how, and to what extent, metastasized cells restore epithelial traits in the lymph node 20 . In the present study, we used the Xenium in situ hybridization platform, which enables spatial gene expression profiling at the single-cell resolution, to track tumour cells in vivo during the metastasis process. The purpose of the study was to clarify the spatial evolution of tumor cells in the extremely early stages of breast cancer lymph nodes metastasis. Results Multi-omics analysis in each breast cancer region Before analyzing spatial patterns at the primary site and metastatic lymph node, we first examined the molecular features of the primary cancer via multiregional bulk analysis. The patient, an 80-year-old female (BRC-26), underwent surgery without neoadjuvant therapy. Her clinical diagnosis was stage II invasive breast carcinoma of the HER2 type, i.e., ER-negative, PgR-negative, and HER2-positive (for full clinicopathologic details, refer to Supplementary Table 1). For high-quality molecular analysis, especially RNA-seq, we freshly froze part of the tumor tissue harvested during surgery. We roughly dissected samples from histologically noncancerous regions (NC regions; NC-1 and NC-2) and cancerous regions (Ca regions; Ca-3 and Ca-4) (Fig. 1a and Supplementary Fig. 1a). Histologically, Ca-4 appeared sparser and less adhesive compared with Ca-3 (Supplementary Fig. 1b). We subjected these tissues to RNA-seq, enzymatic methylation sequencing (EM-seq), and whole-genome sequencing (WGS) (refer to Supplementary Table 2 for statistical analyses). We first performed transcriptome analysis using RNA-seq to compare gene expression between the NC and Ca regions. Of 4,584 differentially expressed genes (DEGs), 2,248 and 2,336 were significantly upregulated and downregulated in the cancer region (Supplementary Fig. 1c). Gene pathway enrichment analysis of these DEGs showed that “cell adhesion” pathways were enriched in NC regions, whereas pathways for “carboxylic acid metabolism,” and “cell cycle” were enriched in Ca regions (Fig. 1b), suggesting progressive malignant transformation from NC to Ca regions 21 . DNA methylation analysis revealed more extensive genome-wide hypomethylation in Ca regions compared with NC regions (Fig. 1c). Furthermore, promoters, intergenic regions, and gene bodies all showed reduced methylation (Supplementary Fig. 1d). Hypomethylation of intergenic and intron regions, especially repetitive sequences, can promote chromosomal instability and mutations, potentially increasing cancer risk 22-24 . These results confirm that breast cells undergo malignant transformation as they progress from NC to Ca regions. Within Ca regions, we observed significant changes in key metastasis-related genes (Fig. 1d). Expression of EMT markers, such as Vimentin 25 , and transcription factors Snail and ZEB2, considered, known EMT regulators 26-27 , increased from Ca-3 to Ca-4. Correspondingly, DNA methylation levels declined from Ca-3 to Ca-4 (Supplementary Fig. 1e), consistent with previous reports that EMT-related genes are epigenetically regulated 28-30 . Similarly, methylation of the THY1 gene, which was highly expressed in Ca-4, decreased from Ca-3 to Ca-4 (Fig. 1e, f). THY1, a GPI-anchored protein involved cell adhesion, migration, and polarity, can suppress lung metastasis when its integrin signaling is inhibited 31-32 . THY1 , which regulates cancer cell migration and invasion, also appears epigenetically controlled. These molecular profiles indicate that malignant transformation from NC to Ca regions, along with DNA methylation differences within Ca regions, may underlie varying EMT phenotypes in tumor cells 33-34 . WGS results further supported cancer progression from NC to Ca regions. We detected mutations in KMT2C (Y987H), the most frequently mutated histone methyltransferase in breast cancer 35 , across all regions (Fig. 1g). PIK3CA (H1047R), a known cancer driver, was mutated in all regions except NC-1. Mutations in transcription factor FOXA1 (R262P) and histone demethylase KDM6A (P281R) were specific to Ca regions (Fig. 1h). Copy number variation (CNV) analysis revealed chromosomal abnormalities in Ca regions (Fig.1i). These findings suggest that tumor development occurred against a background of clonal mammary gland expansion 36 . Overall, multiregional analyses indicated that breast cancer follows a well-characterized molecular trajectory from NC to Ca regions. Importantly, these molecular features aligned closely with histological observations, especially regarding gene regulation. Spatial analysis of the breast cancer at single-cell resolution To investigate the molecular mechanisms of metastasis in detail, we performed spatial gene expression analysis at single-cell resolution on primary tumors and paired axillary metastatic lymph nodes from the patient shown in Fig. 1. In addition to the primary tumor adjacent to the bulk analysis region and lymph node metastases from the same case, we included a tumor-draining lymph node (TdL) as a comparison specimen (Fig. 2a). The TdL was a clinically diagnosed as metastasis-negative by rapid intraoperative evaluation. We applied Xenium spatial gene expression profiling using a custom 380-gene breast cancer panel (refer to Supplementary Table 3 for the gene list). Across the primary tumor, metastatic lymph node and TdL, Xenium detected 1,005,436 cells (Fig. 2b). Using gene expression signatures 37-39 and histological location information, we annotated clusters by major cell types. The luminal cell cluster was subdivided into two clusters, (i) and (ii), with (i) further branching into four subclusters: L1: luminal 1; L2: luminal 2; B: basal; and M: mesenchymal like (Fig. 2c). We identified genes significantly upregulated within each cluster. KRT8 and KRT18 , markers of mammary ductal epithelial cells, were enriched in L1 and L2, whereas KRT5 and KRT14 , myoepithelial markers, were enriched in cluster B (Supplementary Fig. 2). Consistent with expression information, spatial localization of L1, L2, and B clusters matched the lobular architecture of the mammary gland (Fig. 2d). High ANKRD30A expression in L2 and elevated KIT and KRT23 expression in L1 aligned with two previously reported mature luminal cell types 40 : mature luminal cell 2 ( KRT18 +/ ANKRD30A +) and mature luminal cell 1 ( KRT18 +/ KRT23 +/ KIT +). L1, L2, and B clustered on the lobules, whereas M cells were distributed sporadically throughout the tumor, suggesting reduced cell adhesion (Fig. 2e). Consistent with this, the M cluster showed increased expression of EMT markers, including THY1 , an EMT-related gene shown in Fig. 1e. Collectively, these findings support a role for THY1 expression in EMT initiation from normal lobule cells (L1, L2, and B), potentially regulated through DNA methylation. Although cluster (i) encompassed lobular and early EMT-stage cells, cluster (ii) localized pathologically to tumor regions and comprised definitive cancer cells. Strong ANKRD30A expression in cluster (ii) (Fig. 2f) suggested that L2 may be the tumor’s cell of origin 40 . Uniform manifold approximation and projection (UMAP) showed that metastatic lymph node tumor cells formed distinct clusters. Based on high ANKRD30A expression, we defined L2 as the tumor-originating cell and applied pseudotemporal trajectory analysis using Monocle 3 41 . The inferred trajectory showed branching from L2 into (ii)-1, progressing to (ii)-2, and finally differentiating into metastatic tumor cells in lymph nodes (Fig. 2g). Taken together, these findings suggest a progression pathway in which cancer cells originate from L2, acquire EMT potential in the M cluster, and ultimately establish metastases within cluster (ii) in lymph nodes. Spatial analysis of the lymph nodes reveals six transcriptionally distinct EMT states When inspecting spatial gene expression data to further investigate the mechanisms underlying metastasis, we unexpectedly identified a small cluster of only 30 cells within the TdL, a site clinically diagnosed as metastasis-free; these cells spanned a region approximately 200 μm in diameter (Fig. 3a). Gene expression analysis revealed clear KRT19 expression, an epithelial marker absent in normal lymph node cells. Additional hematoxylin and eosin (HE) staining confirmed their morphology was consistent with cancer cells. The identification of these cells in the TdL provided a rare opportunity to examine molecular events at an extremely early stage of metastasis. This form of early dissemination is clinically referred to as isolated tumor cells (ITCs), a term we adopt hereafter. When we visualized ITCs on the UMAP plane, we observed transcriptional heterogeneity (Fig. 3b, left panel), indicating that ITCs comprise several distinct cell types. To better characterize the biology of ITCs and the roles of different cell states in metastasis, we applied nonhierarchical clustering to subclassify tumor cells into six groups (C1–C6) and assigned each ITC to one of these clusters (Fig. 3b, right panel). Owing to the small number of ITCs detected, subsequent analyses focused on the clusters to which each ITC mapped. For each subpopulation associated with an ITC, we identified distinct signaling pathways linking histology, spatial location, and gene expression (Fig. 3c–e and Supplementary Table 4). C6 occupied a position on the UMAP closest to EMT initiation (corresponding to the M cluster in Fig. 2c). Its scattered distribution within the primary tumor suggested low intercellular adhesion. C6 cells expressed high levels of ETM-related markers, such as VIM and CDH3 , as well as stemness markers characteristic of breast cancer, including CD24 low and CD44 high (Fig. 3d). In contrast, epithelial marker expression was low, indicating a sustained mesenchymal phenotype. Gene enrichment analysis showed activation of EMT, stem cell signaling, TGF-β signaling, WNT signaling, immune checkpoint pathways, and antianoikis signaling (Fig. 3e). Signals from the microenvironment, such as TGF-β and Wnt ligands, induce various EMT-related transcription factors 42-43 . EMT also promotes PD-L1 expression, activating immune checkpoint pathways 44-45 . These findings suggest that C6 represents tumor cells in an extreme mesenchymal phenotype, capable of evading anoikis and immune clearance. In contrast, C4 most strongly retained epithelial features and was the most predominant cell population within lymph node metastases, forming the tumor mass (Fig. 3c). C4 cells showed enrichment of cholesterol homeostasis and fatty acid metabolism pathways, with upregulation of fatty acid metabolism–related molecules in metastatic lymph nodes (Fig. 3e, f). A previous mouse-based lymph node metastasis study demonstrated that a metabolic shift toward FAO in tumor cells is required for lymph node colonization by cancer cells 8 . Thus, C4 represents a cancer cell population that has acquired an extreme epithelial phenotype and metabolic reprogramming suited for metastatic growth. Clusters C1, C2, C3, and C5 represented hybrid epithelial–mesenchymal (EM) cell populations that coexpressed epithelial and mesenchymal markers (Fig. 3d). The degree of coexpression differs among clusters, each cluster exhibited distinct features. C5 showed the strongest enrichment for the G2M checkpoint pathway (Fig. 3g). Prior studies have indicated that breast cancers with G2M activity display higher proliferative activity, increased MYC pathway activation, earlier metastasis, and worse survival 46 , suggesting that C5 may represent a more aggressive hybrid EMT subpopulation. C3 was enriched for VEGF signaling and glycolysis, likely reflecting a hypoxic environment in the tumor core that promotes angiogenesis and glycolysis 47 . C2, localized at the tumor periphery, i.e., the invasive front, displayed activation of matrix metalloproteinases (MMPs) and canonical Wnt signaling. Cells at tumor margins exhibit EM plasticity and migratory behavior 48 , with MMP signaling 42 and the Wnt/β catenin pathway 49 known to promote invasion. Canonical Wnt signaling was also enriched in C2 and hybrid EMT subpopulations (C3 and C6), consistent with reports that Wnt7A/B maintains the hybrid EM state 12,48 . C1 showed activation of retinoid metabolism and lipid transport pathways. Retinoic acid activates retinoic acid receptors and retinoid X receptors, driving expression of fatty acid metabolism genes 50 . Considering that C1 gene expression aligns closely with C4 (the most significantly proliferating cluster in lymph node metastases, involving elevated fatty acid metabolism), we speculate that C1 may act as a precursor population to C4. Collectively, these findings suggest that hybrid EM cells undergo functional adaptation within primary tumor, with C6 as a possible EMT origin, C2 as an intermediate, C3 and C5 representing divergent aggressive states, and C1 transitioning toward the epithelialized, proliferative C4. The important thing is that each hybrid EMT subtype has different functional capabilities and that these cells are present in primary tumour from the very early stages of metastasis (Fig. 3h). Analysis of the tumor microenvironment in metastasized lymph node The major bottleneck for cancer cells infiltrating a lymph node is the new tumor microenvironment (TME) 51 . To address this, we compared interactions between tumor cells and their surrounding TME in the TdL and the metastatic lymph node (Fig. 4a). Notably, we detected distinct interaction patterns specific to each TME (Fig. 4b). One interaction detected uniquely in the TdL was the CD45–MRC1 ligand–receptor interaction, occurring between CD45 on T cells and MRC1 on macrophages (Fig. 4c). MRC1, an endocytosis receptor belonging to the C-type lectin family, is expressed on cells like dendritic cells, macrophages, and endothelial cells 52 . It has been shown to inhibit CD45 activity on T cells via direct interaction, leading to upregulation of the immune checkpoint protein CTLA4 and induction of antigen-specific T cell tolerance 53-55 . Therefore, immune escape may begin through this mechanism in the TdL. Metastatic lymph nodes uniquely exhibited coinhibitory ligand–receptor interactions, such as VEGFA–KDR and CD86–CTLA4, consistent with the promotion of angiogenesis and the establishment of a T cell–exhausted environment already formed at this stage (Fig. 4d). Among interactions involving collagen, cancer-associated fibroblasts (CAFs) exhibited the highest number of interactions (Supplementary Fig. 3). Representative interactions in metastatic lymph nodes included COL1A1–CD44 and COL4A1–CD44, both known to promote tumor proliferation and progression 56-57 . In the chemokine signaling category, the CXCL16–CXCR6 interaction was particularly prominent in metastatic lymph nodes (Supplementary Fig. 3). High CXCL16 expression is associated with histological malignancy and is implicated in tumor progression and metastasis via activation of the CXCL16–CXCR6 axis 58 . These TME-derived signals appear to constitute necessary conditions for successful metastasis. Transcriptional profiling of hybrid EMT and its association with patient outcomes Prompted by the preceding analyses, we sought to determine whether the transcriptional profiles we identified are associated with patient prognosis in a broader clinical context. To this end, we reanalyzed a previously reported scRNA-seq dataset from metastatic lymph nodes of patients with breast cancer 4 (Fig. 5a). Further clustering of cancer cells from the dataset revealed five clusters (sc8,sc4,sc2,sc15 and sc10) (Fig. 5b). According to publicly available information, sc2,sc4 and sc8 were metastatic lymph node-derived samples, and sc15, sc10 was derived from the primary tumor. A metabolic shift toward fatty acids were observed in metastatic lymph node-derived samples(sc2,sc4 and sc8)(Fig. 5c). These observations are consistent with our results shown in Figure 3 and support our findings that a metabolic shift toward fatty acid metabolism is necessary for the formation of metastatic colonies in lymph nodes. To investigate whether molecular subtype signatures derived from our Xenium data are associated with clinical outcomes, we evaluated their prognostic relevance using the METABRIC dataset 59 . Interestingly, the subtypes associated with poor prognosis were not clusters that formed metastatic colonies in lymph node (e.g., C4) or those showing strong mesenchymal features (C6) but rather clusters C5 and C3, which exhibited signatures of G2M cell cycle progression and increased glycolytic activity (Supplementary Fig. 4). To validate this finding, we mapped the hybrid EM clusters identified in our Xenium analysis (C1–C6) to the corresponding cancer cell subclusters (sc2, sc4, sc8, sc10, and sc15) in the public scRNA-seq dataset. This revealed a one-to-one correspondence between C5 and sc8 and between C6 and sc15 (Fig. 5d). Consistent with our findings, only sc8 (corresponding to the C5 signature in our Xenium data) was significantly associated with poor prognosis [log-rank P = 4e-09, hazard ratio (HR) = 1.536, 95% confidence interval (CI) = 1.330–1.772] (Fig. 5e). In contrast, sc15, which displayed the highest mesenchymal programming and corresponded to the C6 signature in our Xenium data, was associated with a more favorable prognosis (log-rank P = 0.01, HR = 0.836, 95% CI = 0.727–0.961). Discussion In this study, we tracked the entire metastatic cascade of breast cancer at single-cell resolution and identified the relationship between spatial evolution of tumour cells and dynamic heterogeneity within the EMT spectrum. Our capture of early-stage metastatic events, including exceptionally early dissemination to lymph nodes, provides a rare in vivo snapshot of EMT dynamics during initial colonization. Thus, we were able to uncover distinct subpopulations with specific biological roles and their potential clinical relevance. To our knowledge, this is the first report to detail such fine-grained EMT profiles in vivo during early metastasis. Through detailed analysis of pioneering metastatic tumour cells within metastatic cascade, we identified six distinct types of hybrid epithelial–mesenchymal (EM) cells that play key roles in the early stages of metastasis. The subpopulation that underwent a lipid metabolic shift successfully formed colonies in lymph nodes. This is consistent with previous studies indicating that FAO supports lymph node metastasis in breast cancer 8 . Interestingly, the subpopulation responsible for forming lymph node colonies that underwent a lipid metabolic shift, as well as the subpopulation with the most mesenchymal features, was not associated with poor clinical prognosis. This supports previous evidence that mesenchymal-like phenotypes may correlate with favorable outcomes in patients with breast cancer 60 . Instead, poor prognosis was linked to transcriptional programs associated with MYC and E2F signaling, along with increased aerobic glycolysis and G2/M cell cycle activity. In particular, MYC reprograms metabolism toward enhanced glycolytic flux 61 , and breast cancers with enriched G2M signaling have been associated with elevated MYC pathway activity, increased proliferation, early distant metastasis, and worse survival 46 . Collectively, these findings suggest that distinct EMT subpopulations drive different aspects of tumor progression: some mediate niche colonization, whereas others may drive aggressive metastasis and recurrence, implying the need for tailored therapeutic strategies 11 . A key question concerns how hybrid EM cells disseminate in lymph nodes: do metastases arise via the sequential model, where fully EMT-induced M cells colonize distant sites by MET, or via the cooperative model, which mesenchymal cells support Epithelial cells, which serve as metastasis-initiating cells 62 ? Both models are considered plausible based on prior studies. In the present study, the observation of multiple EMT states in ITCs at an extremely early metastatic stage (Fig. 3a), along with similar mixtures of hybrid EM cells in the primary tumor (Fig. 3h), supports the cooperative metastasis model, where different cell states collaborate to initiate metastasis. We also examined the early interactions between metastatic tumor cells and the lymph node microenvironment. Comparisons between the TdL and metastatic lymph nodes showed that a coordinated, immunosuppressive environment forms around tumor cells in newly colonized niches. These environments are likely shaped, in part, by the metabolic phenotypes of the tumor cells themselves. For example, lactic acid, a byproduct of glycolysis, has been shown to promote alternative splicing in T cells and enhance CTLA-4 expression in a Foxp3-dependent manner 63 . The plasticity of EMT observed in this study may contribute to the ecosystem that supports metastasis outgrowth, including immunomodulatory effects 64 . Despite the novel insights in this study, it has several limitations. First, our in-depth analysis was conducted on tissue from a single patient, limiting the generalizability of our conclusions. Second, owing to the limited number of cells available, we were unable to perform molecular or biochemical validations. However, by integrating our findings with prior knowledge and publicly available single-cell RNA-seq datasets, we were able to infer the likely roles of key subpopulations and validate their clinical relevance. Finally, we could not definitively distinguish whether the observed metastasis occurred via hematogenous or lymphatic dissemination. Future studies should also aim to elucidate whether distinct molecular features exist between these metastatic pathways. In summary, our findings highlight therapeutic vulnerabilities embedded within the EMT spectrum, particularly involving metabolic plasticity and niche adaptation. These insights offer promising avenues for precision diagnostics and the development of targeted therapies aimed at halting metastasis at its earliest and most tractable stage. Materials and Methods Clinical sample This study was approved by the Clinical Ethics Committee of St. Marianna University School of Medicine (approval number: 2297-i103). Informed consent was obtained from the patient. The surgical specimen was from an 80-year-old female patient, who received surgery for breast cancer (BRC-26). The patient had not received neoadjuvant or adjuvant chemotherapy (detailed clinicopathologic findings are provided in Supplementary Table 1). The same specimen was used in our previous study 65 . For multiregional, multi-omics analyses, fresh frozen samples of the primary tumor and adjacent normal tissue were used. Two nontumor regions (NC1 and NC2) and two tumor regions (Ca3, a solid tumor region, and Ca4, a sparse tumor region), as well as normal breast tissue from the same specimen (used as the control), were microdissected using the AVENIO MilliSect system (Roche, Pleasanton, CA) and subjected to multi-omics analysis, including RNA-seq, WGS, and EM-seq. For spatial analyses via Xenium, formalin-fixed paraffin-embedded (FFPE) samples of the primary tumor, lymph node metastases, and TdL were used. RNA-seq RNA-seq was conducted as previously described 65 . Briefly, total RNA was extracted from frozen normal breast tissue using the RNeasy Micro Kit (Qiagen). RNA-seq libraries were prepared using the SMART-Seq Stranded Kit (Takara Bio). Paired-end 150-bp sequencing was performed on the NovaSeq 6000 system (Illumina). RNA-seq data analysis RNA-seq data from the normal tissue obtained in this study and four tumor regions from our previous study 65 were analyzed. Adapter trimming was performed using fastp (v0.23.2) 66 . Reads mapping to rRNA were removed using Bowtie 2 (v2.3.4.3) 67 . The retained reads were aligned to the human reference genome GRCh38.p12 using STAR (v2.7.5c) 68 . Gene-level read counts were obtained using featureCounts (v2.0.2) 69 , and RPKM values were calculated. DEGs between histological classes were identified using DESeq2 (v1.42.1) 70 applying the Wald test with an adjusted P-value cutoff of <0.1. Gene enrichment analysis was performed using Metascape 71 . Preparation of EM-seq and WGS libraries EM-seq and WGS library preparations were conducted as described previously 72 . Briefly, gDNA was extracted from normal tissue using NucleoSpin Tissue XS (MACHEREY-NAGEL), and 100 ng of this gDNA was fragmented using the M220 Focused-ultrasonicator (Covaris). Adapter ligation was conducted using the NEBNext Enzymatic Methyl-seq Kit (New England BioLabs). Half of the adapter-ligated DNA was used for WGS library preparation (five cycles PCR amplification), and the other half was employed for EM-seq library preparation, which included TET oxidation, APOBEC conversion, and PCR amplification (six cycles). Libraries were sequenced as paired-end 150-bp reads on the NovaSeq 6000 system (Illumina). EM-seq data analysis EM-seq data from the normal tissue obtained in this study and four tumor regions from our previous study were analyzed following established protocols 65 . Briefly, adapter trimming was conducted using Trim Galore (v0.6.4_dev; https://github.com/FelixKrueger/TrimGalore). Trimmed reads were aligned to the human reference genome using Bismark (v0.22.1) 73 . Duplicate reads were removed using deduplicate_bismark, and methylation information in a CpG context was extracted via bismark_methylation_extractor. Genome-wide methylation profiles, including patterns in CpG islands, CpG shores, and promoter regions, were visualized as described previously 74 . Differentially methylated regions (DMRs) among NC-1, NC-2, Ca-3, and Ca-4 and the broader NC and Ca regions were identified using metilene v0.2-8 75 , as per our previous study 65 . Promotor methylation was analyzed among NC-1, NC-2, Ca-3, and Ca-4 spanning 1 kb upstream to 500 bp downstream of transcription start sites, focusing on DMRs with methylation differences >10%. Visualisation was performed using IGV 76 . Point mutation detection WGS reads were aligned to the human reference genome using BWA-MEM (v0.7.17) 77 with default settings. Duplicate reads were marked using Picard MarkDuplicates (v2.23.8) (https://broadinstitute.github.io/picard/). Somatic mutations were called using GATK Mutect2 and filtered with FilterMutectCalls (v4.1.3.0) 78 . Variants were annotated via ANNOVAR 79 . CNV detection CNVs in NC-1, NC-2, Ca-3, Ca-4, and normal samples were detected using FACETS (v0.6.2), as described previously 74 . Copy number gains and losses were defined as ≥4 copies and ≤1 copy, respectively. Spatial heatmaps of gene expression levels and average methylation rates Spatial heatmaps of gene expression levels and average methylation rates were generated using a custom workflow implemented on Nikon’s pilot analysis platform. Specifically, microscopic images (Olympus BX53, 2× objective) with substantial overlap were aligned using similarity transformation based on scale-invariant feature transform–detected feature point pairs. This alignment produced a composite image of the entire specimen with a common coordinate system. Manually dissected region contours were aligned to this coordinate system. Regions of interest were color-coded by gene expression levels or average methylation rates. The average methylation rate was calculated as the mean frequency of methylation at CpG sites covered by at least five reads within the predefined promoter regions. In situ gene expression analysis via Xenium Spatial subclonal analysis at single-cell resolution was performed using the Xenium Slides & Sample Prep Reagents (10× Genomics), a predesigned human breast panel, and a custom panel, as previously described 80 . Briefly, FFPE tissue sections (5 µm thickness) of the primary tumor, lymph node metastases, and TdL were mounted onto Xenium slides. Following deparaffinization and decrosslinking, probe hybridization, probe circularization, and rolling circle amplification were performed. Detection of amplified probes was conducted using the Xenium Analyzer (10× Genomics). In total, 380 target genes, comprising 280 from the predesigned panel and 100 from the custom panel, are listed in Supplementary Table 1. Computational processing of Xenium in situ expression data and analysis Raw output files from the Xenium Analyzer were processed using Seurat (v5.0.2) 81 . Data normalization was performed using the SCTransform method. Clustering analysis and UMAP visualization were conducted using the first 30 principal components, with the former performed via the FindNeighbors and FindClustersing functions. The same parameters were applied to subclustering analyses (Figs. 2c and 3). For each cluster relevant to downstream analysis, DEGs were identified using the FindMarkers() function in Seurat with the default parameters (only.pos = TRUE, min.pct = 0.25, logfc.threshold = 0.25). Cluster annotation was based on previously published gene expression signatures 37-39 and spatial histological context. Trajectory analysis Trajectory and pseudotime analysis of luminal cells was performed using Monocle 3 (v1.3.4) 41 . Preprocessing was performed using the preprocess_cds function, with the number of dimensions set to 100. Dimensionality reduction and clustering were implemented using the reduce_dimension and cluster_cells functions, respectively. Principal graph learning was conducted with the learn_graph function. Based on spatial and Xenium data, mature luminal cell 2 was designated as the trajectory root. Pseudotime ordering was visualized using the plot_cells function. Identifying enriched gene signatures in Xenium clusters To investigate the functional characteristics of Xenium clusters (C1–C6 in Fig. 3e), pathway and gene set enrichment analyses were performed using Metascape 71 based on the top 50 DEGs per cluster (Supplementary Table 5). Enrichment was assessed across GO categories, including Biological Process, Cellular Component, and Molecular Function, as well as other biological pathways. Analysis of intercellular communication networks Intercellular communication within metastatic lymph nodes and the TdL was analyzed using the CellChat package (v1.6.1) 82 , which quantitatively identifies cell–cell interactions and communications. Following the official pipeline (https://github.com/jinworks/CellChat), analyses were performed using default settings. Statistically significant ligand–receptor interactions (P < 0.05) were extracted. Visualization of cell–cell interactions was achieved using circle plots, bar plots, heatmaps, and violin plots. Published dataset processing We downloaded a publicly available single-cell RNA-seq dataset (GSE180286 from Guan et al. 4 ). Data integration and clustering were performed in Seurat (v5.0.2) 81 , focusing on HER2-positive breast tumors and two lymph node samples (GSM5457205, GSM5457206, and GSM5457207). Louvain clustering at a resolution of 0.1 yielded nine clusters: one T cell cluster, three B cell clusters, one NK cell cluster, one CAF1 cluster, and one macrophage cluster. Clusters of the same cell type were merged (Fig. 5a). Subclustering was performed on the cancer cell populations (Fig. 5b). Clinical validation using METABRIC data Clinical validation was conducted using the METABRIC dataset 59 . METABRIC gene transcriptome data, as well as clinical and sample level metadata, were downloaded from cBioPortal (https://www.cbioportal.org). Analysis was performed on 1,980 breast cancer cases with available survival outcomes. For each cancer cluster defined in Fig. 5c, the top 100 DEGs (Supplementary Table 6) were used to score individual METABRIC samples. ROC curve analysis was used to define optimal cutoffs for stratifying patients into high and low score groups. Disease-free survival was analyzed using the survival (v3.5-7) and survminer (v0.4.9) R packages. HRs and 95% CIs were calculated using Cox proportional hazards models. Statistical significance was evaluated via the log-rank test. Statistical Analysis Statistical methods and tests are detailed in the figure legends. All analyses were performed using R (v4.3.2) or Python. Declarations Data availability The Xenium and bulk multi-omics data supporting the findings of this study are available at the DDBJ Japanese Genotype-phenotype Archive (https://gr-sharingdbs.dbcls.jp) under accession number JGAD000946. The single-cell RNA-seq dataset from Guan et al. 4 was downloaded from Gene Expression Omnibus (GSE180286). METABRIC transcriptome and clinical data were downloaded from cBioPortal (https://www.cbioportal.org). Additional data are available from the corresponding authors upon reasonable request. Acknowledgements: We thank K. Imamura, K. Abe, M. Satake, J. Zenkou, E. Sekimori, R. Fujinaga, and A. Gouda for their technical assistance. The authors would like to thank Enago (www.enago.jp) for the English language review. Funding: This work was supported by the Japan Agency for Medical Research and Development (AMED) P-PROMOTE, grant number JP21ck0106700 (to S.N.), the Japan Society for the Promotion of Science Grants-in-Aid for Scientific Research (KAKENHI Grants) Number JP24K11738 (Platform for Advanced Genome Science), as well as by Nikon Corporation. The supercomputing resource was provided by Human Genome Center, the University of Tokyo (http://sc.hgc.jp/shirokane.html). Computational resources were also supplied by Kashiwa-no-ha Omics Gate (https://www.kog.or.jp/en/server.html) Author contributions: S.N. and K.K collected the human samples, performed experiments and computational analyses, and generated all figures and tables. A.M., T.I., T.O., and K.T. aided in human sample collection. T.O. generated Fig. 2a. E.I. performed computational analyses. M.K. and M.N. assisted in the Xenium experiments. A.H., I.M. and J.K. conducted the pathological review. A.M. and K.T. developed viewer software. A.K. and A.S. offered advice and reviewed the manuscript. S.N., M.S., and Y.S. designed the project and wrote the manuscript. All authors read and approved the final manuscript. Competing interests: A.M. and K.T. are employees of Nikon Corporation. References Veronesi, U., et al. A randomized comparison of sentinel-node biopsy with routine axillary dissection in breast cancer. N Engl J Med . 349 , 546-553 (2003). Turner, R. R., Ollila, D. W., Krasne, D. L. & Giuliano, A. E. Histopathologic validation of the sentinel lymph node hypothesis for breast carcinoma. Ann Surg . 226 , 271-276; discussion 276-278 (1997). Nagasawa, S., et al. Genomic profiling reveals heterogeneous populations of ductal carcinoma in situ of the breast. Commun Biol . 4 , 438 (2021). Xu, K., et al. Microenvironment components and spatially resolved single-cell transcriptome atlas of breast cancer metastatic axillary lymph nodes. Acta Biochim Biophys Sin (Shanghai) . 54 , 1336-1348 (2022). Pesce, C. & Morrow, M. The need for lymph node dissection in nonmetastatic breast cancer. Annu Rev Med . 64 , 119-129 (2013). de Boer, M., van Dijck, J. A., Bult, P., Borm, G. F. & Tjan-Heijnen, V. C. Breast cancer prognosis and occult lymph node metastases, isolated tumor cells, and micrometastases. J Natl Cancer Inst . 102 , 410-425 (2010). Hinshaw, D. C. & Shevde, L. A. The Tumor Microenvironment Innately Modulates Cancer Progression. Cancer Res . 79 , 4557-4566 (2019). Lee, C. K., et al. Tumor metastasis to lymph nodes requires YAP-dependent metabolic adaptation. Science . 363 , 644-649 (2019). Cohen, M. R. Check concentrations of medications ordered for children. Nursing . 16 , 8 (1986). Nieto, M. A., Huang, R. Y., Jackson, R. A. & Thiery, J. P. Emt: 2016. Cell . 166 , 21-45 (2016). Puisieux, A., Brabletz, T. & Caramel, J. Oncogenic roles of EMT-inducing transcription factors. Nat Cell Biol . 16 , 488-494 (2014). Rubtsova, S. N., Zhitnyak, I. Y. & Gloushankova, N. A. Dual role of E-cadherin in cancer cells. Tissue Barriers . 10 , 2005420 (2022). Fischer, K. R., et al. Epithelial-to-mesenchymal transition is not required for lung metastasis but contributes to chemoresistance. Nature . 527 , 472-476 (2015). Ye, X., et al. Upholding a role for EMT in breast cancer metastasis. Nature . 547 , E1-E3 (2017). Diepenbruck, M. & Christofori, G. Epithelial-mesenchymal transition (EMT) and metastasis: yes, no, maybe? Curr Opin Cell Biol . 43 , 7-13 (2016). Aiello, N. M., et al. Upholding a role for EMT in pancreatic cancer metastasis. Nature . 547 , E7-E8 (2017). Pastushenko, I., et al. Identification of the tumour transition states occurring during EMT. Nature . 556 , 463-468 (2018). Kroger, C., et al. Acquisition of a hybrid E/M state is essential for tumorigenicity of basal breast cancer cells. Proc Natl Acad Sci U S A . 116 , 7353-7362 (2019). Tsai, J. H., Donaher, J. L., Murphy, D. A., Chau, S. & Yang, J. Spatiotemporal regulation of epithelial-mesenchymal transition is essential for squamous cell carcinoma metastasis. Cancer Cell . 22 , 725-736 (2012). Haerinck, J., Goossens, S. & Berx, G. The epithelial-mesenchymal plasticity landscape: principles of design and mechanisms of regulation. Nat Rev Genet . 24 , 590-609 (2023). Bian, X., et al. Carboxylic acid metabolism in cancer: Mechanisms, microenvironment interactions, and therapeutic opportunities. Metabolism . 171 , 156334 (2025). Vymetalkova, V., Vodicka, P., Vodenkova, S., Alonso, S. & Schneider-Stock, R. DNA methylation and chromatin modifiers in colorectal cancer. Mol Aspects Med . 69 , 73-92 (2019). Ehrlich, M. DNA hypomethylation in cancer cells. Epigenomics . 1 , 239-259 (2009). Rodriguez, J., et al. Chromosomal instability correlates with genome-wide DNA demethylation in human primary colorectal cancers. Cancer Res . 66 , 8462-9468 (2006). Usman, S., et al. Vimentin Is at the Heart of Epithelial Mesenchymal Transition (EMT) Mediated Metastasis. Cancers (Basel) . 13 (2021). Wang, Y., Shi, J., Chai, K., Ying, X. & Zhou, B. P. The Role of Snail in EMT and Tumorigenesis. Curr Cancer Drug Targets . 13 , 963-972 (2013). DaSilva-Arnold, S. C., et al. ZEB2, a master regulator of the epithelial-mesenchymal transition, mediates trophoblast differentiation. Mol Hum Reprod . 25 , 61-75 (2019). Sekimata, M., Murakami-Sekimata, A. & Homma, Y. CpG methylation prevents YY1-mediated transcriptional activation of the vimentin promoter. Biochem Biophys Res Commun . 414 , 767-772 (2011). Chen, Y., Wang, K., Qian, C. N. & Leach, R. DNA methylation is associated with transcription of Snail and Slug genes. Biochem Biophys Res Commun . 430 , 1083-1090 (2013). Xie, H., et al. Significance of ZEB2 in the immune microenvironment of colon cancer. Front Genet . 13 , 995333 (2022). Valdivia, A., Avalos, A. M. & Leyton, L. Thy-1 (CD90)-regulated cell adhesion and migration of mesenchymal cells: insights into adhesomes, mechanical forces, and signaling pathways. Front Cell Dev Biol . 11 , 1221306 (2023). Brenet, M., et al. Thy-1 (CD90)-Induced Metastatic Cancer Cell Migration and Invasion Are beta3 Integrin-Dependent and Involve a Ca(2+)/P2X7 Receptor Signaling Axis. Front Cell Dev Biol . 8 , 592442 (2020). van Zijl, F., et al. Hepatic tumor-stroma crosstalk guides epithelial to mesenchymal transition at the tumor edge. Oncogene . 28 , 4022-4033 (2009). Sethi, S., Macoska, J., Chen, W. & Sarkar, F. H. Molecular signature of epithelial-mesenchymal transition (EMT) in human prostate cancer bone metastasis. Am J Transl Res . 3 , 90-99 (2010). Tinsley, E., Bredin, P., Toomey, S., Hennessy, B. T. & Furney, S. J. KMT2C and KMT2D aberrations in breast cancer. Trends Cancer . 10 , 519-530 (2024). Nishimura, T., et al. Evolutionary histories of breast cancer and related clones. Nature . 620 , 607-614 (2023). Aran, D., Hu, Z. & Butte, A. J. xCell: digitally portraying the tissue cellular heterogeneity landscape. Genome Biol . 18 , 220 (2017). Lim, E., et al. Aberrant luminal progenitors as the candidate target population for basal tumor development in BRCA1 mutation carriers. Nat Med . 15 , 907-913 (2009). Pliner, H. A., Shendure, J. & Trapnell, C. Supervised classification enables rapid annotation of cell atlases. Nat Methods . 16 , 983-986 (2019). Nguyen, Q. H., et al. Profiling human breast epithelial cells using single cell RNA sequencing identifies cell diversity. Nat Commun . 9 , 2028 (2018). Qiu, X., et al. Reversed graph embedding resolves complex single-cell trajectories. Nat Methods . 14 , 979-982 (2017). Lamouille, S., Xu, J. & Derynck, R. Molecular mechanisms of epithelial-mesenchymal transition. Nat Rev Mol Cell Biol . 15 , 178-196 (2014). Zeng, Y. A. & Nusse, R. Wnt proteins are self-renewal factors for mammary stem cells and promote their long-term expansion in culture. Cell Stem Cell . 6 , 568-577 (2010). Zhang, Y. C., Zhang, Y. T., Wang, Y., Zhao, Y. & He, L. J. What role does PDL1 play in EMT changes in tumors and fibrosis? Front Immunol . 14 , 1226038 (2023). Imodoye, S. O. & Adedokun, K. A. EMT-induced immune evasion: connecting the dots from mechanisms to therapy. Clin Exp Med . 23 , 4265-4287 (2023). Oshi, M., et al. G2M Cell Cycle Pathway Score as a Prognostic Biomarker of Metastasis in Estrogen Receptor (ER)-Positive Breast Cancer. Int J Mol Sci . 21 (2020). Graziano, F., et al. Glycolysis gene expression analysis and selective metabolic advantage in the clinical progression of colorectal cancer. Pharmacogenomics J . 17 , 258-264 (2017). Sacco, J. L. & Gomez, E. W. Epithelial-Mesenchymal Plasticity and Epigenetic Heterogeneity in Cancer. Cancers (Basel) . 16 (2024). Teeuwssen, M. & Fodde, R. Cell Heterogeneity and Phenotypic Plasticity in Metastasis Formation: The Case of Colon Cancer. Cancers (Basel) . 11 (2019). Keller, H., et al. Fatty acids and retinoids control lipid metabolism through activation of peroxisome proliferator-activated receptor-retinoid X receptor heterodimers. Proc Natl Acad Sci U S A . 90 , 2160-2164 (1993). Galimberti, V., et al. Axillary dissection versus no axillary dissection in patients with sentinel-node micrometastases (IBCSG 23-01): a phase 3 randomised controlled trial. Lancet Oncol . 14 , 297-305 (2013). Martinez-Pomares, L. The mannose receptor. J Leukoc Biol . 92 , 1177-1186 (2012). Schuette, V., et al. Mannose receptor induces T-cell tolerance via inhibition of CD45 and up-regulation of CTLA-4. Proc Natl Acad Sci U S A . 113 , 10649-10654 (2016). Luca, M. E., et al. Mannosylated PLP(139-151) induces peptide-specific tolerance to experimental autoimmune encephalomyelitis. J Neuroimmunol . 160 , 178-187 (2005). Kel, J., et al. Soluble mannosylated myelin peptide inhibits the encephalitogenicity of autoreactive T cells during experimental autoimmune encephalomyelitis. Am J Pathol . 170 , 272-280 (2007). Choi, J. H., et al. Single-cell transcriptome profiling of the stepwise progression of head and neck cancer. Nat Commun . 14 , 1055 (2023). Wang, T., et al. COL4A1 promotes the growth and metastasis of hepatocellular carcinoma cells by activating FAK-Src signaling. J Exp Clin Cancer Res . 39 , 148 (2020). Han, J., et al. CXCL16 Promotes Gastric Cancer Tumorigenesis via ADAM10-Dependent CXCL16/CXCR6 Axis and Activates Akt and MAPK Signaling Pathways. Int J Biol Sci . 17 , 2841-2852 (2021). Curtis, C., et al. The genomic and transcriptomic architecture of 2,000 breast tumours reveals novel subgroups. Nature . 486 , 346-352 (2012). Grosse-Wilde, A., et al. Loss of inter-cellular cooperation by complete epithelial-mesenchymal transition supports favorable outcomes in basal breast cancer patients. Oncotarget . 9 , 20018-20033 (2018). Dong, Y., Tu, R., Liu, H. & Qing, G. Regulation of cancer cell metabolism: oncogenic MYC in the driver's seat. Signal Transduct Target Ther . 5 , 124 (2020). Tsuji, T., Ibaragi, S. & Hu, G. F. Epithelial-mesenchymal transition and cell cooperativity in metastasis. Cancer Res . 69 , 7135-7139 (2009). Ding, R., et al. Lactate modulates RNA splicing to promote CTLA-4 expression in tumor-infiltrating regulatory T cells. Immunity . 57 , 528-540 e526 (2024). Reticker-Flynn, N. E., et al. Lymph node colonization induces tumor-immune tolerance to promote distant metastasis. Cell . 185 , 1924-1942 e1923 (2022). Keisuke Kunigo, et al. Deep long-read methylation analysis using targeted nanopore sequencing and application to local region analysis of clinical specimen, 2025). Chen, S., Zhou, Y., Chen, Y. & Gu, J. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics . 34 , i884-i890 (2018). Langmead, B. & Salzberg, S. L. Fast gapped-read alignment with Bowtie 2. Nat Methods . 9 , 357-359 (2012). Dobin, A., et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics . 29 , 15-21 (2013). Liao, Y., Smyth, G. K. & Shi, W. featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics . 30 , 923-930 (2014). Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol . 15 , 550 (2014). Zhou, Y., et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun . 10 , 1523 (2019). Sakamoto, Y., et al. Long-read whole-genome methylation patterning using enzymatic base conversion and nanopore sequencing. Nucleic Acids Res . 49 , e81 (2021). Krueger, F. & Andrews, S. R. Bismark: a flexible aligner and methylation caller for Bisulfite-Seq applications. Bioinformatics . 27 , 1571-1572 (2011). Haga, Y., et al. Whole-genome sequencing reveals the molecular implications of the stepwise progression of lung adenocarcinoma. Nat Commun . 14 , 8375 (2023). Juhling, F., et al. metilene: fast and sensitive calling of differentially methylated regions from bisulfite sequencing data. Genome Res . 26 , 256-262 (2016). Thorvaldsdóttir, H., Robinson, J. T. & Mesirov, J. P. Integrative Genomics Viewer (IGV): high-performance genomics data visualization and exploration. Brief Bioinform . 14 , 178-192 (2013). Li, H. & Durbin, R. Fast and accurate long-read alignment with Burrows-Wheeler transform. Bioinformatics . 26 , 589-595 (2010). McKenna, A., et al. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res . 20 , 1297-1303 (2010). Wang, K., Li, M. & Hakonarson, H. ANNOVAR: functional annotation of genetic variants from high-throughput sequencing data. Nucleic Acids Res . 38 , e164 (2010). Takano, Y., et al. Spatially resolved gene expression profiling of tumor microenvironment reveals key steps of lung adenocarcinoma development. Nat Commun . 15 , 10637 (2024). Hao, Y., et al. Dictionary learning for integrative, multimodal and scalable single-cell analysis. Nat Biotechnol . 42 , 293-304 (2024). Jin, S., et al. Inference and analysis of cell-cell communication using CellChat. Nat Commun . 12 , 1088 (2021). Additional Declarations Competing interest reported. A.M. and K.T. are employees of Nikon Corporation. Supplementary Files SupplementaryFiguresPDF.pdf SupplementaryTables.xlsx Cite Share Download PDF Status: Published Journal Publication published 23 Jan, 2026 Read the published version in npj Breast Cancer → Version 1 posted Editorial decision: Revision requested 22 Nov, 2025 Reviews received at journal 11 Nov, 2025 Reviews received at journal 03 Oct, 2025 Reviewers agreed at journal 19 Sep, 2025 Reviewers agreed at journal 15 Sep, 2025 Reviewers invited by journal 14 Sep, 2025 Editor assigned by journal 20 Aug, 2025 Submission checks completed at journal 20 Aug, 2025 First submitted to journal 14 Aug, 2025 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-7371816","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":517713491,"identity":"d220beab-6115-40c8-ac7e-e38592de550b","order_by":0,"name":"Satoi Nagasawa","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAABPElEQVRIie2QMUvDQBTH3xFIlitZL0Tar3DhICCt+FVSArqkuDlJiAQySbO2oN9BEIKbkYO6nLrGLUWoS4QERAoFsZoGhx7OgvlN/3uPH/feA2hp+YvoQZOUFDAArR/ppkgkCmmaoDqNgoJfFXCagKlEkdAL7m+rtxO/R+9Elb9GwIA4Wl4K3gWN57B7vaVQPnbNnRm3LsXoyrqIwF4r6HSScQb4gIIhthUQ1CRqipJ0lJidCAZAjqoQl3wYgAdgRJLBPLYiH/5+8li8rGrFQbWiF1IFUs82qkgZJpmnKp3NYCHO1gqR/0IzYffRmLs32YIZ5w+EqThH04k4ZCpZ0FSySy8+Y0/Ld39vGrvzsjgeWLHmQFnO+l1dd+fPxvbFvlDwTyag1iH8DtwIpApaSop+3arkSktLS8u/4hPXbm6okXdtxgAAAABJRU5ErkJggg==","orcid":"","institution":"The University of Tokyo","correspondingAuthor":true,"prefix":"","firstName":"Satoi","middleName":"","lastName":"Nagasawa","suffix":""},{"id":517713492,"identity":"3fbade8f-4eda-4e6c-8215-3591adc17a0d","order_by":1,"name":"Keiko Kajiya","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Keiko","middleName":"","lastName":"Kajiya","suffix":""},{"id":517713493,"identity":"e2ebe71d-ad24-44b0-9081-0f7f808a0f1b","order_by":2,"name":"Erina Ishikawa","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Erina","middleName":"","lastName":"Ishikawa","suffix":""},{"id":517713494,"identity":"7a81aeb8-abc1-44ee-ac6e-225c37fef87d","order_by":3,"name":"Akinori Kanai","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Akinori","middleName":"","lastName":"Kanai","suffix":""},{"id":517713495,"identity":"4eb26b6f-fd85-4b10-ad78-18f0d449a6fe","order_by":4,"name":"Ayako Suzuki","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Ayako","middleName":"","lastName":"Suzuki","suffix":""},{"id":517713496,"identity":"d35046ff-b459-41df-83bb-089320bfcdb3","order_by":5,"name":"Ai Motoyoshi","email":"","orcid":"","institution":"St. Marianna University School of Medicine","correspondingAuthor":false,"prefix":"","firstName":"Ai","middleName":"","lastName":"Motoyoshi","suffix":""},{"id":517713497,"identity":"2fd8fa72-41cc-4248-8017-3f022c653182","order_by":6,"name":"Tsuguo Iwatani","email":"","orcid":"","institution":"St. Marianna University School of Medicine","correspondingAuthor":false,"prefix":"","firstName":"Tsuguo","middleName":"","lastName":"Iwatani","suffix":""},{"id":517713498,"identity":"bf1e8f08-7045-4542-91fb-139d3fcf63d0","order_by":7,"name":"Manabu Kubota","email":"","orcid":"","institution":"St. Marianna University School of Medicine","correspondingAuthor":false,"prefix":"","firstName":"Manabu","middleName":"","lastName":"Kubota","suffix":""},{"id":517713500,"identity":"81037304-cad8-44ed-843c-856124e6d0c3","order_by":8,"name":"Masaru Nakamura","email":"","orcid":"","institution":"St. Marianna University School of Medicine","correspondingAuthor":false,"prefix":"","firstName":"Masaru","middleName":"","lastName":"Nakamura","suffix":""},{"id":517713505,"identity":"0f0c95e1-c0f5-48f3-82d9-51703e55c9e4","order_by":9,"name":"Tatsuya Onishi","email":"","orcid":"","institution":"National Cancer Center Hospital East","correspondingAuthor":false,"prefix":"","firstName":"Tatsuya","middleName":"","lastName":"Onishi","suffix":""},{"id":517713506,"identity":"93d8e0d3-9af3-4471-a8de-656927760349","order_by":10,"name":"Akiyoshi Hoshino","email":"","orcid":"","institution":"Kitasato University Kitasato Institute Hospital","correspondingAuthor":false,"prefix":"","firstName":"Akiyoshi","middleName":"","lastName":"Hoshino","suffix":""},{"id":517713507,"identity":"0a261853-9103-4ad4-b732-714ec58b746c","order_by":11,"name":"Ichiro Maeda","email":"","orcid":"","institution":"Kitasato University Kitasato Institute Hospital","correspondingAuthor":false,"prefix":"","firstName":"Ichiro","middleName":"","lastName":"Maeda","suffix":""},{"id":517713508,"identity":"25e23cb4-b85c-416f-884a-bcf4a37702f8","order_by":12,"name":"Akihiko Morozumi","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Akihiko","middleName":"","lastName":"Morozumi","suffix":""},{"id":517713509,"identity":"c255e7cf-2909-4ca1-968d-cb58683eb058","order_by":13,"name":"Kenji Takatsuka","email":"","orcid":"","institution":"Nikon Corporation","correspondingAuthor":false,"prefix":"","firstName":"Kenji","middleName":"","lastName":"Takatsuka","suffix":""},{"id":517713510,"identity":"e8dda654-b29d-401e-8d7c-0ff9d2c2da3c","order_by":14,"name":"Junki Koike","email":"","orcid":"","institution":"St. Marianna University School of Medicine","correspondingAuthor":false,"prefix":"","firstName":"Junki","middleName":"","lastName":"Koike","suffix":""},{"id":517713511,"identity":"75e20608-d987-4c1d-8a4e-b7a686b84b8c","order_by":15,"name":"Masahide Seki","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Masahide","middleName":"","lastName":"Seki","suffix":""},{"id":517713512,"identity":"2d749016-74d3-4ff8-b225-b59ce9b00dac","order_by":16,"name":"Koichiro Tsugawa","email":"","orcid":"","institution":"St. Marianna University School of Medicine","correspondingAuthor":false,"prefix":"","firstName":"Koichiro","middleName":"","lastName":"Tsugawa","suffix":""},{"id":517713513,"identity":"a9658588-044c-4806-9296-97bbbd52333c","order_by":17,"name":"Yutaka Suzuki","email":"","orcid":"","institution":"The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Yutaka","middleName":"","lastName":"Suzuki","suffix":""}],"badges":[],"createdAt":"2025-08-14 08:53:22","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-7371816/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-7371816/v1","draftVersion":[],"editorialEvents":[{"content":"https://doi.org/10.1038/s41523-026-00897-1","type":"published","date":"2026-01-23T15:56:58+00:00"}],"editorialNote":"","failedWorkflow":false,"files":[{"id":91960289,"identity":"bae7a7bf-8e7c-4e42-911b-9b62d9543d8b","added_by":"auto","created_at":"2025-09-23 07:46:01","extension":"jpg","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":853779,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eMulti-omics analysis across regions of breast cancer tissue.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ea)\u003c/strong\u003eHE images of breast cancer tissue with regions sampled for multi-omics analysis. Two noncancerous regions (NC-1 and NC-2; pink lines) and two cancerous regions (Ca-3, a solid tumor region; Ca-4, a sparse tumor region; light blue lines) were dissected (also refer to Supplementary Fig. 1). Scale bar: 1 mm.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eb) \u003c/strong\u003eGene enrichment analysis using DEGs from noncancerous (upper panel) and cancerous regions (lower panel). Bar graphs show the top 10 enriched clusters; color scale indicates −log10 adjusted P-values.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ec)\u003c/strong\u003e Violin plots showing genome-wide (50 kb window) methylation rates measured by nanoEM across each region. Red line indicate mean.P-values were calculated using the two-sided Wilcoxon rank-sum test with Bonferroni correction. *P \u0026lt; 0.05, **P \u0026lt; 0.01, ***P \u0026lt; 0.001, ****P \u0026lt; 0.0001.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ed) \u003c/strong\u003eHeatmaps showing Z-score–normalized gene expression levels (from TPM values) and promoter CpG methylation levels of selected genes involved in cancer metastasis in Ca-3 and Ca-4. Promoter regions were defined as 1 kb upstream to 500 bp downstream of the transcription start site. Genes with the most variable promoter methylation across NC-1, NC-2, Ca-3, Ca-4, and normal breast tissue were selected.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ee) \u003c/strong\u003eSpatial visualization of THY1 gene expression (left) and promotor methylation (right) using our custom viewer. Color scale represents Z-scores (normalized from TPM) for gene expression and average methylation rate for methylation.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ef)\u003c/strong\u003e Methylation pattern of THY1 detected by EM-seq visualized using the IGV; the middle panel shows a heatmap of methylation rates across regions.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eg)\u003c/strong\u003e Mutation status of breast cancer-related genes across regions. Categories include MAPK-PI3K signaling, cell cycle regulation, DNA damage and apoptosis, chromatin remodeling, cell adhesion, transcription regulation, and clinically relevant genes.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eh) \u003c/strong\u003eRepresentative mutations in KMT2C (Y987H) and KDM6A (P281R) across regions.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ei) \u003c/strong\u003eGenome-wide copy number profiles by chromosome for each region, with cancer-related genes showing gains (red) or losses (blue) indicated in the inset.\u003c/p\u003e","description":"","filename":"Figure1.jpg","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/5135a1efe20dc19adba9959c.jpg"},{"id":91958584,"identity":"5a414de5-c52d-4b61-8cda-35188496ea8f","added_by":"auto","created_at":"2025-09-23 07:38:01","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":19918626,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eSpatial analysis of the breast cancer at single-cell resolution.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ea)\u003c/strong\u003eWorkflow of the study. Primary tumor (P), metastatic lymph node (ML), tumor-draining lymph node (TdL) samples were collected and processed for Xenium analysis. Middle panel: HE staining; right panel: annotated major cell types from b).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eb) \u003c/strong\u003eUMAP of 1,005,436 cells analyzed via Xenium Explore, annotated by cell type using canonical markers and histological location. Right panel shows subpopulation of luminal cells (n = 245,385 cells).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ec)\u003c/strong\u003e UMAP visualization of luminal cells ; right panel shows subclustered labeled as mature luminal 1 (L1), mature luminal 2 (L2), basal (B), and mesenchymal like (M) cells.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ed)\u003c/strong\u003e Histological location and gene expression in clusters representing normal lobules. Violin plot shows log-normalized expression levels of markers for mature luminal (ANKRD30A and KIT), basal (KRT23) cells. Scale bar: 100 μm.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ee)\u003c/strong\u003e Histological location of M clusters. Violin plot shows log-normalized expression of DEGs for M clusters; bubble plot shows EMT-related gene expression. Dot size represents the percentage of cells expressing. Color indicates mean expression relative to other subsets. Scale bar: 500 μm.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ef)\u003c/strong\u003e Histological location in cancer cluster(ii)v and expression of L2 marker ANKRD30A.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eg) \u003c/strong\u003eMonocle 3 trajectory analysis of luminal cell differentiation, revealing two main trajectories. The UMAP axes correspond to b). vCells are color-coded by pseudotime. Red winding lines in the objects indicate putative developmental trajectories.\u003c/p\u003e","description":"","filename":"Figure2.png","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/ca0bb8641af08961baafb81d.png"},{"id":91958579,"identity":"ff6815b0-12e0-4e46-9c44-52ac3ee8c2a2","added_by":"auto","created_at":"2025-09-23 07:38:01","extension":"jpg","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":1511341,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eSpatial analysis of lymph nodes reveals six transcriptionally distinct EMT states.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ea) \u003c/strong\u003eIsolated tumor cells (ITCs) within the TdL. Enlarged view (upper panel), KRT19 expression (middle panel), and HE staining (lower panel). Scale bar: 100 μm.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eb) \u003c/strong\u003eUMAP placement of ITCs and reclustering of tumor cells, cell numbers per cluster shown in footnotes.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ec) \u003c/strong\u003eHistological localization of cancer cell clusters.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ed)\u003c/strong\u003e Bubble plot showing expression of related to EMT transcription factors, adhesion molecules, and stemness markers across clusters; violin plot shows representative marker expression.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ee)\u003c/strong\u003e Pathway enrichment using top 50 DEGs (Supplementary Table 5) in clusters C1–C6; heatmap shows enrichment scores (gray: not significant).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ef) \u003c/strong\u003eXenium spatial plot of FASN and SCD transcripts. Color shows log-normalized expression levels.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eg) \u003c/strong\u003eViolin plot of G2M signature scores across cancer clusters using Seurat; P-values were determined via the Wilcoxon rank-sum test. *P \u0026lt; 0.05, **P \u0026lt; 0.01, ***P \u0026lt; 0.001, and ****P \u0026lt; 0.0001.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eh) \u003c/strong\u003eRepresentation of cancer clusters within the primary tumor; cell counts shown in footnotes.\u003c/p\u003e","description":"","filename":"Figure3.jpg","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/50e5ecee2b47456e8b8a7b80.jpg"},{"id":91958578,"identity":"0c074962-decd-4f02-ae11-64b64b58a9c4","added_by":"auto","created_at":"2025-09-23 07:38:01","extension":"jpg","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":510711,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eAnalysis of the TME in metastatic lymph nodes.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ea)\u003c/strong\u003eROI analyzed for cell–cell communication. Pie chart shows proportions of TME components. Scale bar: 100 μm.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eb)\u003c/strong\u003e Overview of intercellular interactions. Arrows and edge show directionality; color indicates cell type; circle size is proportional to cell count; edge thickness indicates the number of interactions between populations. Bar chart ranks signaling pathways by differential information flow between metastatic lymph nodes and the TdL (dark red shows enrichment in metastatic lymph nodes; dark green shows enrichment in the TdL).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ec)\u003c/strong\u003e Heatmaps showing CD45 signaling interactions in metastatic lymph nodes; color scale represents communication probabilities; empty squares show zero probability. Violin plot shows gene expression levels.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ed)\u003c/strong\u003e VEGF and CD86 signaling interactions in the TdL, shown as heatmaps and violin plots, as in c).\u003c/p\u003e","description":"","filename":"Figure4.jpg","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/ff43ec455010b6fd03876c9b.jpg"},{"id":91960290,"identity":"45a7f031-d41c-4850-b3f1-ef503418a858","added_by":"auto","created_at":"2025-09-23 07:46:01","extension":"jpg","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":702778,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eTranscriptional profiling of hybrid EMT and its association with patient outcomes.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ea)\u003c/strong\u003e UMAP annotated by dataset (left) and major cell types (right). Datasets: primary tumor (GSM5457205); lymph nodes (GSM5457206 and GSM5457207).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eb)\u003c/strong\u003e UMAP of reclustered tumor cells.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ec)\u003c/strong\u003e Dot plot compares metabolic pathway scores from MSigDB C2 gene sets across subclusters (sc2–sc15); color scale shows −log10-adjusted P-value.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ed)\u003c/strong\u003e Correlation of public tumor clusters (sc2, sc4, sc8, sc10, and sc15) with Xenium clusters (C1–C6). The top 50 DEGs (Supplementary Table 5) from Xenium clusters were used to generate signature scores. Violin and feature plots show Xenium signature scores. The UMAP axes correspond to b). C5 and C6 correspond to sc8 and sc15, respectively.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ee)\u003c/strong\u003e Survival analysis using the METABRIC cohort. Signature scores based on top 100 DEGs (Supplementary Table 6) for sc8 and sc15 were computed using UCell. Kaplan–Meier survival plots show disease-free survival stratified by high vs. low signature scores. P-values from log-rank test; adjusted HRs from multivariate Cox proportional hazards model, with the low-scores as the reference.\u003c/p\u003e","description":"","filename":"Figure5.jpg","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/3d92996d399b6a49595b44ba.jpg"},{"id":101153387,"identity":"8a64d901-5790-466f-89db-b8cb9e2fb11c","added_by":"auto","created_at":"2026-01-26 16:15:06","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":22466874,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/f5a88481-3d98-423a-8877-3e69221d7f60.pdf"},{"id":91958582,"identity":"4e40d6d7-ea81-402b-953e-bfcb2684aeb9","added_by":"auto","created_at":"2025-09-23 07:38:01","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"supplement","size":4148780,"visible":true,"origin":"","legend":"","description":"","filename":"SupplementaryFiguresPDF.pdf","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/2ff858b38e69edd94439698b.pdf"},{"id":91958577,"identity":"64bd0ad4-25a1-43a1-894b-ffafad78094d","added_by":"auto","created_at":"2025-09-23 07:38:01","extension":"xlsx","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":640574,"visible":true,"origin":"","legend":"","description":"","filename":"SupplementaryTables.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7371816/v1/87a7ac8d63e8f0944e52e351.xlsx"}],"financialInterests":"Competing interest reported. A.M. and K.T. are employees of Nikon Corporation.","formattedTitle":"Spatial Gene Expression Analysis Reveals Drivers of Extremely Early Lymph Node Metastasis in Breast Cancer","fulltext":[{"header":"Introduction","content":"\u003cp\u003eAxillary lymph nodes are considered the “first stop” in breast cancer metastasis \u003csup\u003e1-2\u003c/sup\u003e. The presence of metastatic tumor cells in these nodes is a well-recognized marker of poor prognosis \u003csup\u003e3-7\u003c/sup\u003e. Tumor cell profiles in metastases have been studied extensively. For example, Lee et al. used a mouse model and found that a shift toward fatty acid oxidation (FAO), driven by the selective activation of the transcription coactivator YAP in lymph node metastases, is essential for their establishment \u003csup\u003e8\u003c/sup\u003e. Liu et al. used single-cell RNA sequencing (RNA-seq) and spatial transcriptomics to compare the microenvironments of primary and lymph node metastatic tumors, reporting that in lymph node metastases, compared with primary tumors, T cell activation, cytotoxicity, and proliferation are markedly suppressed while dendritic cells show reduced ability to prime and activate T cells \u003csup\u003e9\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003eAlthough the requirements for lymph node metastasis and the niche environment are becoming clearer, how these molecular characteristics relate to cellular plasticity, possibly involving EMT and MET, remains uncertain. Notably, almost no studies have reported the first metastatic event in the lymph nodes. This is because it is difficult to track the trajectory of metastasis in humans, and it is extremely difficult to observe transient and plastic EMT states in vivo.\u003c/p\u003e\n\u003cp\u003eSubstantial evidence supports the role of epithelial-to-mesenchymal transition (EMT) and mesenchymal-to-epithelial transition (MET) in metastasis \u003csup\u003e10-11\u003c/sup\u003e.\u0026nbsp;EMT is a cellular process in which epithelial cells acquire mesenchymal cell characteristics, causing dramatic changes in tissue structure and function. Cells take on a fibroblast-like morphology and increase their motility (ability to move) and invasiveness (ability to invade surrounding tissues). MET is a process that reverses the changes in cellular phenotype induced by EMT. In cancer cells, suppression of EMT-Transcription factors\u0026nbsp;like\u0026nbsp;\u003cem\u003eTWIST1\u003c/em\u003e is necessary to promote MET, which is required for proliferation after metastasis \u003csup\u003e12\u003c/sup\u003e.Despite extensive debate, it remains unclear how EMT or MET contributes to metastatic potential and whether it is a necessary condition for metastasis \u003csup\u003e10,13-16\u003c/sup\u003e. A major challenge in addressing this question lies in EMT’s broad and evolving definition. EMT was once considered a binary switch between epithelial and mesenchymal states. However, several \u003cem\u003ein vitro\u003c/em\u003e studies have shown that EMT progresses gradually, forming hybrid subpopulations that occupy intermediate states between epithelial and mesenchymal identities \u003csup\u003e17-18\u003c/sup\u003e. Interestingly, tumor cell subpopulations nearest to the mesenchymal state do not necessarily exhibit the highest metastatic potential \u003csup\u003e19\u003c/sup\u003e. Indeed, high metastatic potential often lies in the hybrid state. Conversely, evidence suggests that metastatic potential and the reverse process, MET, are not correlated, indicating that EMT and MET may need to be considered separately \u003csup\u003e17\u003c/sup\u003e. To understand the molecular events driving metastasis, it is essential to examine how EMT progresses in the primary tumor by dissecting its intermediate states and how, and to what extent, metastasized cells restore epithelial traits\u0026nbsp; in the lymph node \u003csup\u003e20\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003eIn the present study, we used the Xenium in situ hybridization platform, which enables spatial gene expression profiling at the single-cell resolution, to track tumour cells in vivo during the metastasis process. The purpose of the study was to clarify the spatial evolution of tumor cells in the extremely early stages of breast cancer lymph nodes metastasis.\u003c/p\u003e"},{"header":"Results","content":"\u003cp\u003e\u003cstrong\u003eMulti-omics analysis in each breast cancer region\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eBefore analyzing spatial patterns at the primary site and metastatic lymph node, we first examined the molecular features of the primary cancer via multiregional bulk analysis. The patient, an 80-year-old female (BRC-26), underwent surgery without neoadjuvant therapy. Her clinical diagnosis was stage II invasive breast carcinoma of the HER2 type, i.e., ER-negative, PgR-negative, and HER2-positive (for full clinicopathologic details, refer to Supplementary Table 1). For high-quality molecular analysis, especially RNA-seq, we freshly froze part of the tumor tissue harvested during surgery. We roughly dissected samples from histologically noncancerous regions (NC regions; NC-1 and NC-2) and cancerous regions (Ca regions; Ca-3 and Ca-4) (Fig. 1a and Supplementary Fig. 1a). Histologically, Ca-4 appeared sparser and less adhesive compared with Ca-3 (Supplementary Fig. 1b). We subjected these tissues to RNA-seq, enzymatic methylation sequencing (EM-seq), and whole-genome sequencing (WGS) (refer to Supplementary Table 2 for statistical analyses).\u003c/p\u003e\n\u003cp\u003eWe first performed transcriptome analysis using RNA-seq to compare gene expression between the NC and Ca regions. Of 4,584 differentially expressed genes (DEGs), 2,248 and 2,336 were significantly upregulated and downregulated in the cancer region (Supplementary Fig. 1c). Gene pathway enrichment analysis of these DEGs showed that \u0026ldquo;cell adhesion\u0026rdquo; pathways were enriched in NC regions, whereas pathways for \u0026ldquo;carboxylic acid metabolism,\u0026rdquo; and \u0026ldquo;cell cycle\u0026rdquo; were enriched in Ca regions (Fig. 1b), suggesting progressive malignant transformation from NC to Ca regions \u003csup\u003e21\u003c/sup\u003e. DNA methylation analysis revealed more extensive genome-wide hypomethylation in Ca regions compared with NC regions (Fig. 1c). Furthermore, promoters, intergenic regions, and gene bodies all showed reduced methylation (Supplementary Fig. 1d). Hypomethylation of intergenic and intron regions, especially repetitive sequences, can promote chromosomal instability and mutations, potentially increasing cancer risk \u003csup\u003e22-24\u003c/sup\u003e. These results confirm that breast cells undergo malignant transformation as they progress from NC to Ca regions.\u003c/p\u003e\n\u003cp\u003eWithin Ca regions, we observed significant changes in key metastasis-related genes (Fig. 1d). Expression of EMT markers, such as Vimentin \u003csup\u003e25\u003c/sup\u003e, and transcription factors Snail and ZEB2, considered, known EMT regulators \u003csup\u003e26-27\u003c/sup\u003e, increased from Ca-3 to Ca-4. Correspondingly, DNA methylation levels declined from Ca-3 to Ca-4 (Supplementary Fig. 1e), consistent with previous reports that EMT-related genes are epigenetically regulated \u003csup\u003e28-30\u003c/sup\u003e. Similarly, methylation of the \u003cem\u003eTHY1\u003c/em\u003e gene, which was highly expressed in Ca-4, decreased from Ca-3 to Ca-4 (Fig. 1e, f). THY1, a GPI-anchored protein involved cell adhesion, migration, and polarity, can suppress lung metastasis when its integrin signaling is inhibited \u003csup\u003e31-32\u003c/sup\u003e. \u003cem\u003eTHY1\u003c/em\u003e, which regulates cancer cell migration and invasion, also appears epigenetically controlled. These molecular profiles indicate that malignant transformation from NC to Ca regions, along with DNA methylation differences within Ca regions, may underlie varying EMT phenotypes in tumor cells \u003csup\u003e33-34\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003eWGS results further supported cancer progression from NC to Ca regions. We detected mutations in KMT2C (Y987H), the most frequently mutated histone methyltransferase in breast cancer \u003csup\u003e35\u003c/sup\u003e, across all regions (Fig. 1g). PIK3CA (H1047R), a known cancer driver, was mutated in all regions except NC-1. Mutations in transcription factor FOXA1 (R262P) and histone demethylase KDM6A (P281R) were specific to Ca regions (Fig. 1h). Copy number variation (CNV) analysis revealed chromosomal abnormalities in Ca regions (Fig.1i). These findings suggest that tumor development occurred against a background of clonal mammary gland expansion \u003csup\u003e36\u003c/sup\u003e. Overall, multiregional analyses indicated that breast cancer follows a well-characterized molecular trajectory from NC to Ca regions. Importantly, these molecular features aligned closely with histological observations, especially regarding gene regulation.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eSpatial analysis of the breast cancer at single-cell resolution\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo investigate the molecular mechanisms of metastasis in detail, we performed spatial gene expression analysis at single-cell resolution on primary tumors and paired axillary metastatic lymph nodes from the patient shown in Fig. 1.\u0026nbsp;In addition to the primary tumor adjacent to the bulk analysis region and lymph node metastases from the same case, we included a tumor-draining lymph node (TdL) as a comparison specimen (Fig. 2a). The\u0026nbsp;TdL was a clinically diagnosed as metastasis-negative by rapid intraoperative evaluation.\u003c/p\u003e\n\u003cp\u003eWe applied Xenium spatial gene expression profiling using a custom 380-gene breast cancer panel (refer to\u0026nbsp;Supplementary Table 3\u0026nbsp;for the gene list). Across the primary tumor, metastatic lymph node and TdL, Xenium detected 1,005,436 cells\u0026nbsp;(Fig.\u0026nbsp;2b). Using gene expression signatures \u003csup\u003e37-39\u003c/sup\u003e and histological location information, we annotated clusters by major cell types. The luminal cell cluster was subdivided into two clusters, (i) and (ii), with (i) further branching into four subclusters: L1: luminal 1; L2: luminal 2; B: basal; and M: mesenchymal like (Fig.\u0026nbsp;2c).\u003c/p\u003e\n\u003cp\u003eWe identified genes significantly upregulated within each cluster. \u003cem\u003eKRT8\u003c/em\u003e and \u003cem\u003eKRT18\u003c/em\u003e, markers of mammary ductal epithelial cells, were enriched in L1 and L2, whereas \u003cem\u003eKRT5\u003c/em\u003e and \u003cem\u003eKRT14\u003c/em\u003e, myoepithelial markers, were enriched in cluster B (Supplementary Fig. 2). Consistent with expression information, spatial localization of L1, L2, and B clusters matched the lobular architecture of the mammary gland (Fig. 2d). High \u003cem\u003eANKRD30A\u003c/em\u003e expression in L2 and elevated \u003cem\u003eKIT\u003c/em\u003e and \u003cem\u003eKRT23\u003c/em\u003e expression in L1 aligned with two previously reported mature luminal cell types \u003csup\u003e40\u003c/sup\u003e: mature luminal cell 2 (\u003cem\u003eKRT18\u003c/em\u003e+/\u003cem\u003eANKRD30A\u003c/em\u003e+) and mature luminal cell 1 (\u003cem\u003eKRT18\u003c/em\u003e+/\u003cem\u003eKRT23\u003c/em\u003e+/\u003cem\u003eKIT\u003c/em\u003e+). L1, L2, and B clustered on the lobules, whereas M cells were distributed sporadically throughout the tumor, suggesting reduced cell adhesion (Fig. 2e). Consistent with this, the M cluster showed increased expression of EMT markers, including \u003cem\u003eTHY1\u003c/em\u003e, an EMT-related gene shown in Fig. 1e. Collectively, these findings support a role for \u003cem\u003eTHY1\u003c/em\u003e expression in EMT initiation from normal lobule cells (L1, L2, and B), potentially regulated through DNA methylation.\u003c/p\u003e\n\u003cp\u003eAlthough cluster (i) encompassed lobular and early EMT-stage cells, cluster (ii) localized pathologically to tumor regions and comprised definitive cancer cells. Strong \u003cem\u003eANKRD30A\u003c/em\u003e expression in cluster (ii) (Fig. 2f) suggested that L2 may be the tumor\u0026rsquo;s cell of origin \u003csup\u003e40\u003c/sup\u003e. Uniform manifold approximation and projection (UMAP) showed that metastatic lymph node tumor cells formed distinct clusters. Based on high \u003cem\u003eANKRD30A\u003c/em\u003e expression, we defined L2 as the tumor-originating cell and applied pseudotemporal trajectory analysis using Monocle 3 \u003csup\u003e41\u003c/sup\u003e. The inferred trajectory showed branching from L2 into (ii)-1, progressing to (ii)-2, and finally differentiating into metastatic tumor cells in lymph nodes (Fig. 2g). Taken together, these findings suggest a progression pathway in which cancer cells originate from L2, acquire EMT potential in the M cluster, and ultimately establish metastases within cluster (ii) in lymph nodes.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eSpatial analysis of the lymph nodes reveals six transcriptionally distinct EMT states\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eWhen inspecting spatial gene expression data to further investigate the mechanisms underlying metastasis, we unexpectedly identified a small cluster of only 30 cells within the TdL, a site clinically diagnosed as metastasis-free; these cells spanned a region approximately 200 \u0026mu;m in diameter (Fig.\u0026nbsp;3a).\u0026nbsp;Gene expression analysis revealed clear \u003cem\u003eKRT19\u003c/em\u003e expression, an epithelial marker absent in normal lymph node cells. Additional hematoxylin and eosin (HE) staining confirmed their morphology was consistent with cancer cells. The identification of these cells in the TdL provided a rare opportunity to examine molecular events at an extremely early stage of metastasis. This form of early dissemination is clinically referred to as isolated tumor cells (ITCs), a term we adopt hereafter.\u003c/p\u003e\n\u003cp\u003eWhen we visualized ITCs on the UMAP plane, we observed transcriptional heterogeneity (Fig.\u0026nbsp;3b, left panel), indicating that ITCs comprise several distinct cell types. To better characterize the biology of ITCs and the roles of different cell states in metastasis, we applied nonhierarchical clustering to subclassify tumor cells into six groups (C1\u0026ndash;C6) and assigned each ITC to one of these clusters (Fig.\u0026nbsp;3b, right panel). Owing to the small number of ITCs detected, subsequent analyses focused on the clusters to which each ITC mapped. For each subpopulation associated with an ITC, we identified distinct signaling pathways linking histology, spatial location, and gene expression (Fig.\u0026nbsp;3c\u0026ndash;e\u0026nbsp;and\u0026nbsp;Supplementary Table 4).\u003c/p\u003e\n\u003cp\u003eC6 occupied a position on the UMAP closest to EMT initiation (corresponding to the M cluster in Fig.\u0026nbsp;2c). Its scattered distribution within the primary tumor suggested low intercellular adhesion. C6 cells expressed high levels of ETM-related markers, such as \u003cem\u003eVIM\u003c/em\u003e and \u003cem\u003eCDH3\u003c/em\u003e, as well as stemness markers characteristic of breast cancer, including CD24\u003csup\u003elow\u003c/sup\u003e and CD44\u003csup\u003ehigh\u0026nbsp;\u003c/sup\u003e(Fig. 3d).\u0026nbsp;In contrast, epithelial marker expression was low, indicating a sustained mesenchymal phenotype. Gene enrichment analysis showed activation of EMT, stem cell signaling, TGF-\u0026beta; signaling, WNT signaling, immune checkpoint pathways, and antianoikis signaling (Fig.\u0026nbsp;3e). Signals from the microenvironment, such as TGF-\u0026beta; and Wnt ligands, induce various EMT-related transcription factors \u003csup\u003e42-43\u003c/sup\u003e. EMT also promotes \u003cem\u003ePD-L1\u003c/em\u003e expression, activating immune checkpoint pathways \u003csup\u003e44-45\u003c/sup\u003e. These findings suggest that C6 represents tumor cells in an extreme mesenchymal phenotype, capable of evading anoikis and immune clearance. In contrast, C4 most strongly retained epithelial features and was the most predominant cell population within lymph node metastases, forming the tumor mass\u0026nbsp;(Fig.\u0026nbsp;3c). C4 cells showed enrichment of cholesterol homeostasis and fatty acid metabolism pathways, with upregulation of fatty acid metabolism\u0026ndash;related molecules in metastatic lymph nodes (Fig.\u0026nbsp;3e, f). A previous mouse-based lymph node metastasis study demonstrated that a metabolic shift toward FAO in tumor cells is required for lymph node colonization by cancer cells \u003csup\u003e8\u003c/sup\u003e. Thus, C4 represents a cancer cell population that has acquired an extreme epithelial phenotype and metabolic reprogramming suited for metastatic growth.\u003c/p\u003e\n\u003cp\u003eClusters C1, C2, C3, and C5 represented hybrid epithelial\u0026ndash;mesenchymal (EM) cell populations that coexpressed epithelial and mesenchymal markers\u0026nbsp;(Fig. 3d).\u0026nbsp;The degree of coexpression differs among clusters, each cluster exhibited distinct features. C5 showed the strongest enrichment for the G2M checkpoint pathway (Fig. 3g). Prior studies have indicated that breast cancers with G2M activity display higher proliferative activity, increased MYC pathway activation, earlier metastasis, and worse survival \u003csup\u003e46\u003c/sup\u003e, suggesting that C5 may represent a more aggressive hybrid EMT subpopulation.\u0026nbsp;C3 was enriched for VEGF signaling and glycolysis, likely reflecting a hypoxic environment in the tumor core that promotes angiogenesis and glycolysis \u003csup\u003e47\u003c/sup\u003e.\u0026nbsp;C2, localized at the tumor periphery, i.e., the invasive front, displayed activation of\u0026nbsp;matrix metalloproteinases (MMPs) and canonical Wnt signaling.\u0026nbsp;Cells at\u0026nbsp;tumor margins exhibit EM plasticity and migratory behavior \u003csup\u003e48\u003c/sup\u003e, with MMP signaling \u003csup\u003e42\u003c/sup\u003e and the Wnt/\u0026beta; catenin pathway \u003csup\u003e49\u003c/sup\u003e known to promote invasion.\u0026nbsp;Canonical Wnt signaling was also enriched in C2 and hybrid EMT subpopulations (C3 and C6), consistent with reports that Wnt7A/B maintains the hybrid EM state \u003csup\u003e12,48\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003eC1 showed activation of retinoid metabolism and lipid transport pathways. Retinoic acid activates retinoic acid receptors and retinoid X receptors, driving expression of fatty acid metabolism genes \u003csup\u003e50\u003c/sup\u003e. Considering that C1 gene expression aligns closely with C4 (the most significantly proliferating cluster in lymph node metastases, involving elevated fatty acid metabolism), we speculate that C1 may act as a precursor population to C4.\u003c/p\u003e\n\u003cp\u003eCollectively, these findings suggest that hybrid EM cells undergo functional adaptation within primary tumor, with C6 as a possible EMT origin, C2 as an intermediate, C3 and C5 representing divergent aggressive states, and C1 transitioning toward the epithelialized, proliferative C4. The important thing is that each hybrid EMT subtype has different functional capabilities and that these cells are present in primary tumour from the very early stages of metastasis (Fig. 3h).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAnalysis of the tumor microenvironment in metastasized lymph node\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe major bottleneck for cancer cells infiltrating a lymph node is the new tumor microenvironment (TME) \u003csup\u003e51\u003c/sup\u003e. To address this, we compared interactions between tumor cells and their surrounding TME in the TdL and the metastatic lymph node (Fig. 4a). Notably, we detected distinct interaction patterns specific to each TME (Fig. 4b). One interaction detected uniquely in the TdL was the CD45\u0026ndash;MRC1 ligand\u0026ndash;receptor interaction, occurring between CD45 on T cells and MRC1 on macrophages (Fig. 4c). MRC1, an endocytosis receptor belonging to the C-type lectin family, is expressed on cells like dendritic cells, macrophages, and endothelial cells \u003csup\u003e52\u003c/sup\u003e. It has been shown to inhibit CD45 activity on T cells via direct interaction, leading to upregulation of the immune checkpoint protein CTLA4 and induction of antigen-specific T cell tolerance \u003csup\u003e53-55\u003c/sup\u003e. Therefore, immune escape may begin through this mechanism in the TdL.\u003c/p\u003e\n\u003cp\u003eMetastatic lymph nodes uniquely exhibited coinhibitory ligand\u0026ndash;receptor interactions, such as VEGFA\u0026ndash;KDR and CD86\u0026ndash;CTLA4, consistent with the promotion of angiogenesis and the establishment of a T cell\u0026ndash;exhausted environment already formed at this stage (Fig. 4d). Among interactions involving collagen, cancer-associated fibroblasts (CAFs) exhibited the highest number of interactions (Supplementary Fig. 3). Representative interactions in metastatic lymph nodes included COL1A1\u0026ndash;CD44 and COL4A1\u0026ndash;CD44, both known to promote tumor proliferation and progression \u003csup\u003e56-57\u003c/sup\u003e. In the chemokine signaling category, the CXCL16\u0026ndash;CXCR6 interaction was particularly prominent in metastatic lymph nodes (Supplementary Fig. 3). High CXCL16 expression is associated with histological malignancy and is implicated in tumor progression and metastasis via activation of the CXCL16\u0026ndash;CXCR6 axis \u003csup\u003e58\u003c/sup\u003e. These TME-derived signals appear to constitute necessary conditions for successful metastasis.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTranscriptional profiling of hybrid EMT and its association with patient outcomes\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003ePrompted by the preceding analyses, we sought to determine whether the transcriptional profiles we identified are associated with patient prognosis in a broader clinical context. To this end, we reanalyzed a previously reported scRNA-seq dataset from metastatic lymph nodes of patients with breast cancer \u003csup\u003e4\u003c/sup\u003e (Fig.\u0026nbsp;5a). Further clustering of cancer cells from the dataset revealed five clusters (sc8,sc4,sc2,sc15 and sc10) (Fig.\u0026nbsp;5b).\u0026nbsp;According to publicly available information, sc2,sc4 and sc8 were metastatic lymph node-derived samples, and sc15, sc10 was derived from the primary tumor. A metabolic shift toward fatty acids were observed in metastatic lymph node-derived samples(sc2,sc4 and sc8)(Fig. 5c). These observations are consistent with our results shown in Figure 3 and support our findings that a metabolic shift toward fatty acid metabolism is necessary for the formation of metastatic colonies in lymph nodes.\u0026nbsp;To investigate whether molecular subtype signatures derived from our Xenium data are associated with clinical outcomes, we evaluated their prognostic relevance using the METABRIC dataset \u003csup\u003e59\u003c/sup\u003e. Interestingly, the subtypes associated with poor prognosis were not clusters that formed metastatic colonies in lymph node (e.g., C4) or those showing strong mesenchymal features (C6) but rather clusters C5 and C3, which exhibited signatures of G2M cell cycle progression and increased glycolytic activity (Supplementary Fig. 4).\u003c/p\u003e\n\u003cp\u003eTo validate this finding, we mapped the hybrid EM clusters identified in our Xenium analysis (C1\u0026ndash;C6) to the corresponding cancer cell subclusters (sc2, sc4, sc8, sc10, and sc15) in the public scRNA-seq dataset. This revealed a one-to-one correspondence between C5 and sc8 and between C6 and sc15 (Fig. 5d). Consistent with our findings, only sc8 (corresponding to the C5 signature in our Xenium data) was significantly associated with poor prognosis [log-rank P = 4e-09, hazard ratio (HR) = 1.536, 95% confidence interval (CI) = 1.330\u0026ndash;1.772] (Fig. 5e). In contrast, sc15, which displayed the highest mesenchymal programming and corresponded to the C6 signature in our Xenium data, was associated with a more favorable prognosis (log-rank P = 0.01, HR = 0.836, 95% CI = 0.727\u0026ndash;0.961).\u003c/p\u003e"},{"header":"Discussion","content":"\u003cp\u003eIn this study, we tracked the entire metastatic cascade of breast cancer at single-cell resolution and identified the relationship between spatial evolution of tumour cells \u0026nbsp;and dynamic heterogeneity within the EMT spectrum. Our capture of early-stage metastatic events, including exceptionally early dissemination to lymph nodes, provides a rare in vivo snapshot of EMT dynamics during initial colonization. Thus, we were able to uncover distinct subpopulations with specific biological roles and their potential clinical relevance. To our knowledge, this is the first report to detail such fine-grained EMT profiles in vivo during early metastasis.\u003c/p\u003e\n\u003cp\u003eThrough detailed analysis of pioneering metastatic tumour cells within metastatic cascade, we identified six distinct types of hybrid epithelial–mesenchymal (EM) cells that play key roles in the early stages of metastasis. The subpopulation that underwent a lipid metabolic shift successfully formed colonies in lymph nodes. This is consistent with previous studies indicating that FAO supports lymph node metastasis in breast cancer \u003csup\u003e8\u003c/sup\u003e . Interestingly, the subpopulation responsible for forming lymph node colonies that underwent a lipid metabolic shift, as well as the subpopulation with the most mesenchymal features, was not associated with poor clinical prognosis. This supports previous evidence that mesenchymal-like phenotypes may correlate with favorable outcomes in patients with breast cancer \u003csup\u003e60\u003c/sup\u003e.\u0026nbsp;Instead,\u0026nbsp;poor prognosis was linked to transcriptional programs associated with MYC and E2F signaling, along with increased aerobic glycolysis and G2/M cell cycle activity.\u0026nbsp;In particular,\u0026nbsp;MYC reprograms metabolism toward enhanced glycolytic flux \u003csup\u003e61\u003c/sup\u003e, and breast cancers with enriched G2M signaling have been associated with elevated MYC pathway activity, increased proliferation, early distant metastasis, and worse survival \u003csup\u003e46\u003c/sup\u003e.\u0026nbsp;Collectively, these findings suggest that distinct EMT subpopulations drive different aspects of tumor progression: some mediate niche colonization, whereas others may drive aggressive metastasis and recurrence, implying the need for tailored therapeutic strategies \u003csup\u003e11\u003c/sup\u003e. A key question concerns how hybrid EM cells disseminate in lymph nodes: do metastases arise via the sequential model, where fully EMT-induced M cells colonize distant sites by MET, or via the cooperative model, which mesenchymal cells support Epithelial cells, which serve as metastasis-initiating cells \u003csup\u003e62\u003c/sup\u003e? Both models are considered plausible based on prior studies. In the present study, the observation of multiple EMT states in ITCs at an extremely early metastatic stage (Fig. 3a), along with similar mixtures of hybrid EM cells in the primary tumor (Fig. 3h), supports the cooperative metastasis model, where different cell states collaborate to initiate metastasis.\u003c/p\u003e\n\u003cp\u003eWe also examined the early interactions between metastatic tumor cells and the lymph node microenvironment.\u0026nbsp;Comparisons between the TdL and metastatic lymph nodes showed that a coordinated, immunosuppressive environment forms around tumor cells in newly colonized niches.\u0026nbsp;These environments are likely shaped, in part, by the metabolic phenotypes of the tumor cells themselves. For example,\u0026nbsp;lactic acid, a byproduct of glycolysis, has been shown to promote alternative splicing in T cells and enhance CTLA-4 expression in a Foxp3-dependent manner \u003csup\u003e63\u003c/sup\u003e.\u0026nbsp;The plasticity of EMT observed in this study may contribute to the ecosystem that supports metastasis outgrowth, including immunomodulatory effects \u003csup\u003e64\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003eDespite the novel insights in this study, it has several limitations.\u0026nbsp;First, our in-depth analysis was conducted on tissue from a single patient, limiting the generalizability of our conclusions.\u0026nbsp;Second, owing to the limited number of cells available, we were unable to perform molecular or biochemical validations.\u0026nbsp;However, by integrating our findings with prior knowledge and publicly available single-cell RNA-seq datasets, we were able to infer the likely roles of key subpopulations and validate their clinical relevance.\u0026nbsp;Finally, we could not definitively distinguish whether the observed metastasis occurred via hematogenous or lymphatic dissemination.\u0026nbsp;Future studies should also aim to elucidate whether distinct molecular features exist between these metastatic pathways.\u003c/p\u003e\n\u003cp\u003eIn summary, our findings highlight therapeutic vulnerabilities embedded within the EMT spectrum, particularly involving metabolic plasticity and niche adaptation. These insights offer promising avenues for precision diagnostics and the development of targeted therapies aimed at halting metastasis at its earliest and most tractable stage.\u003c/p\u003e"},{"header":"Materials and Methods","content":"\u003cp\u003e\u003cstrong\u003eClinical sample\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThis study was approved by the Clinical Ethics Committee of St. Marianna University School of Medicine (approval number: 2297-i103). Informed consent was obtained from the patient. The surgical specimen was from an 80-year-old female patient, who received surgery for breast cancer (BRC-26). The patient had not received neoadjuvant or adjuvant chemotherapy (detailed clinicopathologic findings are provided in\u0026nbsp;Supplementary Table 1). The same specimen was used in our previous study \u003csup\u003e65\u003c/sup\u003e. For multiregional, multi-omics analyses, fresh frozen samples of the primary tumor and adjacent normal tissue were used. Two nontumor regions (NC1 and NC2) and two tumor regions (Ca3, a solid tumor region, and Ca4, a sparse tumor region), as well as normal breast tissue from the same specimen (used as the control), were microdissected using the AVENIO MilliSect system (Roche, Pleasanton, CA) and subjected to multi-omics analysis, including RNA-seq, WGS, and EM-seq. For spatial analyses via Xenium, formalin-fixed paraffin-embedded (FFPE) samples of the primary tumor, lymph node metastases, and TdL were used.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eRNA-seq\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eRNA-seq was conducted as previously described \u003csup\u003e65\u003c/sup\u003e. Briefly, total RNA was extracted from frozen normal breast tissue using the RNeasy Micro Kit (Qiagen). RNA-seq libraries were prepared using the SMART-Seq Stranded Kit (Takara Bio). Paired-end 150-bp sequencing was performed on the NovaSeq 6000 system (Illumina).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eRNA-seq data analysis\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eRNA-seq data from the normal tissue obtained in this study and four tumor regions from our previous study \u003csup\u003e65\u003c/sup\u003e were analyzed. Adapter trimming was performed using fastp (v0.23.2) \u003csup\u003e66\u003c/sup\u003e. Reads mapping to rRNA were removed using Bowtie 2 (v2.3.4.3) \u003csup\u003e67\u003c/sup\u003e. The retained reads were aligned to the human reference genome GRCh38.p12 using STAR (v2.7.5c) \u003csup\u003e68\u003c/sup\u003e. Gene-level read counts were obtained using featureCounts (v2.0.2) \u003csup\u003e69\u003c/sup\u003e, and RPKM values were calculated. DEGs between histological classes were identified using DESeq2 (v1.42.1) \u003csup\u003e70\u003c/sup\u003e applying the Wald test with an adjusted P-value cutoff of \u0026lt;0.1. Gene enrichment analysis was performed using Metascape \u003csup\u003e71\u003c/sup\u003e .\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ePreparation of EM-seq and WGS libraries\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eEM-seq and WGS library preparations were conducted as described previously \u003csup\u003e72\u003c/sup\u003e. Briefly, gDNA was extracted from normal tissue using NucleoSpin Tissue XS (MACHEREY-NAGEL), and 100 ng of this gDNA was fragmented using the M220 Focused-ultrasonicator (Covaris). Adapter ligation was conducted using the NEBNext Enzymatic Methyl-seq Kit (New England BioLabs). Half of the adapter-ligated DNA was used for WGS library preparation (five cycles PCR amplification), and the other half was employed for EM-seq library preparation, which included TET oxidation, APOBEC conversion, and PCR amplification (six cycles). Libraries were sequenced as paired-end 150-bp reads on the NovaSeq 6000 system (Illumina).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eEM-seq data analysis\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eEM-seq data from the normal tissue obtained in this study and four tumor regions from our previous study were analyzed following established protocols \u003csup\u003e65\u003c/sup\u003e. Briefly, adapter trimming was conducted using Trim Galore (v0.6.4_dev; https://github.com/FelixKrueger/TrimGalore). Trimmed reads were aligned to the human reference genome using Bismark (v0.22.1) \u003csup\u003e73\u003c/sup\u003e. Duplicate reads were removed using deduplicate_bismark, and methylation information in a CpG context was extracted via bismark_methylation_extractor. Genome-wide methylation profiles, including patterns in CpG islands, CpG shores, and promoter regions, were visualized as described previously \u003csup\u003e74\u003c/sup\u003e. Differentially methylated regions (DMRs) among NC-1, NC-2, Ca-3, and Ca-4 and the broader NC and Ca regions were identified using metilene v0.2-8 \u003csup\u003e75\u003c/sup\u003e, as per our previous study \u003csup\u003e65\u003c/sup\u003e. Promotor methylation was analyzed among NC-1, NC-2, Ca-3, and Ca-4 spanning 1 kb upstream to 500 bp downstream of transcription start sites, focusing on DMRs with methylation differences \u0026gt;10%. Visualisation was performed using IGV \u003csup\u003e76\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ePoint mutation detection\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eWGS reads were aligned to the human reference genome using BWA-MEM (v0.7.17) \u003csup\u003e77\u003c/sup\u003e with default settings. Duplicate reads were marked using Picard MarkDuplicates (v2.23.8) (https://broadinstitute.github.io/picard/). Somatic mutations were called using GATK Mutect2 and filtered with FilterMutectCalls (v4.1.3.0) \u003csup\u003e78\u003c/sup\u003e. Variants were annotated via ANNOVAR \u003csup\u003e79\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCNV detection\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eCNVs in NC-1, NC-2, Ca-3, Ca-4, and normal samples were detected using FACETS (v0.6.2), as described previously \u003csup\u003e74\u003c/sup\u003e. Copy number gains and losses were defined as ≥4 copies and ≤1 copy, respectively.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eSpatial heatmaps of gene expression levels and average methylation rates\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eSpatial heatmaps of gene expression levels and average methylation rates were generated using a custom workflow implemented on Nikon’s pilot analysis platform. Specifically, microscopic images (Olympus BX53, 2× objective) with substantial overlap were aligned using similarity transformation based on scale-invariant feature transform–detected feature point pairs. This alignment produced a composite image of the entire specimen with a common coordinate system. Manually dissected region contours were aligned to this coordinate system. Regions of interest were color-coded by gene expression levels or average methylation rates. The average methylation rate was calculated as the mean frequency of methylation at CpG sites covered by at least five reads within the predefined promoter regions.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003e\u003cem\u003eIn situ\u003c/em\u003e\u003c/strong\u003e\u003cstrong\u003e\u0026nbsp;gene expression analysis via Xenium\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eSpatial subclonal analysis at single-cell resolution was performed using the Xenium Slides \u0026amp; Sample Prep Reagents (10× Genomics), a predesigned human breast panel, and a custom panel, as previously described \u003csup\u003e80\u003c/sup\u003e. Briefly, FFPE tissue sections (5 µm thickness) of the primary tumor, lymph node metastases, and TdL were mounted onto Xenium slides. Following deparaffinization and decrosslinking, probe hybridization, probe circularization, and rolling circle amplification were performed. Detection of amplified probes was conducted using the Xenium Analyzer (10× Genomics). In total, 380 target genes, comprising 280 from the predesigned panel and 100 from the custom panel, are listed in Supplementary Table 1.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eComputational processing of Xenium \u003cem\u003ein situ\u003c/em\u003e expression data and analysis\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eRaw output files from the Xenium Analyzer were processed using Seurat (v5.0.2) \u003csup\u003e81\u003c/sup\u003e. Data normalization was performed using the SCTransform method. Clustering analysis and UMAP visualization were conducted using the first 30 principal components, with the former performed via the FindNeighbors and FindClustersing functions. The same parameters were applied to subclustering analyses (Figs. 2c and 3). For each cluster relevant to downstream analysis, DEGs were identified using the FindMarkers() function in Seurat with the default parameters (only.pos = TRUE, min.pct = 0.25, logfc.threshold = 0.25). Cluster annotation was based on previously published gene expression signatures \u003csup\u003e37-39\u003c/sup\u003e and spatial histological context.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eTrajectory analysis\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTrajectory and pseudotime analysis of luminal cells was performed using Monocle 3 (v1.3.4) \u003csup\u003e41\u003c/sup\u003e. Preprocessing was performed using the preprocess_cds function, with the number of dimensions set to 100. Dimensionality reduction and clustering were implemented using the reduce_dimension and cluster_cells functions, respectively. Principal graph learning was conducted with the learn_graph function. Based on spatial and Xenium data, mature luminal cell 2 was designated as the trajectory root. Pseudotime ordering was visualized using the plot_cells function.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eIdentifying enriched gene signatures in Xenium clusters\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo investigate the functional characteristics of Xenium clusters (C1–C6 in Fig. 3e), pathway and gene set enrichment analyses were performed using Metascape \u003csup\u003e71\u003c/sup\u003e based on the top 50 DEGs per cluster (Supplementary Table 5). Enrichment was assessed across GO categories, including Biological Process, Cellular Component, and Molecular Function, as well as other biological pathways.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAnalysis of intercellular communication networks\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eIntercellular communication within metastatic lymph nodes and the TdL was analyzed using the CellChat package (v1.6.1) \u003csup\u003e82\u003c/sup\u003e, which quantitatively identifies cell–cell interactions and communications. Following the official pipeline (https://github.com/jinworks/CellChat), analyses were performed using default settings. Statistically significant ligand–receptor interactions (P \u0026lt; 0.05) were extracted. Visualization of cell–cell interactions was achieved using circle plots, bar plots, heatmaps, and violin plots.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003ePublished dataset processing\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eWe downloaded a publicly available single-cell RNA-seq dataset (GSE180286 from Guan et al.\u0026nbsp;\u003csup\u003e4\u003c/sup\u003e). Data integration and clustering were performed in Seurat (v5.0.2) \u003csup\u003e81\u003c/sup\u003e, focusing on HER2-positive breast tumors and two lymph node samples (GSM5457205, GSM5457206, and GSM5457207). Louvain clustering at a resolution of 0.1 yielded nine clusters: one T cell cluster, three B cell clusters, one NK cell cluster, one CAF1 cluster, and one macrophage cluster. Clusters of the same cell type were merged (Fig. 5a). Subclustering was performed on the cancer cell populations (Fig. 5b).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eClinical validation using METABRIC data\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eClinical validation was conducted using the METABRIC dataset \u003csup\u003e59\u003c/sup\u003e. METABRIC gene transcriptome data, as well as clinical and sample level metadata, were downloaded from cBioPortal (https://www.cbioportal.org). Analysis was performed on 1,980 breast cancer cases with available survival outcomes. For each cancer cluster defined in Fig. 5c, the top 100 DEGs (Supplementary Table 6) were used to score individual METABRIC samples. ROC curve analysis was used to define optimal cutoffs for stratifying patients into high and low score groups. Disease-free survival was analyzed using the survival (v3.5-7) and survminer (v0.4.9) R packages. HRs and 95% CIs were calculated using Cox proportional hazards models. Statistical significance was evaluated via the log-rank test.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eStatistical Analysis\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eStatistical methods and tests are detailed in the figure legends. All analyses were performed using R (v4.3.2) or Python.\u003c/p\u003e"},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eData availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe Xenium and bulk multi-omics data supporting the findings of this study are available at the DDBJ Japanese Genotype-phenotype Archive (https://gr-sharingdbs.dbcls.jp) under accession number JGAD000946. The single-cell RNA-seq dataset from Guan et al. \u003csup\u003e4\u003c/sup\u003e was downloaded from Gene Expression Omnibus (GSE180286). METABRIC transcriptome and clinical data were downloaded from cBioPortal (https://www.cbioportal.org). Additional data are available from the corresponding authors upon reasonable request.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAcknowledgements:\u0026nbsp;\u003c/strong\u003eWe thank K. Imamura, K. Abe, M. Satake, J. Zenkou, E. Sekimori, R. Fujinaga, and A. Gouda for their technical assistance. The authors would like to thank Enago (www.enago.jp) for the English language review.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eFunding:\u0026nbsp;\u003c/strong\u003eThis work was supported by the Japan Agency for Medical Research and Development (AMED) P-PROMOTE, grant number JP21ck0106700 (to S.N.), the Japan Society for the Promotion of Science Grants-in-Aid for Scientific Research (KAKENHI Grants) Number JP24K11738 (Platform for Advanced Genome Science), as well as by Nikon Corporation. The supercomputing resource was provided by Human Genome Center, the University of Tokyo (http://sc.hgc.jp/shirokane.html). Computational resources were also supplied by Kashiwa-no-ha Omics Gate (https://www.kog.or.jp/en/server.html)\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAuthor contributions:\u0026nbsp;\u003c/strong\u003eS.N. and K.K collected the human samples, performed experiments and computational analyses, and generated all figures and tables.\u003c/p\u003e\n\u003cp\u003eA.M., T.I., T.O., and K.T. aided in human sample collection.\u003c/p\u003e\n\u003cp\u003eT.O. generated Fig. 2a.\u003c/p\u003e\n\u003cp\u003eE.I. performed computational analyses.\u003c/p\u003e\n\u003cp\u003eM.K. and M.N. assisted in the Xenium experiments.\u003c/p\u003e\n\u003cp\u003eA.H., I.M. and J.K. conducted the pathological review.\u003c/p\u003e\n\u003cp\u003eA.M. and K.T. developed viewer software.\u003c/p\u003e\n\u003cp\u003eA.K. and A.S. offered advice and reviewed the manuscript.\u003c/p\u003e\n\u003cp\u003eS.N., M.S., and Y.S. designed the project and wrote the manuscript.\u003c/p\u003e\n\u003cp\u003eAll authors read and approved the final manuscript.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCompeting interests:\u0026nbsp;\u003c/strong\u003eA.M. and K.T. are employees of Nikon Corporation.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\n\u003cli\u003eVeronesi, U., et al. A randomized comparison of sentinel-node biopsy with routine axillary dissection in breast cancer. \u003cem\u003eN Engl J Med\u003c/em\u003e. \u003cstrong\u003e349\u003c/strong\u003e, 546-553 (2003).\u003c/li\u003e\n\u003cli\u003eTurner, R. R., Ollila, D. W., Krasne, D. L. \u0026amp; Giuliano, A. E. Histopathologic validation of the sentinel lymph node hypothesis for breast carcinoma. \u003cem\u003eAnn Surg\u003c/em\u003e. \u003cstrong\u003e226\u003c/strong\u003e, 271-276; discussion 276-278 (1997).\u003c/li\u003e\n\u003cli\u003eNagasawa, S., et al. Genomic profiling reveals heterogeneous populations of ductal carcinoma in situ of the breast. \u003cem\u003eCommun Biol\u003c/em\u003e. \u003cstrong\u003e4\u003c/strong\u003e, 438 (2021).\u003c/li\u003e\n\u003cli\u003eXu, K., et al. Microenvironment components and spatially resolved single-cell transcriptome atlas of breast cancer metastatic axillary lymph nodes. \u003cem\u003eActa Biochim Biophys Sin (Shanghai)\u003c/em\u003e. \u003cstrong\u003e54\u003c/strong\u003e, 1336-1348 (2022).\u003c/li\u003e\n\u003cli\u003ePesce, C. \u0026amp; Morrow, M. The need for lymph node dissection in nonmetastatic breast cancer. \u003cem\u003eAnnu Rev Med\u003c/em\u003e. \u003cstrong\u003e64\u003c/strong\u003e, 119-129 (2013).\u003c/li\u003e\n\u003cli\u003ede Boer, M., van Dijck, J. A., Bult, P., Borm, G. F. \u0026amp; Tjan-Heijnen, V. C. Breast cancer prognosis and occult lymph node metastases, isolated tumor cells, and micrometastases. \u003cem\u003eJ Natl Cancer Inst\u003c/em\u003e. \u003cstrong\u003e102\u003c/strong\u003e, 410-425 (2010).\u003c/li\u003e\n\u003cli\u003eHinshaw, D. C. \u0026amp; Shevde, L. A. The Tumor Microenvironment Innately Modulates Cancer Progression. \u003cem\u003eCancer Res\u003c/em\u003e. \u003cstrong\u003e79\u003c/strong\u003e, 4557-4566 (2019).\u003c/li\u003e\n\u003cli\u003eLee, C. K., et al. Tumor metastasis to lymph nodes requires YAP-dependent metabolic adaptation. \u003cem\u003eScience\u003c/em\u003e. \u003cstrong\u003e363\u003c/strong\u003e, 644-649 (2019).\u003c/li\u003e\n\u003cli\u003eCohen, M. R. Check concentrations of medications ordered for children. \u003cem\u003eNursing\u003c/em\u003e. \u003cstrong\u003e16\u003c/strong\u003e, 8 (1986).\u003c/li\u003e\n\u003cli\u003eNieto, M. A., Huang, R. Y., Jackson, R. A. \u0026amp; Thiery, J. P. Emt: 2016. \u003cem\u003eCell\u003c/em\u003e. \u003cstrong\u003e166\u003c/strong\u003e, 21-45 (2016).\u003c/li\u003e\n\u003cli\u003ePuisieux, A., Brabletz, T. \u0026amp; Caramel, J. Oncogenic roles of EMT-inducing transcription factors. \u003cem\u003eNat Cell Biol\u003c/em\u003e. \u003cstrong\u003e16\u003c/strong\u003e, 488-494 (2014).\u003c/li\u003e\n\u003cli\u003eRubtsova, S. N., Zhitnyak, I. Y. \u0026amp; Gloushankova, N. A. Dual role of E-cadherin in cancer cells. \u003cem\u003eTissue Barriers\u003c/em\u003e. \u003cstrong\u003e10\u003c/strong\u003e, 2005420 (2022).\u003c/li\u003e\n\u003cli\u003eFischer, K. R., et al. Epithelial-to-mesenchymal transition is not required for lung metastasis but contributes to chemoresistance. \u003cem\u003eNature\u003c/em\u003e. \u003cstrong\u003e527\u003c/strong\u003e, 472-476 (2015).\u003c/li\u003e\n\u003cli\u003eYe, X., et al. Upholding a role for EMT in breast cancer metastasis. \u003cem\u003eNature\u003c/em\u003e. \u003cstrong\u003e547\u003c/strong\u003e, E1-E3 (2017).\u003c/li\u003e\n\u003cli\u003eDiepenbruck, M. \u0026amp; Christofori, G. Epithelial-mesenchymal transition (EMT) and metastasis: yes, no, maybe? \u003cem\u003eCurr Opin Cell Biol\u003c/em\u003e. \u003cstrong\u003e43\u003c/strong\u003e, 7-13 (2016).\u003c/li\u003e\n\u003cli\u003eAiello, N. M., et al. Upholding a role for EMT in pancreatic cancer metastasis. \u003cem\u003eNature\u003c/em\u003e. \u003cstrong\u003e547\u003c/strong\u003e, E7-E8 (2017).\u003c/li\u003e\n\u003cli\u003ePastushenko, I., et al. Identification of the tumour transition states occurring during EMT. \u003cem\u003eNature\u003c/em\u003e. \u003cstrong\u003e556\u003c/strong\u003e, 463-468 (2018).\u003c/li\u003e\n\u003cli\u003eKroger, C., et al. Acquisition of a hybrid E/M state is essential for tumorigenicity of basal breast cancer cells. \u003cem\u003eProc Natl Acad Sci U S A\u003c/em\u003e. \u003cstrong\u003e116\u003c/strong\u003e, 7353-7362 (2019).\u003c/li\u003e\n\u003cli\u003eTsai, J. H., Donaher, J. L., Murphy, D. A., Chau, S. \u0026amp; Yang, J. Spatiotemporal regulation of epithelial-mesenchymal transition is essential for squamous cell carcinoma metastasis. \u003cem\u003eCancer Cell\u003c/em\u003e. \u003cstrong\u003e22\u003c/strong\u003e, 725-736 (2012).\u003c/li\u003e\n\u003cli\u003eHaerinck, J., Goossens, S. \u0026amp; Berx, G. The epithelial-mesenchymal plasticity landscape: principles of design and mechanisms of regulation. \u003cem\u003eNat Rev Genet\u003c/em\u003e. \u003cstrong\u003e24\u003c/strong\u003e, 590-609 (2023).\u003c/li\u003e\n\u003cli\u003eBian, X., et al. Carboxylic acid metabolism in cancer: Mechanisms, microenvironment interactions, and therapeutic opportunities. \u003cem\u003eMetabolism\u003c/em\u003e. \u003cstrong\u003e171\u003c/strong\u003e, 156334 (2025).\u003c/li\u003e\n\u003cli\u003eVymetalkova, V., Vodicka, P., Vodenkova, S., Alonso, S. \u0026amp; Schneider-Stock, R. DNA methylation and chromatin modifiers in colorectal cancer. \u003cem\u003eMol Aspects Med\u003c/em\u003e. \u003cstrong\u003e69\u003c/strong\u003e, 73-92 (2019).\u003c/li\u003e\n\u003cli\u003eEhrlich, M. DNA hypomethylation in cancer cells. \u003cem\u003eEpigenomics\u003c/em\u003e. \u003cstrong\u003e1\u003c/strong\u003e, 239-259 (2009).\u003c/li\u003e\n\u003cli\u003eRodriguez, J., et al. Chromosomal instability correlates with genome-wide DNA demethylation in human primary colorectal cancers. \u003cem\u003eCancer Res\u003c/em\u003e. \u003cstrong\u003e66\u003c/strong\u003e, 8462-9468 (2006).\u003c/li\u003e\n\u003cli\u003eUsman, S., et al. Vimentin Is at the Heart of Epithelial Mesenchymal Transition (EMT) Mediated Metastasis. \u003cem\u003eCancers (Basel)\u003c/em\u003e. \u003cstrong\u003e13\u003c/strong\u003e (2021).\u003c/li\u003e\n\u003cli\u003eWang, Y., Shi, J., Chai, K., Ying, X. \u0026amp; Zhou, B. P. The Role of Snail in EMT and Tumorigenesis. \u003cem\u003eCurr Cancer Drug Targets\u003c/em\u003e. \u003cstrong\u003e13\u003c/strong\u003e, 963-972 (2013).\u003c/li\u003e\n\u003cli\u003eDaSilva-Arnold, S. C., et al. ZEB2, a master regulator of the epithelial-mesenchymal transition, mediates trophoblast differentiation. \u003cem\u003eMol Hum Reprod\u003c/em\u003e. \u003cstrong\u003e25\u003c/strong\u003e, 61-75 (2019).\u003c/li\u003e\n\u003cli\u003eSekimata, M., Murakami-Sekimata, A. \u0026amp; Homma, Y. CpG methylation prevents YY1-mediated transcriptional activation of the vimentin promoter. \u003cem\u003eBiochem Biophys Res Commun\u003c/em\u003e. \u003cstrong\u003e414\u003c/strong\u003e, 767-772 (2011).\u003c/li\u003e\n\u003cli\u003eChen, Y., Wang, K., Qian, C. N. \u0026amp; Leach, R. DNA methylation is associated with transcription of Snail and Slug genes. \u003cem\u003eBiochem Biophys Res Commun\u003c/em\u003e. \u003cstrong\u003e430\u003c/strong\u003e, 1083-1090 (2013).\u003c/li\u003e\n\u003cli\u003eXie, H., et al. Significance of ZEB2 in the immune microenvironment of colon cancer. \u003cem\u003eFront Genet\u003c/em\u003e. \u003cstrong\u003e13\u003c/strong\u003e, 995333 (2022).\u003c/li\u003e\n\u003cli\u003eValdivia, A., Avalos, A. M. \u0026amp; Leyton, L. Thy-1 (CD90)-regulated cell adhesion and migration of mesenchymal cells: insights into adhesomes, mechanical forces, and signaling pathways. \u003cem\u003eFront Cell Dev Biol\u003c/em\u003e. \u003cstrong\u003e11\u003c/strong\u003e, 1221306 (2023).\u003c/li\u003e\n\u003cli\u003eBrenet, M., et al. Thy-1 (CD90)-Induced Metastatic Cancer Cell Migration and Invasion Are beta3 Integrin-Dependent and Involve a Ca(2+)/P2X7 Receptor Signaling Axis. \u003cem\u003eFront Cell Dev Biol\u003c/em\u003e. \u003cstrong\u003e8\u003c/strong\u003e, 592442 (2020).\u003c/li\u003e\n\u003cli\u003evan Zijl, F., et al. Hepatic tumor-stroma crosstalk guides epithelial to mesenchymal transition at the tumor edge. \u003cem\u003eOncogene\u003c/em\u003e. \u003cstrong\u003e28\u003c/strong\u003e, 4022-4033 (2009).\u003c/li\u003e\n\u003cli\u003eSethi, S., Macoska, J., Chen, W. \u0026amp; Sarkar, F. H. Molecular signature of epithelial-mesenchymal transition (EMT) in human prostate cancer bone metastasis. \u003cem\u003eAm J Transl Res\u003c/em\u003e. \u003cstrong\u003e3\u003c/strong\u003e, 90-99 (2010).\u003c/li\u003e\n\u003cli\u003eTinsley, E., Bredin, P., Toomey, S., Hennessy, B. T. \u0026amp; Furney, S. J. KMT2C and KMT2D aberrations in breast cancer. \u003cem\u003eTrends Cancer\u003c/em\u003e. \u003cstrong\u003e10\u003c/strong\u003e, 519-530 (2024).\u003c/li\u003e\n\u003cli\u003eNishimura, T., et al. Evolutionary histories of breast cancer and related clones. \u003cem\u003eNature\u003c/em\u003e. \u003cstrong\u003e620\u003c/strong\u003e, 607-614 (2023).\u003c/li\u003e\n\u003cli\u003eAran, D., Hu, Z. \u0026amp; Butte, A. J. xCell: digitally portraying the tissue cellular heterogeneity landscape. \u003cem\u003eGenome Biol\u003c/em\u003e. \u003cstrong\u003e18\u003c/strong\u003e, 220 (2017).\u003c/li\u003e\n\u003cli\u003eLim, E., et al. Aberrant luminal progenitors as the candidate target population for basal tumor development in BRCA1 mutation carriers. \u003cem\u003eNat Med\u003c/em\u003e. \u003cstrong\u003e15\u003c/strong\u003e, 907-913 (2009).\u003c/li\u003e\n\u003cli\u003ePliner, H. A., Shendure, J. \u0026amp; Trapnell, C. Supervised classification enables rapid annotation of cell atlases. \u003cem\u003eNat Methods\u003c/em\u003e. \u003cstrong\u003e16\u003c/strong\u003e, 983-986 (2019).\u003c/li\u003e\n\u003cli\u003eNguyen, Q. H., et al. Profiling human breast epithelial cells using single cell RNA sequencing identifies cell diversity. \u003cem\u003eNat Commun\u003c/em\u003e. \u003cstrong\u003e9\u003c/strong\u003e, 2028 (2018).\u003c/li\u003e\n\u003cli\u003eQiu, X., et al. Reversed graph embedding resolves complex single-cell trajectories. \u003cem\u003eNat Methods\u003c/em\u003e. \u003cstrong\u003e14\u003c/strong\u003e, 979-982 (2017).\u003c/li\u003e\n\u003cli\u003eLamouille, S., Xu, J. \u0026amp; Derynck, R. Molecular mechanisms of epithelial-mesenchymal transition. \u003cem\u003eNat Rev Mol Cell Biol\u003c/em\u003e. \u003cstrong\u003e15\u003c/strong\u003e, 178-196 (2014).\u003c/li\u003e\n\u003cli\u003eZeng, Y. A. \u0026amp; Nusse, R. Wnt proteins are self-renewal factors for mammary stem cells and promote their long-term expansion in culture. \u003cem\u003eCell Stem Cell\u003c/em\u003e. \u003cstrong\u003e6\u003c/strong\u003e, 568-577 (2010).\u003c/li\u003e\n\u003cli\u003eZhang, Y. C., Zhang, Y. T., Wang, Y., Zhao, Y. \u0026amp; He, L. J. What role does PDL1 play in EMT changes in tumors and fibrosis? \u003cem\u003eFront Immunol\u003c/em\u003e. \u003cstrong\u003e14\u003c/strong\u003e, 1226038 (2023).\u003c/li\u003e\n\u003cli\u003eImodoye, S. O. \u0026amp; Adedokun, K. A. EMT-induced immune evasion: connecting the dots from mechanisms to therapy. \u003cem\u003eClin Exp Med\u003c/em\u003e. \u003cstrong\u003e23\u003c/strong\u003e, 4265-4287 (2023).\u003c/li\u003e\n\u003cli\u003eOshi, M., et al. G2M Cell Cycle Pathway Score as a Prognostic Biomarker of Metastasis in Estrogen Receptor (ER)-Positive Breast Cancer. \u003cem\u003eInt J Mol Sci\u003c/em\u003e. \u003cstrong\u003e21\u003c/strong\u003e (2020).\u003c/li\u003e\n\u003cli\u003eGraziano, F., et al. Glycolysis gene expression analysis and selective metabolic advantage in the clinical progression of colorectal cancer. \u003cem\u003ePharmacogenomics J\u003c/em\u003e. \u003cstrong\u003e17\u003c/strong\u003e, 258-264 (2017).\u003c/li\u003e\n\u003cli\u003eSacco, J. L. \u0026amp; Gomez, E. W. Epithelial-Mesenchymal Plasticity and Epigenetic Heterogeneity in Cancer. \u003cem\u003eCancers (Basel)\u003c/em\u003e. \u003cstrong\u003e16\u003c/strong\u003e (2024).\u003c/li\u003e\n\u003cli\u003eTeeuwssen, M. \u0026amp; Fodde, R. Cell Heterogeneity and Phenotypic Plasticity in Metastasis Formation: The Case of Colon Cancer. \u003cem\u003eCancers (Basel)\u003c/em\u003e. \u003cstrong\u003e11\u003c/strong\u003e (2019).\u003c/li\u003e\n\u003cli\u003eKeller, H., et al. Fatty acids and retinoids control lipid metabolism through activation of peroxisome proliferator-activated receptor-retinoid X receptor heterodimers. \u003cem\u003eProc Natl Acad Sci U S A\u003c/em\u003e. \u003cstrong\u003e90\u003c/strong\u003e, 2160-2164 (1993).\u003c/li\u003e\n\u003cli\u003eGalimberti, V., et al. Axillary dissection versus no axillary dissection in patients with sentinel-node micrometastases (IBCSG 23-01): a phase 3 randomised controlled trial. \u003cem\u003eLancet Oncol\u003c/em\u003e. \u003cstrong\u003e14\u003c/strong\u003e, 297-305 (2013).\u003c/li\u003e\n\u003cli\u003eMartinez-Pomares, L. The mannose receptor. \u003cem\u003eJ Leukoc Biol\u003c/em\u003e. \u003cstrong\u003e92\u003c/strong\u003e, 1177-1186 (2012).\u003c/li\u003e\n\u003cli\u003eSchuette, V., et al. Mannose receptor induces T-cell tolerance via inhibition of CD45 and up-regulation of CTLA-4. \u003cem\u003eProc Natl Acad Sci U S A\u003c/em\u003e. \u003cstrong\u003e113\u003c/strong\u003e, 10649-10654 (2016).\u003c/li\u003e\n\u003cli\u003eLuca, M. E., et al. Mannosylated PLP(139-151) induces peptide-specific tolerance to experimental autoimmune encephalomyelitis. \u003cem\u003eJ Neuroimmunol\u003c/em\u003e. \u003cstrong\u003e160\u003c/strong\u003e, 178-187 (2005).\u003c/li\u003e\n\u003cli\u003eKel, J., et al. Soluble mannosylated myelin peptide inhibits the encephalitogenicity of autoreactive T cells during experimental autoimmune encephalomyelitis. \u003cem\u003eAm J Pathol\u003c/em\u003e. \u003cstrong\u003e170\u003c/strong\u003e, 272-280 (2007).\u003c/li\u003e\n\u003cli\u003eChoi, J. H., et al. Single-cell transcriptome profiling of the stepwise progression of head and neck cancer. \u003cem\u003eNat Commun\u003c/em\u003e. \u003cstrong\u003e14\u003c/strong\u003e, 1055 (2023).\u003c/li\u003e\n\u003cli\u003eWang, T., et al. COL4A1 promotes the growth and metastasis of hepatocellular carcinoma cells by activating FAK-Src signaling. \u003cem\u003eJ Exp Clin Cancer Res\u003c/em\u003e. \u003cstrong\u003e39\u003c/strong\u003e, 148 (2020).\u003c/li\u003e\n\u003cli\u003eHan, J., et al. CXCL16 Promotes Gastric Cancer Tumorigenesis via ADAM10-Dependent CXCL16/CXCR6 Axis and Activates Akt and MAPK Signaling Pathways. \u003cem\u003eInt J Biol Sci\u003c/em\u003e. \u003cstrong\u003e17\u003c/strong\u003e, 2841-2852 (2021).\u003c/li\u003e\n\u003cli\u003eCurtis, C., et al. The genomic and transcriptomic architecture of 2,000 breast tumours reveals novel subgroups. \u003cem\u003eNature\u003c/em\u003e. \u003cstrong\u003e486\u003c/strong\u003e, 346-352 (2012).\u003c/li\u003e\n\u003cli\u003eGrosse-Wilde, A., et al. Loss of inter-cellular cooperation by complete epithelial-mesenchymal transition supports favorable outcomes in basal breast cancer patients. \u003cem\u003eOncotarget\u003c/em\u003e. \u003cstrong\u003e9\u003c/strong\u003e, 20018-20033 (2018).\u003c/li\u003e\n\u003cli\u003eDong, Y., Tu, R., Liu, H. \u0026amp; Qing, G. Regulation of cancer cell metabolism: oncogenic MYC in the driver\u0026apos;s seat. \u003cem\u003eSignal Transduct Target Ther\u003c/em\u003e. \u003cstrong\u003e5\u003c/strong\u003e, 124 (2020).\u003c/li\u003e\n\u003cli\u003eTsuji, T., Ibaragi, S. \u0026amp; Hu, G. F. Epithelial-mesenchymal transition and cell cooperativity in metastasis. \u003cem\u003eCancer Res\u003c/em\u003e. \u003cstrong\u003e69\u003c/strong\u003e, 7135-7139 (2009).\u003c/li\u003e\n\u003cli\u003eDing, R., et al. Lactate modulates RNA splicing to promote CTLA-4 expression in tumor-infiltrating regulatory T cells. \u003cem\u003eImmunity\u003c/em\u003e. \u003cstrong\u003e57\u003c/strong\u003e, 528-540 e526 (2024).\u003c/li\u003e\n\u003cli\u003eReticker-Flynn, N. E., et al. Lymph node colonization induces tumor-immune tolerance to promote distant metastasis. \u003cem\u003eCell\u003c/em\u003e. \u003cstrong\u003e185\u003c/strong\u003e, 1924-1942 e1923 (2022).\u003c/li\u003e\n\u003cli\u003eKeisuke Kunigo, et al. Deep long-read methylation analysis using targeted nanopore sequencing and application to local region analysis of clinical specimen, 2025).\u003c/li\u003e\n\u003cli\u003eChen, S., Zhou, Y., Chen, Y. \u0026amp; Gu, J. fastp: an ultra-fast all-in-one FASTQ preprocessor. \u003cem\u003eBioinformatics\u003c/em\u003e. \u003cstrong\u003e34\u003c/strong\u003e, i884-i890 (2018).\u003c/li\u003e\n\u003cli\u003eLangmead, B. \u0026amp; Salzberg, S. L. Fast gapped-read alignment with Bowtie 2. \u003cem\u003eNat Methods\u003c/em\u003e. \u003cstrong\u003e9\u003c/strong\u003e, 357-359 (2012).\u003c/li\u003e\n\u003cli\u003eDobin, A., et al. STAR: ultrafast universal RNA-seq aligner. \u003cem\u003eBioinformatics\u003c/em\u003e. \u003cstrong\u003e29\u003c/strong\u003e, 15-21 (2013).\u003c/li\u003e\n\u003cli\u003eLiao, Y., Smyth, G. K. \u0026amp; Shi, W. featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. \u003cem\u003eBioinformatics\u003c/em\u003e. \u003cstrong\u003e30\u003c/strong\u003e, 923-930 (2014).\u003c/li\u003e\n\u003cli\u003eLove, M. I., Huber, W. \u0026amp; Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. \u003cem\u003eGenome Biol\u003c/em\u003e. \u003cstrong\u003e15\u003c/strong\u003e, 550 (2014).\u003c/li\u003e\n\u003cli\u003eZhou, Y., et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. \u003cem\u003eNat Commun\u003c/em\u003e. \u003cstrong\u003e10\u003c/strong\u003e, 1523 (2019).\u003c/li\u003e\n\u003cli\u003eSakamoto, Y., et al. Long-read whole-genome methylation patterning using enzymatic base conversion and nanopore sequencing. \u003cem\u003eNucleic Acids Res\u003c/em\u003e. \u003cstrong\u003e49\u003c/strong\u003e, e81 (2021).\u003c/li\u003e\n\u003cli\u003eKrueger, F. \u0026amp; Andrews, S. R. Bismark: a flexible aligner and methylation caller for Bisulfite-Seq applications. \u003cem\u003eBioinformatics\u003c/em\u003e. \u003cstrong\u003e27\u003c/strong\u003e, 1571-1572 (2011).\u003c/li\u003e\n\u003cli\u003eHaga, Y., et al. Whole-genome sequencing reveals the molecular implications of the stepwise progression of lung adenocarcinoma. \u003cem\u003eNat Commun\u003c/em\u003e. \u003cstrong\u003e14\u003c/strong\u003e, 8375 (2023).\u003c/li\u003e\n\u003cli\u003eJuhling, F., et al. metilene: fast and sensitive calling of differentially methylated regions from bisulfite sequencing data. \u003cem\u003eGenome Res\u003c/em\u003e. \u003cstrong\u003e26\u003c/strong\u003e, 256-262 (2016).\u003c/li\u003e\n\u003cli\u003eThorvaldsd\u0026oacute;ttir, H., Robinson, J. T. \u0026amp; Mesirov, J. P. Integrative Genomics Viewer (IGV): high-performance genomics data visualization and exploration. \u003cem\u003eBrief Bioinform\u003c/em\u003e. \u003cstrong\u003e14\u003c/strong\u003e, 178-192 (2013).\u003c/li\u003e\n\u003cli\u003eLi, H. \u0026amp; Durbin, R. Fast and accurate long-read alignment with Burrows-Wheeler transform. \u003cem\u003eBioinformatics\u003c/em\u003e. \u003cstrong\u003e26\u003c/strong\u003e, 589-595 (2010).\u003c/li\u003e\n\u003cli\u003eMcKenna, A., et al. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. \u003cem\u003eGenome Res\u003c/em\u003e. \u003cstrong\u003e20\u003c/strong\u003e, 1297-1303 (2010).\u003c/li\u003e\n\u003cli\u003eWang, K., Li, M. \u0026amp; Hakonarson, H. ANNOVAR: functional annotation of genetic variants from high-throughput sequencing data. \u003cem\u003eNucleic Acids Res\u003c/em\u003e. \u003cstrong\u003e38\u003c/strong\u003e, e164 (2010).\u003c/li\u003e\n\u003cli\u003eTakano, Y., et al. Spatially resolved gene expression profiling of tumor microenvironment reveals key steps of lung adenocarcinoma development. \u003cem\u003eNat Commun\u003c/em\u003e. \u003cstrong\u003e15\u003c/strong\u003e, 10637 (2024).\u003c/li\u003e\n\u003cli\u003eHao, Y., et al. Dictionary learning for integrative, multimodal and scalable single-cell analysis. \u003cem\u003eNat Biotechnol\u003c/em\u003e. \u003cstrong\u003e42\u003c/strong\u003e, 293-304 (2024).\u003c/li\u003e\n\u003cli\u003eJin, S., et al. Inference and analysis of cell-cell communication using CellChat. \u003cem\u003eNat Commun\u003c/em\u003e. \u003cstrong\u003e12\u003c/strong\u003e, 1088 (2021).\u003c/li\u003e\n\u003c/ol\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":"npj-breast-cancer","isNatureJournal":false,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"npjbcancer","sideBox":"Learn more about [npj Breast Cancer](http://www.nature.com/npjbcancer/)","snPcode":"41523","submissionUrl":"https://mts-npjbcancer.nature.com/","title":"npj Breast Cancer","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"NPJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-7371816/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-7371816/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"Lymph node metastasis critically influences breast cancer prognosis; however, the cellular mechanisms underlying the earliest metastatic events remain unclear. In spatial transcriptomic analysis of a patient with breast cancer at single cell resolution, we identified 30 tumor cells representing the initial metastatic seeding in a lymph node. These cells originated from multiple epithelial–mesenchymal (EM) transition status and included six distinct subpopulations with biological significance. Only cells exhibiting metabolic shift toward fatty acid metabolism successfully established lymph node colonies, implicating this shift in metastatic fitness. The tumor microenvironment surrounding these cells showed immunosuppressive and tumor-promoting features, supporting metastasis establishment. Cross-referencing these expression profiles with public datasets revealed that poor prognosis correlated not with fully mesenchymal or metastatic populations, but with hybrid EM cells exhibiting epithelial and mesenchymal traits. These findings highlight the metabolic and phenotypic plasticity of metastatic cells and serve as translational bridges between the spatial evolution of tumor cells in the extremely early stages of lymph nodes metastasis and clinical prognosis in breast cancer.","manuscriptTitle":"Spatial Gene Expression Analysis Reveals Drivers of Extremely Early Lymph Node Metastasis in Breast Cancer","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-09-23 07:37:56","doi":"10.21203/rs.3.rs-7371816/v1","editorialEvents":[{"type":"communityComments","content":0},{"type":"decision","content":"Revision requested","date":"2025-11-22T17:00:08+00:00","index":"","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2025-11-11T16:05:25+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2025-10-03T09:04:01+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"97490610253786960053934344771052714141","date":"2025-09-19T14:27:26+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"121239200669755630732956844372700376141","date":"2025-09-15T05:24:35+00:00","index":"hide","fulltext":""},{"type":"reviewersInvited","content":"","date":"2025-09-14T19:39:36+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2025-08-20T20:18:40+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2025-08-20T13:33:36+00:00","index":"","fulltext":""},{"type":"submitted","content":"npj Breast Cancer","date":"2025-08-14T08:45:52+00:00","index":"","fulltext":""}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"npj-breast-cancer","isNatureJournal":false,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"npjbcancer","sideBox":"Learn more about [npj Breast Cancer](http://www.nature.com/npjbcancer/)","snPcode":"41523","submissionUrl":"https://mts-npjbcancer.nature.com/","title":"npj Breast Cancer","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"NPJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"89314b98-bef9-43fa-bba5-33619b97a4f6","owner":[],"postedDate":"September 23rd, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"published-in-journal","subjectAreas":[{"id":55021349,"name":"Biological sciences/Cancer"},{"id":55021350,"name":"Biological sciences/Computational biology and bioinformatics"},{"id":55021351,"name":"Health sciences/Oncology"}],"tags":[],"updatedAt":"2026-01-26T16:12:52+00:00","versionOfRecord":{"articleIdentity":"rs-7371816","link":"https://doi.org/10.1038/s41523-026-00897-1","journal":{"identity":"npj-breast-cancer","isVorOnly":false,"title":"npj Breast Cancer"},"publishedOn":"2026-01-23 15:56:58","publishedOnDateReadable":"January 23rd, 2026"},"versionCreatedAt":"2025-09-23 07:37:56","video":"","vorDoi":"10.1038/s41523-026-00897-1","vorDoiUrl":"https://doi.org/10.1038/s41523-026-00897-1","workflowStages":[]},"version":"v1","identity":"rs-7371816","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-7371816","identity":"rs-7371816","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}

Text is read by the "Ask this paper" AI Q&A widget below. Extraction quality varies by source — PMC NXML preserves structure cleanly, OA-HTML may include some navigation residue, and OA-PDF can have broken hyphenation. The publisher copy (via DOI) is the canonical version.

My notes (saved in your browser only)

Ask this paper AI returns verbatim quotes from the full text · source: preprint-html

Answers must be backed by verbatim quotes from this paper's full text. Hallucinated quotes are dropped automatically; if no verbatim passage answers the question, we say so. How this works

Citation neighborhood (no data yet)

We don't have any in-corpus citations linked to this paper yet. This is a recent paper (2025) — citers typically take a year or two to land, and the OpenAlex reference graph may still be filling in.

Source provenance

europepmc
last seen: 2026-05-20T01:45:00.602351+00:00
unpaywall
last seen: 2026-05-22T02:00:06.705733+00:00
License: CC-BY-4.0