Profiles of nasopharyngeal carcinoma cells in 10 distinct clusters identified through scRNA-seq analysis

To generate a deep transcriptional atlas of nasopharyngeal carcinomas (NPCs), high-resolution 10 × genomics scRNA-seq was performed on data from 24 primary, 7 peripheral blood mononuclear cell (PBMC), and 7 metastatic nasopharyngeal carcinoma (Fig. 1A). Following high quality control and filtration, a total of 292,298 cells from these tumors were classified into 10 clusters: T cells (177,664 cells), B cells (56,496 cells), Macrophages/Monocytes (14,506 cells), Natural Killer (NK) cells (27,534 cells), Plasma cells (6,347 cells), plasmacytoid Dendritic Cells (pDC, 1,251 cells), Migratory Dendritic Cells (805 cells), Mast cells (855 cells), Cancer-Associated Fibroblasts (CAFs, 807 cells), and Epithelial cells (6,033 cells) (Fig. 1B–C).

Fig. 1

Dissection of the tumor microenvironment in Nasopharyngeal Carcinoma (NPC) with scRNA-seq. A Workflow diagram showing the samples from 24 primary NPC tumors, 7 PBMC and 7 metastatic tumors for scRNA-seq. B UMAP plots of cells from the 38 samples profiled in this study, with each cell color-coded to indicate its associated cell types. C Using Harmony to integrate three sets of scRNA data, primary tumors (blue), PBMC (green), and metastatic tumors (red). D A dotsplot illustrates marker genes exhibiting high expression across 10 clusters, high expression (Red), low expression (blue). E, F FeaturePlot and ViolinPlot depict the distribution and quantification of specific marker genes within each subpopulation, respectively

We used batch correction and data integration to remove batch effects across datasets. Moreover, each group was identified according to marker genes. T cells exhibited significantly higher CD3D, CD3E, CD8A, and GZMK. B cells were responsible for producing antibodies and participating in the humoral immune response, such as MS4A1, CD79A, CD79B, and CD19. Macrophages and monocytes demonstrated high levels of LYZ, CD14, S100A8, C1QA, CD163, and MRC1. Natural Killer (NK) cells expressed GNLY, NKG7, KLRD1, and PRF1. Plasma cells showed increased expression of MZB1, IGHG1, IGHA1, and TNFRSF17. Plasmacytoid dendritic cells (pDCs) were defined by IL3RA, LILRA4, and CLEC4C. Migratory dendritic cells associated with ability to migrate from peripheral tissues to secondary lymphoid organs, such as CCR7, LAMP3, EBI3, and PDCD1LG2. Mast cells were characterized by the presence of TPSAB1, CPA3, MS4A2, and TPSB2. Cancer-associated fibroblasts (CAFs) were identified by COL1A1, COL3A1, DCN, and FAP. Epithelial cells were characterized by the expression of specific markers, such as CD24, KRT19, EPCAM, and KRT18. (Fig. 1D).

We also found that PTPRC was highly expressed in immune cell-related fractions, indicating its importance in these groups. CD3D was present in the majority of T cells, whereas CD8A selectively identified a subset, highlighting the diversity of T cell functions by distinguishing CD4 + and CD8+ cells. CD19 gene distinguished B cell subtypes, and myeloid cells were detected by LYZ, revealing their immunological involvement. Due to their exclusive relationship with GNLY, NK cells were crucial to innate immunity. Plasma, pDCs, Migratory DCs, and Mast cells were specifically defined by IGHG1, IL3RA, LAMP3, and TPSAB1, highlighting their diverse and specialized functions. COL1A1 was mostly identified in fibroblasts, while EPCAM was only found in epithelial cells, demonstrating their vital involvement in tissue integrity. This marker gene distribution mapping illuminates the complex interaction and functional distinctions that distinguish healthy and pathological tissue states (Fig. 1E) and enhances our understanding of cellular architecture. Violin plots show similar patterns and distributions. (Fig. 1F).

Each groups exhibited diverse functional potentiality in different cancer types

We examined nasopharyngeal cancer primary, metastatic, and peripheral blood mononuclear cell heterogeneity in three datasets. The PBMC group had significantly more B and NK cells than the primary and metastatic groups. There were no significant differences in T cell populations. Epithelial cells were observed solely in primary and metastatic samples, with metastatic samples exhibiting a higher proportion. This study illustrates the complicated immunological and epithelial cell patterns of nasopharyngeal cancer phases and components (Fig. 2A).

Fig. 2

Comparison of single-cell data from three groups. A The proportion of cells from three groups (primary NPC, PBMC, and metastatic NPC) among 10 populations. B A scatter plot displaying differential expression across multiple groups, with the top 10 upregulated genes marked for each group. The selection standard include |avg_log2FC|> = 1 and p_val_adj < 0.05. Upregulated genes are depicted in red, while downregulated genes are shown in green. C A violin plot show the top 10 genes with the most significant differences among primary NPC, PBMC, and metastatic NPC. D, E The dotplot illustrating the gene ontology and reactome pathway enrichment analysis conducted on the top 50 differentially expressed genes across the three groups, respectively

Volcano plot showed up-regulated and down-regulated genes across dataset, including the top10 differential genes in each group (Fig. 2B), with these same genes highlighted in the violin diagram (Fig. 2C). The upregulation of HSPA1A and HSPA1B in primary tumor tissues emphasizes their role in protein folding and cellular stress responses. EGR1 and CDKN1A (p21), which regulate cell proliferation, apoptosis, and cycle regulation, are also elevated. NR4A1, DUSP4 and GZMK show cellular differentiation, survival, and immune response, and ITM2C, CD83 and FABP5 show signaling and metabolic pathways. In PBMCs, high expression of ADGRG1, CX3CR1, S1PR5, and KLRD1 indicates substantial immune system interaction with the tumor, whereas SPON2, FGR, FAM65B, and FCGR3A indicate active immune cell signaling and response SNHG series expression increases in metastatic tumors, affecting RNA modification and cell metastatic pathways. GAS5 may reduce cancer, while SEPTIN6, PCED1B-AS1 may help cells migrate and adapt to different microenvironments. Across tumor phases, a complex regulatory network controls tumor formation, immune system, and cellular activity (Fig. 2C).

The functional enrichment analysis of the 50 most variably expressed genes reveals diverse biological activities across primary, PBMC, and metastatic tumor samples. Primary tumor immune activation requires T cell responses and apoptotic signals, therefore viral interactions may help them hide. These malignancies address protein stress and cell adhesion, suggesting ERK1/2 for proliferation, differentiation, invasiveness, and migration. Adhesion molecules and cytokines help myeloid and lymphocytes defend pathogens and detect tumors in PBMCs. Cancer-specific cytoskeletal alterations and autophagy help metastatic tumors move and adapt to stress. Survival, proliferation in new habitats, and immune evasion are promoted by these adaptations. Additionally, TGF-β pathways boost cell differentiation and metastatic properties, showing the complex link between tumor evolution and host immune responses (Fig. 2D-E).

Variation of T cells sub-populations across primary and metastatic tumors

T cells from nasopharyngeal carcinoma were isolated, re-clustered, and divided into 10 populations (Fig. 3A). All T cells in these subpopulations express CD3E. Two CD8+ groups, CD8_T-GZMK and CD8_T-GNLY, express GZMB, a cytotoxic effector marker, suggesting they may destroy tumor cells. CD8_T-GNLY have high level GZMK, GNLY, and PRF1 expression boosted their effector T cell status. However, CD8_T-NELL2, expressing NELL2 and CCR7, resembles naïve T cells, indicating an unstimulated CD8+ T cell cluster. Gene expressions show CD4+ T cell variety. CD4_T-SPDYA and CD4_T-FHIT subsets, expressing CCR7, may constitute naïve CD4+ T cells. Their functions match T cell phenotypes, with TXNIP and TOX2 involvement suggesting metabolic and transcriptional regulation. FOXP3 and IL2RA-high Tregs regulate tumor immune evasion and immunological homeostasis. Treg-RBMS3 elevated expression suggests a regulatory role, possibly necessary for tumor-mediated immunosuppression. In response to the immune response, pro-T cells with high MKI67 and TOP2A expression may proliferate rapidly (Fig. 3B, C).

Fig. 3

T cells undergo dynamic re-clustering within NPC microenvironment. A UMAP plot displaying 10 distinct clusters of T cells, including CD4 T cells, CD8 T cells, Treg cells, and proliferating T cells, each uniquely represented by different colors. B The dotsplot illuminates the expression patterns of marker genes across 10 distinct T cell clusters. C The FeaturePlot visualizes the distribution of specific gene expression within each group, effectively illustrating how characteristic genes are expressed across different cell populations. D A violin plot visually highlights the variation in gene expression levels across groups within T cells, focusing on the top 10 genes that exhibit the most significant differences. E RNA velocity delineates the fate trajectories of T cells within primary and metastatic tumors, with purple denoting the early stages of differentiation and yellow signifying the terminal stages. F The heatmap displays the differential transcriptional regulation among T cell subpopulations. G The heatmap exhibits genes of primary and metastatic tumors as they vary along the differentiation trajectory

An analysis of primary and metastatic cancers revealed distinct biological variations in the top 10 genes that were expressed differently. The immune response and regulation of the microenvironment were primarily influenced by genes expressed in the original tumor. Tumors exhibit the expression of immunoglobulin genes IGHA1 and IGHG1, which suggests the presence of immunological responses mediated by antibodies. Elevated concentrations of CCL4, CCL4L2, CCL3, and CCL3L1 suggest that the tumor is attracting T cells and macrophages to carry out immunological surveillance. The activity of the MAPK signaling pathway is crucial for the growth and survival of cancer cells, and it may be inhibited by DUSP4. The increased expression of GZMB indicates that cytotoxic T lymphocytes specifically attack tumor cells, whereas the transcriptional regulation of BHLHE40 may modify the development of tumor cells and their biological cycles. Protein stability and signal transmission are regulated by PTMS and ubiquitination. Genes that are highly expressed in metastatic tumors enhance the migration, invasion, and ability to adapt to stress in the surrounding microenvironment. BACH2, a regulatory protein that controls gene expression, has the potential to modify the process of cell specialization and the body's defense mechanisms. Elevated levels of LEF1 and PLAC8 expression indicate the presence of tumor cells that are likely to adhere and migrate, leading to metastasis. MT1X, MT1E, and MT2A can aid metastatic cells in managing oxidative stress in unfamiliar settings by sequestering and detoxifying metal ions. Metastatic tumors employ SESN3 and HSPH1 to adjust to harsh circumstances, hence controlling protein quality and folding (Fig. 3D).

The transition of pro-T cells to effector cells in primary tumors of nasopharyngeal cancer involved the maturation of undifferentiated T cells into functional T cells (Fig. 3E). The trajectory is demonstrated through signal transmission, immunological stimulation, inflammatory reactions, and the activation of metabolic genes during the process of cell differentiation. The genes ACTG1 and ITGAE govern cell adhesion and cytoskeleton. PCLAF, IL2RA, and CTLA4 are crucial immune response genes that work in conjunction with FOXP3 and IFI16. NR4A2 promotes immune homeostasis and prevents autoimmune diseases by enhancing immunological tolerance and Tregs. ALOX5AP and AOAH regulate lipid metabolism and inflammation, whereas PPP2R2B and TBC1D4 govern transcription and cell differentiation (Fig. 3G). Treg cells undergo a conversion process from FOXP3 to RBMS3, whereas metastatic cancers undergo differentiation from CD4_T_SPDYA to CD4_T_TOX2. The differentiation process starts with the CD4_T_TXNIP molecule (Fig. 3E). This may facilitate tumor cells in adapting to new surroundings by boosting their ability to survive, multiply, and spread. The genes UBA52 and IL7R are involved in cellular metabolism and proliferation, and tumor-immune interactions, respectively. The genes YPEL5, NR4A2, and ZNF331 have a role in signal transduction and stress response, whereas NKG7 and MT1E are involved in immunological control. These genes showcase the intricate behavior of tumor cells. The genes ATP1B3 and ARHGEF39 influence the movement of cells, whereas CD58 and ICOS enhance the ability of cells to evade the immune system, suggesting strategies for metastasis. The study of RNF125 in protein ubiquitination provides insights into how tumor cells adapt to maintain protein balance across different organs (Fig. 3G). Metastatic tumor cells undergo differentiation by evolving their gene expression, which emphasizes metabolic pathways for growth, modulation of stress response, and the intricate immunological milieu.

Discrepancies are observed when comparing transcriptional regulation across in primary and metastatic tumor groupings. The increase of DDIT3 indicates the presence of metastatic stress, whereas the elevation of NFATC1 supports the control of immune response to promote tumor growth. These cancerous growths increase the expression of CCNT2 and CREB1, necessitating improved regulation of the cell cycle in new environments to facilitate rapid proliferation. The expression of HOXB8 enhances the process of tumor cell differentiation and proliferation. Immunological tolerance and evasion in in situ cancers are affected by elevated levels of immunomodulatory genes, such as FOXP3 and GATA3. HIVEP1 and PBX4 enhance cellular differentiation and maintain structural integrity, whereas CEBPD, FOSL2, KLF13, and MAFF have the potential to exacerbate inflammation and stimulate tumor growth. The findings illustrate the intricate and dynamic interaction of regulatory systems in tumors, encompassing processes such as development, immune evasion, and environmental adaptability (Fig. 3F).

Elucidating heterogeneity and fate differentiation of B cells in primary and metastatic nasopharyngeal carcinoma

We employed unsupervised clustering to classify B cells in nasopharyngeal cancer in order to get insight into their mechanisms. This method detected 8 distinct groups, with showing the presence of typical B cell markers such as MS4A1, CD79A, CD79B, and CD19. Naive B cells had higher levels expression of BACH2, TCL1A, and FCER2, suggesting the presence of immature B cells. The B_CLECL1 and memory B DNAH8 exhibited the expression of CLECL1, whereas the latter demonstrated an upregulation of DNAH8, FCRL4, and CCR1, suggesting the presence of a distinct subtype of memory B cells. The presence of notable IL7R and BCL11B expression in the B_IL7R subgroup indicates a more favorable response to the IL7 signal. Both the proliferative_germinal_center_B and germinal_center_B subsets had elevated levels of AICDA, RGS13, and GCSAM, which suggests that there is active germinal center B cell activity. The former also had elevated MKI67 gene, suggesting ongoing cellular proliferation. The Plasma_cells_FNDC3B exhibited substantial expression of FNDC3B, MZB1, IGHG1, IGHA1, and TNFRSF17, indicating the presence of plasma cell characteristics (Fig. 4A–C).

Fig. 4

B Cells manifest Dynamic Reorganization and Diversity within the Microenvironment of NPC. A UMAP plot reveals 8 clusters of B cells, including naïve B cells, Memory B cells, Germinal Center B cells, and Plasma Cells, each uniquely identified by varied colors. B The dot plot elucidates the expression profiles of marker genes within eight unique B cell clusters. C The FeaturePlot offers a visualization of the expression distribution of particular genes across groups, clearly depicting the manifestation of characteristic genes within diverse cell populations. D A violin plot graphically emphasizes the fluctuations in gene expression levels among groups within B cells, centering on the top 10 genes showcasing the most pronounced differences. E RNA velocity traces the progression trajectories of B cells across primary and metastatic tumors, with purple representing the onset of differentiation and yellow denoting the final stages. F The heatmap delineates the differential transcriptional landscapes among B cell subpopulations, highlighting the distinct regulatory mechanisms and gene expression profiles that define the heterogeneity within these cellular subsets. G The heatmap shows the changes in gene expression related to the trajectory across various B cell subgroups

Immunoglobulin genes such as IGHG2 and IGKC exhibit potent antibody-mediated protection by B cells in the early stages of cancer. This activity highlights the tactics of the B cell brigade in combating tumor invasion. TRIB1 and GAS6 collaborate to control cell viability and promote cellular communication, which is crucial for the proliferation and dissemination of tumors. During the metastatic phase, the tumor exhibits resilience to metastatic stress through the expression of heat shock proteins HSPH1 and HSPE1. This resilience enables the tumor to survive and adapt to new surroundings. The presence of metal ion processing genes MT2A and MT1X contributes to the narrative by demonstrating the requirement for equilibrium in order to prevent cellular demise caused by toxicity. The gene regulation and DNA repair mechanisms of ZNF10 contribute to the survival of tumors, whereas the immune evasion abilities of CD72 and FCER2 complicate the progression of the tumor. The emphasis of TCL1A on cell division and inhibition of programmed cell death, along with NIBAN3 in managing low oxygen levels, demonstrates the ability of tumor cells to adjust to challenging surroundings. (Fig. 4D).

The gene expression in primary tumors progresses from Naive B cells to B-CLEC1L, and finally transitions to Memory B cells. This process involves the control of genes such as WDR76, PCLAF, and CYTOR. The process demonstrates the development and functional change of B cells in the primary tumor microenvironment. The genes involved in this process are responsible such as repairing DNA, controlling transcription, regulating the cell cycle, and modulating the immune response. In contrast, the pattern of gene activity in metastatic cancers is more intricate, progressing from B-CLEC1L to Memory B cells, then through Proliferative germinal center B cells to Germinal center B cells, and finally culminating in the formation of Plasma cells. The key genes involved in this process are DERL3, PRDM1, CD247, and others. These genes primarily play a role in B cell proliferation, differentiation, immunological control, and cell death pathways. This alteration exposes B cells may undergo additional modifications during the metastatic process to enhance the survival and spread of tumor cells (Fig. 4E, G).

During the progression of primary tumors and their spread to other parts of the body, transcriptional regulatory factors play a crucial role by impacting different stages of tumor growth in a distinct manner. During the first stage, ELF3, a transcription factor that is exclusive to epithelial cells, has a vital function in maintaining cell identity while promoting cell growth and specialization. Simultaneously, SPI1 and BCL11B play crucial roles in coordinating the differentiation and activity of immune cells, so regulating the dynamic interaction between the tumor and the immunological defenses of the host. SREBF2 regulates lipid metabolism, providing tumor cells with vital energy and the necessary components for cellular membranes. Moreover, a combination of TGIF2, ZNF394, STAT4, TCF7, and NR2F1 contributes to the development of tumor by influencing several signaling pathways and gene expression patterns. On the other hand, metastatic tumors focus on improving DNA repair, speeding up the cell cycle, and effectively utilizing genetic information. The increased expression of BRCA1 serves as evidence of the improved repair mechanisms that are triggered in response to DNA damage. The E2F family members, such as E2F1, E2F2, E2F4, and E2F8, promote fast growth of cancer cells, ensuring their ability to adapt and withstand the challenging process of metastasis. YBX1 plays a critical role in enhancing the stability of mRNA and the efficiency of translation. This narrative provides a detailed description of the regulatory landscape of transcription, highlighting the significant impact these factors have on tumor growth at every stage (Fig. 4F).

Diversity of myeloid sub-populations in primary and metastatic tumors

In order to observe the effect of myeloid cells on the tumor microenvironment, they were divided into 9 groups (Fig. 5A). The Mac1, Mac2, Mono, and Neutrophil subgroups expressed the CD68, C1QA, and C1QB genes, while the LYZ gene was widely activated. The CD163 and MRC1 genes in Mac2 promote tissue healing and anti-inflammatory response. In addition, the Mono, which expresses more S100A12, VCAN, and FCN1, regulates the immune system and inflammation. Neutrophils expresses the S100A8 gene, indicating their synergy with dendritic cells. Moreover, cDC, and mDC activate T cells via CD1C and CLEC10A. FCER1A, LILRA4, and CLEC4C distinguish the mDC, which modulates the immune system. LILRA4 and CLEC4C are critical for virus defense and autoimmune balance. The migratory DC, with proliferation genes MKI67 and TOP2A and expression of CCR7, LAMP3, EBI3, and PDCD1LG2, shows how myeloid cells shape tumor progression and the immunological environment (Fig. 5B, C).

Fig. 5

Exploring the Diversity and Subgroup Dynamics of Myeloid Cells in NPC Across Primary and Metastatic Tumors. A UMAP display reveals 9 unique subsets of myeloid cells, encompassing macrophages, monocytes, neutrophils, dendritic cells, and proliferative cells, each distinguished by a unique color palette. B The dot plot offers a detailed view of the expression patterns for marker genes within the eight delineated myeloid cell clusters, shedding light on their distinct genetic identities. C The FeaturePlot presents the distribution of specific gene expressions among various groups, effectively illustrating the unique gene signatures across different myeloid cell populations. D A violin plot highlights the variance in gene expression across myeloid cell groups, focusing on the top 10 genes with the most notable differences. E RNA velocity charts the developmental paths of myeloid cells in both primary and metastatic tumors, with purple indicating the early stages of differentiation and yellow marking the advanced stages. F A heatmap outlines the diverse transcriptional regulation within myeloid cell subsets, emphasizing the unique regulatory mechanisms and gene expression patterns that contribute to their heterogeneity. G A heatmap captures the dynamic gene expression changes associated with the developmental trajectories among myeloid cell subgroups, showing the evolution of gene expression as myeloid cells progress through their differentiation pathways

Tumor biology is diverse and multifaceted, as shown by primary and metastatic tumor gene expression profiles. High gene expressions of IL6, CCL2, and SLC25A37 in primary tumors indicate chronic inflammation in the tumor microenvironment, which affects immune system modulation. IL6 and CCL2, key regulators, boost tumor immune evasion and promote a tumor-promoting milieu. Furthermore, PTGS2 expression is closely connected to pain, inflammation, and tumor growth. High levels of IGLC3, TREM1, and SPP1 indicate immune cell activation and tumor cell contacts. MT1E and MT1X upregulation increases oxidative stress resistance, while NR3C1 emphasizes hormone signaling in tumor survival and adaptation. SOX4, a transcription factor, may promote phenotypic alterations and increase invasiveness. Metastatic tumors maintain genomic integrity and cope with stress by upregulating POLB and ITM2C, which are involved in DNA repair and the endoplasmic reticulum stress response. The increased expression of PLD4 and AFF3 genes suggests specialized cellular signaling and lipid metabolism transcriptional control (Fig. 5D).

In primary and metastatic tumors, Mono differentiates into Mac2, proliferative cells mature into mDC, and cDC with neutrophils in the late stage (Fig. 5E). At the start of tumor growth, FCRL1, TESPA1, RND3, and NFKBIA regulate inflammatory responses, cell polarity, and immune evasion. In particular, NFKBIA regulates NF-κB signaling pathway activity for tumor cell survival, whereas CCL3 and L1B genes highlight the importance of inflammatory mediators in the tumor microenvironment. ADAM12 and WNT5B expression affects cell adhesion, migration, and wnt signaling pathways, which are crucial for tumor cell invasiveness and microenvironment interactions. These findings reveal distinct molecular mechanisms for cell differentiation, immunoregulatory signaling, and cell–cell interactions in primary tumors. In metastatic samples, elevated SPIB and POLB expression is associated to immune responses and DNA repair, and LILRA4 and IRF4 may boost immunomodulation. ADAM12 shows extracellular matrix remodeling and cell migration during tumor growth. CCL18 and NFKBIA activation during metastasis may also boost tumor microenvironment immunosuppression and inflammation. These findings illuminate the molecular adaptations to new surroundings, enabling their survival, invasion, and dissemination in metastatic cancers (Fig. 5G).

BX1, HMGA1, and TEAD4 may promote tumor growth in transcriptional regulatory networks. In primary tumors, IRF8 promotes cell proliferation, preserves tumor stem cell features, and regulates cell–cell interactions. IRF7 and NFATC2 expression modulates the immunological landscape and immune evasion techniques during tumor cell interactions with the host's immune system, and spreading ability and adaptation depend on RUNX2, NR3C1, ZEB1, and NFKB1. Metastatic cancers express ZEB1, a key Epithelial-mesenchymal transition (EMT) regulator, to increase migration and invasion. ETS2, CEBPB, and MEF2A alter cell differentiation and cell–cell signaling to facilitate tumor metastasis. These factors regulate immunological modulation, cell cycle, cell fate determination, DNA damage response, and cell migration and invasion in primary and metastatic cancers (Fig. 5F).

Cell–cell communication heterogeneity in primary and metastatic nasopharyngeal carcinoma

Cell communication across 30 tumor cell populations (Fig. 6A) showed substantial differences in CAF, pro-T, and B-IL7R cell contacts in primary versus metastatic tumors (Fig. 6B). This discovery brought chemokine, co-inhibitory, and co-stimulatory molecule communications under scrutiny. The interaction between TNF superfamily members, their receptors, and immunomodulatory molecules is crucial in both tumor phases. These receptor pairs shape the tumor microenvironment, cell signaling, immunological control, and life cycles, as well as tumor cell motility and metastasis. TNFSF14-TNFRSF6B and TNFSF13-TNFRSF17 receptor pairings may boost local immune responses, helping tumor cells survive and proliferate. Additionally, the TNF-TNFRSF interaction can increase inflammation via the NF-κB pathway, promoting tumor growth. ICAM1 and VCAM1 linked to integrins make tumor cells tightly integrated with the milieu, facilitating local tumor propagation (Fig. 6C).

Fig. 6

Cell–cell communications in primary and metastatic NPC. A The network diagram compares cell–cell communications among 30 different groups in primary and metastatic tumors, with the thickness of the lines indicating the weight of interaction between subgroups: thicker lines represent more communication, while thinner lines indicate less. B The heatmap illustrates the volume of cell–cell communication among multiple populations within two tumor groups, with red indicating relatively higher levels of communication and blue denoting relatively lower levels. C The dot plot presents the populations with the greatest quantitative differences in receptor-ligand interactions, namely CAFs, pro-T, and B-IL7R cells, highlighting the marked variations in chemokine, co-inhibitory, and co-stimulatory receptor dynamics among these subgroups

Tumor cells need CXCL12-CXCR4/CCR7 axis activation to proliferate, survive, and metastasis to new settings. The CXCL12-CXCR4 axis drives tumor cell migration to CXCL12-rich locations like lymph nodes, lungs, and bone marrow, facilitating metastasis (Fig. 6C). Thus, cell communication receptors are crucial to initial tumor genesis, maintenance, and metastasis. They enhance tumor cell survival and dissemination by changing signaling pathways and immunological responses, enabling tumor cells to thrive in new surroundings. Understanding these mechanisms provides tumor therapy molecular targets for customized treatment.

Exploring cellular heterogeneity and molecular characteristics of metastatic and non-metastatic primary tumors

To explore the differences between primary tumors with and without metastasis, we analyzed single-cell RNA sequencing data from 2 metastatic primary tumors and 10 non-metastatic primary tumors, categorizing them into 10 distinct clusters (Fig. 7A, B). Analysis of the cell type proportions revealed that metastatic primary tumors have higher proportions of myeloid cells, dendritic cells (DCs), epithelial cells, and plasma cells, whereas non-metastatic primary tumors predominantly contain T cells and B cells (Fig. 7C). The cell types identified include: CD4+ T cells (CD3D, CD3E), CD8+ T cells (CD3D, CD3E, CD8A, CD8B, GZMK), B cells (MS4A1, CD79A, CD79B, CD19), macrophages/monocytes (CD14, S100A8, C1QA, CD163, MRC1), plasma cells (MZB1, IGHG1, IGHA1, TNFRSF17), plasmacytoid DCs (IL3RA, LILRA4, CLEC4C), migratory DCs (LAMP3, EBI3, PDCD1LG2), mast cells (TPSAB1, CPA3, MS4A2, TPSB2), cancer-associated fibroblasts (CAFs) (COL1A1, COL3A1, DCN, FAP), and epithelial cells (KRT19, EPCAM, KRT18) (Figure D). In examining the differential gene expression between metastatic and non-metastatic primary tumors, we focused on the top 10 genes with the largest expression differences (Fig. 7E). Notably, CCL3L1 is highly expressed in metastatic primary tumors and plays a crucial role in regulating immune cell responses to tumors, affecting cell migration and invasion within the tumor microenvironment, potentially promoting or inhibiting tumor metastasis (Fig. 7F). Additionally, we conducted Gene Ontology (GO) enrichment analysis on the top 100 differentially expressed genes for both groups. The findings suggest that metastatic tumor cells are more active in energy and nucleotide metabolism, likely correlating with their enhanced proliferation and survival capabilities in new environments; meanwhile, non-metastatic primary tumor cells are more focused on maintaining protein stability and modulating functions related to immune interactions with the host, which may support their prolonged survival at the primary site (Fig. 7G). In our transcription factor analysis, several key transcription factors, including MYC, SOX2, ETS2, GATA1, and GATA2, are found to be highly expressed in metastatic primary tumors. These factors are involved in regulating essential processes such as cell proliferation, migration, and angiogenesis, and are closely linked to the aggressiveness and metastatic potential of various cancers (Fig. 7H). In summary,these results elucidate the cellular and molecular distinctions between metastatic and non-metastatic primary tumors, revealing potential targets for therapeutic intervention.

Fig. 7

Comparative analysis of cellular and molecular profiles bewteen metastatic and non-metastatic primary tumors. A The umap plot displaying 10 distinct clusters, identifying various cell types such as CD4+ T cells, CD8+ T cells, B cells, Macrophages/Monocytes, plasma cells, pDCs, migratory DCs, mast cells, cancer-associated fibroblasts (CAFs), and epithelial cells. B Integration of single-cell RNA sequencing data, with metastatic primary tumors represented in red and non-metastatic primary tumors in green. C Proportions of cell populations from metastatic and non-metastatic primary tumors across the identified 10 clusters. D A dot plot illustrating the expression levels of marker genes across the 10 clusters, with high expression shown in red and low expression in blue. E A violin plot highlighting the top 10 genes exhibiting the most significant expression differences between metastatic and non-metastatic primary tumors. F A FeaturePlot visualizing the distribution of CCL3L1 expression within the two groups, high expression (Red), low expression (blue). G A heatmap displaying the differential transcriptional regulation between the two tumor types