ScRNA-seq reveals the correlation between M2 phenotype of tumor-associated macrophages and lymph node metastasis of breast cancer

The process of lymphatic metastasis was proved to be associated with podoplanin-expressing macrophages in breast cancer (BC). This study aimed to investigate the role of the M2 phenotype of tumor-associated macrophages and mine the key M2 macrophages-related genes for lymph node metastasis in BC. We downloaded the GSE158399 dataset from the Gene Expression Omnibus (GEO) database, which includes transcriptomic profiles of individual cells from primary tumors, negative lymph nodes (NLNs), and positive lymph nodes (PLNs) of breast cancer patients. The cell subsets were identified by clustering analysis after quality control of the scRNA-seq using Seurat. The activation and migration capability of M2 macrophages were evaluated with R package “GSVA”. The key M2 macrophages-related genes were screened from the differential expressed genes (DEGs) and M2 macrophages activation and migration gene sets collected from MSigDB database. Our analysis identified three main cell types in primary tumors, NLNs, and PLNs: basal cells, luminal cells, and immune cell subsets. The further cell type classification of immune cell subsets indicated M2 macrophages accumulation in NLs and PLs. The GSVA enrichment scores for activation and migration capability were increased significantly in M2 macrophages from primary tumors than NLNs and PLNs (p-value < 0.001). Seven M2 macrophages activation-related and 15 M2 macrophages migration-related genes were significantly up-regulated in primary tumors than NLNs and PLNs. The proportion and GSVA enrichment scores for activation and migration of M2 macrophages may be potential markers for lymph node metastasis in breast cancer. Our study demonstrated that twenty-two up-regulated mRNA may be possible therapeutic targets for lymph node metastasis in breast cancer.


Introduction
Breast cancer (BC) is the most frequently diagnosed neoplasm and the leading cause of cancer mortality among women [1].In 2022, there will be estimated 290,560 new cases and 42,780 cancer deaths due to BC in the USA [2].BC represents the most important cancer-related cause of disease burden worldwide, especially in developed countries [3].Despite the continuous improvement of BC treatments in recent years, the survival of advanced BC is still not ideal.Lymph node metastasis is one of the main pathways of tumor metastasis, especially for breast tissue with abundant lymphatic vessels and lymphatic network [4].The metastatic cancer cells in the lymph nodes can directly enter other lymphatic vessels and then spread widely in the body [5].Lymph node metastasis is of great significance to BC's prognosis, indicating the worsened prognosis [6].Lymph node metastases is a complex biological process affected by a complex gene regulatory network and various growth factors, involving tumor movement, vascular invasion, and clonogenicity in the microenvironment [7].In the process of lymph node metastases, lymph nodes provide a supportive environment for tumor cells of a specific genetic background, supporting their clonal growth and further distant metastasis [8].
The process of lymphatic metastasis in breast cancer has been found to be associated with a specific type of macrophages expressing podoplanin [9].Macrophages display two different phenotypes in response to different environmental stimuli: M1 and M2 macrophages.M2 macrophages can produce large amounts of cytokines that cause Th2-type immune responses, suppressing immune function in the tumor microenvironment, inducing angiogenesis, and supporting tumor growth and metastasis [10].In the tumor microenvironment, the majority of tumor-associated macrophages (TAMs) are of the M2 phenotype [11].Increased numbers of M2 macrophages have been significantly correlated with lymph node metastasis, larger tumor size, poor differentiation of cancer cells, and an elevated risk of recurrence in BC [12].The experimental data of Watari et al. showed that M2 macrophages correlated with lymph node metastasis in highly metastatic cancer [13].Therefore, it is crucial to further investigate the biological characteristics of M2 macrophages in the context of breast cancer lymph node metastasis.Understanding their role and mechanisms in promoting metastasis could potentially lead to the development of targeted therapies.
With the development of sequencing technology, singlecell RNA sequencing (scRNA-seq) has brought new data at the cellular level to researchers, enriching the tools for single-cell analysis of tumors.The scRNA-seq enables researchers to explore the different biological properties of individual cells in complex tissues and understand the response of cell subpopulations to environmental elements.It has provided new solutions to various omics-related problems in the life sciences, and the related research has become increasingly popular.In this study, we collected a set of single-cell sequencing data (GSE158399) to investigate the role of M2 phenotype of tumor-associated macrophages and mine the key M2 macrophages-related genes for lymph node metastasis in breast cancer.

Data preparation and quality control
The single-cell sequencing data GSE158399 of lymph node metastasis in breast cancer was obtained from the GEO database (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE158399).The paired single-cell data of three sources in a female patient with luminal B subtype were detected in the sequencing platform of HiSeq X Ten, including primary tumors (GSM4798908), NLNs (GSM4798910), and PLNs (GSM4798909).Table 1 shows the detailed cell and gene number statistics.
The number of detected cells and captured genes in the primary tumors, NLNs, and PLNs were counted.The intersection of the detected genes in the three tissues was collected, and the number of genes in every single cell of each sample was calculated.The threshold for quality control was determined by the overall distribution.The single cell and genes with low quality were eliminated using."merge.SCTAssay" function in the R package "Seurat".

ScRNA-seq clustering analysis
To reduce genomic instability caused by the single cell quality, batch effect and retain most of the gene expression information, the R package "Seurat" was used to perform a series of data processing and analysis on the expression data of single cells.Based on this R package, we performed PCA dimensionality reduction, high-variant gene selection, and K-Nearest Neighbor clustering, and finally divided single cells into more detailed cell subsets and analyzed them.

GSVA-gene set enrichment analysis
Relevant gene sets were collected through the MSigDB (https://www.gsea-msigdb.org/gsea/msigdb)database, and the R package "GSVA" was used to score the functional activity of M2 macrophages from the three samples using the Wilcoxon rank sum test compares differences among them.

Screening and functional annotation of differentially expressed genes (DEGs)
The mean expression of each gene in the three samples was counted and the FC value was calculated.Subsequently, the Wilcoxon rank-sum test was performed on the single cells from the two types of samples, and the corresponding test pvalue and the p-value corrected by the FDR algorithm were recorded, and genes with FDR < 0.001 and |log2FC| > 1 were retained.Finally, the DAVID (https://david-d.ncifcrf.gov/) database was used to analyze the KEGG signaling pathway enrichment of the above differentially expressed genes to understand the biological processes, signaling pathways, molecular functions, and cellular localization.The number of genes > 10 and p ≤ 0.05 were considered statistically significant.

Data quality control results
High-quality, usable single-cell data were obtained after quality control processing.24235 common genes were finally obtained after the intersection of the captured genes in the three tissues (Table 1).According to the overall distribution of the number of captured genes in single cells of the three samples, it is considered that the number of genes captured by cells in GSM4798908/GSM4798909/ GSM4798910 is less than 600/200/100 genes are lower quality cells (Fig. 1A).Genes expressed in fewer than 15 cells were considered lower-quality genes based on the overall distribution of cell numbers expressed in each gene (Fig. 1B).After the above screening, 33857 single cell expression data were obtained by using the "merge.SCTAssay" function in the R package "Seurat".Through the overall distribution of single cells from the three source samples, the proportion of mitochondrial genes in the integrated single-cell data was counted, and cells with mitochondrial gene expression accounting for more than 20% were excluded [14].Ultimately, 17,746 high-quality genes and 32,079 single-cell samples were obtained (Fig. 1C).

Major cell types in different tissues
In order to obtain a more detailed analysis of single cells, The 32,079 single cells were divided into 26 cell subsets using the R package "Seurat".By comparing the tissue sources of the single-cell samples, it was found that the single cells in the GSM4798909 and GSM4798910 samples were generally immune leukocytes.In contrast, a small number of cells in the primary tumor GSM4798908 were immune leukocytes (Figs.2A and 2B).On this basis, three main cell types present in breast tissue were identified according to gene markers expressed explicitly in cell subsets: Basal cells (VIM + , ITGB1/CD29 + ), Luminal cells (KRT19 + ), and immune cells Cell subsets (PTPRC/CD45 + ) [15][16][17].In addition, a small group of mixed cells expressed both Basal cells and Luminal cell markers (Figs. 2B and 2C).Since KRT14 and KRT5 were not detected in the dataset, they could not be used as identification criteria for Basal cells.

M2 macrophages accumulate in PLNs
Cell type division of the immune cell subset Cluster 13 present in single cells derived from the GSM4798908 sample.The results showed that there were a small number of T cells (62) and NK cells (14) in this cell population (850 cells in Cluster 13).Still, most of them were macrophages that specifically expressed CD68, especially M2 macrophages that specifically expressed FCGR2B, and a few were M1 macrophages specifically expressing MRC1 (Figs. 3A and  3B).Second, M2-type macrophages were significantly more numerous than M1-type in GSM4798909 and GSM4798910 tissues (Fig. 3C).In lymph node metastasis positive tissue GSM4798909-derived single cells, the number of M2-type macrophages was 6.7708 times that of M1-type macrophages.In single cells derived from GSM4798908 and GSM4798910, the number of M2 macrophages is 3.4416 and 5.8000 times that of M1 macrophages, respectively (Table 2).
The ability of M2 macrophages to migrate and activate is significantly enhanced in orthotopic tumor tissue For this study, we focused on the ability of M2 macrophages to migrate and activate cellular immune responses.For this study, we focused on the ability of M2 macrophages to migrate and activate cellular immune responses.We used the "GOBP_MACROPHAGE_ACTIVATION_INVOLVED_ IN_IMMUNE_RESPONSE" gene set from the MSigDB database to assess macrophage activation, and the "GOBP_MACROPHAGE_MIGRATION" gene set to evaluate macrophage migration.The results demonstrated that M2 macrophages in orthotopic tumor tissue exhibited significantly enhanced activation ability compared to M2 macrophages in NLNs and PLNs.The statistical analysis yielded test p-values of less than 2.2e-16, indicating a highly significant difference (Fig. 4A).Furthermore, the migration ability of M2 macrophages in orthotopic tumor tissue was found to be significantly higher than that of M2 macrophages in NLNs and PLNs.The statistical analysis resulted in a test p-value of less than 2.2e-16, indicating a highly significant difference (Fig. 4B).These findings suggest that M2 macrophages in orthotopic tumor tissue possess a greater propensity to migrate and activate cellular immune responses compared to M2 macrophages in the nearby and peripheral lymph nodes.The enhanced migratory and activation abilities of M2 macrophages in the tumor microenvironment may contribute to the regulation of immune responses and potentially influence tumor progression and immune surveillance.

Differential analysis of key genes related to M2 macrophages
According to the above signal pathway results, the promotion of lymph node metastasis by M2 macrophages may be related to biological characteristics such as phagosomes and immunodeficiency.Therefore, through the statistics of relevant genes in "gobp_macrophage_activation_invoved_ in_immune_response", we found that there are seven significantly up-regulated genes, namely SUCNR1 (Succinate Receptor 1), TREM2 (Triggering Receptor Expressed On Myeloid Cells 2), TYROBP (Tyrosine Kinase Binding Protein), GRN (Granulin Precursor), HAVCR2 (Hepatitis A Virus Cellular Receptor 2), IFI35 (Interferon Induced Protein 35), NMI (N-Myc And STAT Interactor).In contrast, the SUCNR1 gene was not expressed in the two lymph node-derived samples, The TREM2 gene was not expressed in lymph node metastasis-negative samples (Table 5).To visualize the differential expression of these seven upregulated genes, a heat map was generated (Fig. 6A).

Discussion
Lymph node metastasis is a common problem in breast cancer, which seriously affects the survival and prognosis of patients.Therefore, it has become increasingly urgent to understand its pathogenesis and find a specific and sensitive treatment method.The mobility and invasiveness of cells are the keys to tumor metastasis, and the movement and migration of tumor cells play a leading role in the whole process of tumor metastasis.The study conducted by Watari et al. [18] stated that inflammatory stimuli could help to establish the tumor microenvironment, and the growth, invasion, and matastases could be induced through the activation of macrophages.In our study, we observed a fundamental distinction between samples obtained from primary tumors and lymph nodes.The single-cell analysis of lymph node samples revealed a predominant presence of immune leukocytes, whereas immune leukocytes constituted only a small portion of the primary tumor samples.Further classification of cell types based on single-cell data demonstrated that the majority of cells in the primary tumors were macrophages expressing the CD68 marker, particularly M2 macrophages expressing FCGR2B.Additionally, Mahmoud et al. [19] discovered that higher numbers of CD68+ macrophages were associated with lower overall survival rates in breast cancer patients.This finding underscores the potential significance of M2 macrophages in lymph node metastasis and its clinical implications.
Tumor-associated macrophages (TAMs) play a key role in the growth of breast cancer.The role of TAMs in cancer pathogenesis depends on their phenotypic and functional polarization [20].Tumor cells can disguise as normal cells to deceive macrophages and inhibit the phagocytosis of tumor cells by macrophages [21].Studies have shown that the high infiltration of TAMs is essential in promoting breast cancer cell metastasis.The infiltration density of TAMs in primary breast tumor tissue is significantly higher than in adjacent tissue.At the same time, the infiltration of TAMs in primary breast tumor tissue increases with tumor stage and size [22].Previous studies have shown that infiltration of macrophages in the tumor microenvironment supported tumor growth, angiogenesis, metastasis/invasion, inflammation, and immunosuppression by secreting factors that promote tumor progression [23].TAMs are heterogeneous populations with different subsets performing different functions [24].TAMs primarily exhibit a phenotype resembling alternatively activated M2 macrophages, which possess anti-inflammatory properties and promote tumor growth.M2-type macrophages largely contribute to the aggressiveness of malignant tumors [18].TAMs located in primary tumors are mostly associated with poor prognosis in cancer patients.In this study, we found that M2 macrophages mainly aggregated in breast cancer lymph node metastasis-positive samples, and M2 macrophages in orthotopic tumor tissue had stronger migration and activation ability.This phenomenon suggests to us that M2 macrophages play an important role in promoting lymph node metastasis in patients.Depletion or phenotypic reversal (M2 to M1) of TAMs has been shown to halt tumor progression in mouse models of breast cancer [25].In order to reveal the biological role of M2 macrophages in the process of lymph node metastasis, we analyzed the differences of M2 macrophages from three samples and annotated their functions.The results showed that M2 macrophages may help tumor cells complete metastasis through phagosomes, and M2 macrophages may also have immune function defects.New clinical research suggested that an immunodeficiency of macrophages 22 may cause Crohn's disease (CD).The reason for granulomas is the reduced ability to clear invading bacteria induced by a weakened attraction of granulocytes to the intestinal wall.So far, there are few reports on whether M2 macrophages promote tumor cell metastasis through phagosomes and immunodeficiency.
Based on the above analysis, we analyzed the difference between the gene sets of "gobp_macrophage_activation_ involvedin_immune_response" and "gobp_macrophage_ migration" in M2 macrophages, respectively.The previous study analyzed the data in GSE158399 and found the elevated expression level of CXCL14 in PLs, which might be valuable in predicting the prognosis of breast cancer with lymph node metastasis [23].In this study, we found that SUCNR1, C3aR1, C5aR1, MMP14, THBS 1, TSP-1, and TREM2, which are abnormally expressed in tumor tissue in situ, may be potentially associated with lymph node metastasis.
Keiran et al. [26] discovered that the activation of SUCNR1 plays a crucial role in promoting the antiinflammatory phenotype of macrophages, thus contributing to the anti-inflammatory response.It is demonstrated that Succinate and its receptor SUCNR1 can suppress immune responses.Moreover, the receptor SUCNR1 can also suppresse immune responses, and SUCNR1 deficiency in macrophages could lead to enhanced inflammatory responses [27].C3AR1 was down-regulated in osteosarcoma tissues and cells, and its overexpression inhibited the proliferation, migration, and invasion of osteosarcoma cells and induced apoptosis [28].These findings may appear contradictory to our conclusions, but they emphasize that the interaction between C3AR1 and cell migration can have divergent roles in different types of cancer.In addition, C5α can also promote the malignant development of HBc-positive hepatocellular carcinoma through C5AR1 [29].Matrix metalloproteinases (MMPs) are key factors in extracellular matrix remodeling and cell migration during tumor metastasis.In particular, MMP-14, a membraneanchored MMP, is closely involved in these processes.The study found that Bladers I and IV were associated with cell migration.Bladers IV is required for MMP-14 homodimerization.The interaction between MMP-14 and CD44 leads to phosphorylation of the EGF receptor and downstream activation of MAPK and PI3K signaling pathways involved in cell migration [30].Reports provide evidence [31] that THBS1 derived from oral squamous cell carcinoma (OSCC) exosomes is involved in the polarization of macrophages towards an M1-like phenotype.
In contrast, conditioned medium from exosomes induced M1-like TAMs and significantly promoted the malignant migration of OSCCs.Perturbation of macrophage migration inhibitory factor expression in mouse melanoma suppresses tumor formation by up-regulating Thrombospondin-1 (TSP-1) [32].TSP-1 is a secreted protein that inhibits angiogenesis, modulates anti-tumor immunity, stimulates tumor cell migration, and modulates extracellular proteases and growth factors in the tumor microenvironment.Furthermore, in polyomavirus middle T antigen (Pyt) transgenic mice, TSP-1 in the mammary tumor microenvironment inhibits angiogenesis and tumor growth, yet promotes lung metastasis in Pyt transgenic mice [33].Research has shown that targeting TREM2 on tumorassociated macrophages enhances immunotherapy [34].
However, there are potential limitations in this study.Firstly, the primary breast cancer tissues, PLs, and NLs in the GSE158399 dataset were from one patient with Luminal B subtype breast cancer, which may restrict the generalizability of the findings.Secondly, the current mainstream tissue sequencing methods may not effectively detect the key genes associated with M2 macrophages.Therefore, there could be additional genes involved in the M2 macrophage phenotype that were not identified in this study.Lastly, further investigations are needed to validate and expand upon the interesting findings presented here.In conclusion, this study found that M2 macrophages play an important role in promoting lymph node metastasis in breast cancer patients, possibly through immune activation and phagosomes to help tumor cells complete metastasis, thus exerting their biological function.In addition, we also discovered key genes that M2 macrophages significantly upregulated in immune response and cell migration to discover key molecular targets regulating breast cancer metastasis and providing new strategies for the prevention and treatment of breast cancer metastasis.These results may have important implications for understanding the mechanism of lymph node metastasis in breast cancer.the article and revising it critically for important intellectual content.All authors reviewed the manuscript and approved the version to be published.

FIGURE 1 .
FIGURE 1. Distribution of gene numbers captured from single-cell.

FIGURE 2 .
FIGURE 2. The sample source and cell subgroups in scRNA-seq data after clustering.

FIGURE 3 .
FIGURE 3. Macrophage distribution in single cells from three sample sources.

FIGURE 4 .
FIGURE 4. Differences in the activation capacity and cell migration ability of single-cell M2 macrophages.

FIGURE 5 .
FIGURE 5.The KEGG pathway enriched in the significantly DEGs.

FIGURE 6 .
FIGURE 6. Heatmap of macrophage immune responses and migration related DEGs.

TABLE 3 The
KEGG pathway enriched in the top 200 significantly highly expressed genes

TABLE 4 99
significantly down-regulated genes involved in the KEGG pathway

TABLE 6
Differential expression statistics of 15 macrophage migration-related genes