The experimental design and informed consent procedures were approved by the Ethics Committee of Beijing Tongren Hospital, affiliated with Capital Medical University (Ethics Approval Number: TREKY2020-021). Written informed consent was obtained from all patients participating in this study, and all procedures were conducted in accordance with the Declaration of Helsinki. Surgical procedures were performed on ACC patients who met the surgical indications, and corresponding regions of tumor center and adjacent non-tumorous tissues were collected. In this study, we collected 32 primary tumor samples and 47 lung metastasis samples, including three matched pairs from the same patients. The diagnosis of ACC was histopathologically confirmed by at least two pathologists. Histological grading was conducted based on the pathological subtypes proposed by Szanto et al. in 1984 (26): grade I for predominantly cribriform or tubular patterns, grade II for mixed patterns with cribriform and tubular components with less than 30% solid areas, and grade III for predominantly solid components exceeding 30%. The diagnostic criteria for high-grade transformation were based on the characteristic histological features proposed by Seethala et al. (27).

Appropriate tissue samples were collected into corresponding numbered grinding tubes, and 1.5 mL of TRIzol lysis reagent was added. RNA sequencing (RNA-seq) was performed at the Beijing Genomics Institute. Total RNA was extracted using an RNA extraction kit and reverse-transcribed into cDNA using a reverse transcription kit. The concentration, RNA integrity number, 28S/18S ratio, and fragment size of total RNA were determined using the Agilent 2100 Bioanalyzer. Samples meeting quality standards were used for library construction and sequencing on an MGISEQ2000RS platform. Subsequent data analysis, visualization, and mining were performed using the Dr. Tom multi-omics data mining system (https://biosys.bgi.com). Differential gene expression analysis between groups was conducted using DESeq2 (28) with criteria of Fold Change ≥ 1 and Adjusted P value ≤ 0.05. Differentially expressed genes were functionally classified based on Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) annotations, as well as official classifications. KEGG enrichment analysis was performed using the phyper function in R software, and GO enrichment analysis was carried out using the TermFinder package. Genes with Qvalue ≤ 0.05 were considered significantly enriched in the candidate gene set.

Single-sample gene set enrichment analysis (ssGSEA) was performed using the GSVA package in R software (version 4.3.1) (29). Marker genes characterizing immune cell types were obtained from Bindea et al. (30) and Charoentong et al. (31), while angiogenesis marker genes were sourced from Masiero et al. (32). The marker gene set for MHC class I antigen presenting machinery (APM) included HLA-A, HLA-B, HLA-C, B2M, TAP1, TAP2, and TAPBP (33). The overall immune infiltration score (IIS) was defined as the average of the standardized values of innate and adaptive immune scores (33). The T-cell infiltration score (TIS) was defined as the average of the standardized values of nine T-cell subtypes (33). Cytolytic activity (CYT), an indicator of immune cell cytotoxic activity, was defined as the geometric mean of the transcriptional expression levels of the effector genes GZMA and PRF1 (34). Tumor Immune Dysfunction and Exclusion (TIDE) is a computational method based on transcriptomics that simulates two major mechanisms of tumor immune evasion: T-cell dysfunction and T-cell exclusion. It can predict cancer patient response to immune checkpoint inhibitors (35). Batch-corrected normalized data were input into the TIDE website (http://tide.dfci.harvard.edu) for calculation. Additionally, the CIBERSORT algorithm, based on deconvolution, was used for relative quantification analysis of immune cell proportions (36).

The previously generated single-cell sequencing dataset used in this study has been uploaded to the public repository Gene Expression Omnibus (GEO) under the accession number GSE216852. This dataset includes single-cell 3’-RNA sequencing data from one primary ACC patient and one lung metastasis ACC patient. The Seurat R package was utilized to analyze the ScRNA-seq data according to standard analysis procedures. After log normalization and dimensionality reduction of the differential gene expression data of macrophage populations, the cells were divided into six clusters. Uniform Manifold Approximation and Projection (UMAP) plots, violin plots, and bubble plot analysis results were generated using the Monocle 2.0 package.

Paraffin-embedded tissue sections (4 μm thick) were deparaffinized in fresh xylene and subjected to antigen retrieval. Endogenous peroxidase activity was quenched by covering the sections with 3% hydrogen peroxide for 15 minutes. Tissue sections were blocked with 10% goat serum at 37°C for 1 hour, then incubated overnight at 4°C with primary antibodies against CD8 (ZSGB-Bio, ZA-0508), CD4 (ZSGB-Bio, ZA-0519), and CD68 (PTM-BIO, PTM-5130). Detection of horseradish peroxidase (HRP) activity was performed using the PV-6000D immunohistochemistry kit (ZSGB-Bio). Sections were counterstained with hematoxylin, differentiated with hydrochloric acid, dehydrated, and mounted. Rabbit or mouse monoclonal IgG was used as a negative control. Images were captured using the Pannoramic^® 250 FLASH scanner. The number of positive cells per unit area was quantified using Image J software (version 13.0.6, National Institutes of Health, USA).

Data were analyzed using GraphPad Prism software (version 9.5.1, GraphPad Software, CA, USA) and R software (version 4.3.1). The Mann-Whitney U test was used to compare the median distributions between two sample groups. Spearman correlation analysis was employed to assess the correlations between different indices. Lung metastasis-free survival was estimated using the Kaplan-Meier method, and statistical differences were determined using the Log-rank test. Statistical significance was indicated as follows: P<0.05 (*), P<0.01 (**), P<0.001 (***), P<0.0001 (****).

To elucidate the heterogeneity in the tumor immune microenvironment (TIME) between primary tumors (PT) and lung metastases (L-MET) of ACC, we conducted a comprehensive analysis of both by transcriptome sequencing, single-cell sequencing, and immunohistochemistry ( Figure 1A ). The basic clinical characteristics of the cohort are summarized in Supplementary Table S1 . Transcriptome sequencing was performed on 25 primary tumor samples and 11 matched adjacent non-tumor primary tissues, as well as 34 lung metastasis samples and 19 matched normal lung tissues. Single-sample gene set enrichment analysis was conducted based on the characteristic gene sets of 28 immune cell types. The heatmap displayed the abundance of various immune cells across all samples, revealing distinct differences between different sample types ( Figure 1B ). Adjacent normal lung tissues exhibited a richer immune contexture due to the enrichment of T cells, B cells, and innate cells ( Figures 1C–E ). Both lung metastases and primary tumors displayed relatively low immune infiltration environments. Additionally, the APM and CYT score of tumor samples were lower than those of adjacent normal tissues, accompanied by elevated expression of the immune checkpoint TIGIT ( Figures 1F–H ). Overall, the abundance and functions of immune cells in the TIME of primary and metastatic lesions were significantly lower than those in adjacent normal tissues.

Next, we compared the differences in immune cell types, functional activation, and inhibition levels between primary tumors and lung metastases, revealing a highly distinct immune landscape ( Figures 2A, B ; Supplementary Figure S1A ). The IIS was higher in L-MET than in PT ( Figure 2C ). Although the TIS showed no significant difference ( Figure 2D ), the composition of specific immune lineages varied ( Figure 2B ). Effector memory CD4 T cells were significantly more abundant in PT (P<0.0001), while effector memory CD8 T cells were significantly less abundant in PT compared to L-MET (P<0.01). Certain immune cell infiltrations, such as activated B cells, effector memory CD8 T cells, and Th1 cells, were markedly higher in lung metastases, suggesting a stronger adaptive immune response. Conversely, immunosuppressive immature dendritic cells and macrophages were also more prevalent in L-MET and positively correlated with activated CD8 T cells, natural killer (NK) cells, and natural killer T (NKT) cells ( Figure 2B ; Supplementary Figure S1B ). The anti-tumor immune response partially depends on the antigen-presenting capacity. Our results showed that the APM score was significantly higher in L-MET compared to PT (P<0.05) ( Figure 2E ). The expression of HLA class I molecules and HLA class II molecules was significantly higher in L-MET ( Supplementary Figures S2A–G ). However, there was no substantial difference in CYT score between PT and L-MET (P=0.15) ( Figure 2F ), suggesting additional immunosuppressive mechanisms may exist in lung metastases. We further analyzed T-cell dysfunction and T-cell exclusion scores, finding that both scores were significantly higher in L-MET (P<0.05) ( Figures 2G, H ). The expression of the immune checkpoint molecule CTLA4 was also higher (P<0.05). CD274 (encoding PD-L1) expression was slightly higher in lung metastases, though not statistically significant (P=0.0504).

In summary, our study reveals heterogeneity in the tumor immune microenvironment between primary head and neck ACC and its lung metastases. The differences in immune cell composition, immune responses function, and inhibitory mechanisms between lung metastases and primary tumors are substantial. These differences should be carefully considered when designing preclinical drug development studies targeting lung metastases.

Our transcriptome sequencing analysis revealed that the overall immune cell infiltration is higher in lung metastases compared to primary tumors. However, the TIME also exhibits elevated immune suppression mechanisms in L-MET. To further validate the RNA-seq results and explore the specific cell types exerting immunosuppressive effects, we performed immunohistochemical validation on 30 tumor samples. Consistent with the sequencing results, the number of CD8+ and CD68+ positive immune cells infiltrating lung metastases was significantly higher than in primary tumors (P<0.05 and P<0.01, respectively) ( Figures 3A, B ). We found a remarkable enrichment of macrophages in the stroma of lung metastases ( Figure 3A ), which were positively correlated with CD8+ T cells (Spearman r=0.4588, P<0.01) ( Figure 3C ). Tumor-associated macrophages typically have immunosuppressive and pro-tumorigenic roles in many cancers (37). These findings suggest that macrophages may be a critical component of the immunosuppressive microenvironment in lung metastases of ACC.

Tumor immune phenotypes are often correlated with malignant progression and response to immunotherapy (38, 39). Based on the spatial distribution of T-cell infiltration, TIME can be classified into three types: immune-inflamed (I-I), immune-excluded (I-E), and immune-desert (I-D) (38). We examined the TIME classification in PT and L-MET, finding that lung metastases predominantly exhibited the I-I phenotype ( Figure 3D ), which is consistent with our previous findings (40). CD68+ immune cells were predominant across different TIME classifications ( Figure 3E ). The distribution of different histological grades in PT and L-MET samples was relatively uniform ( Figure 3F ), with higher densities of CD68+ immune cells in Grade II and Grade III tissues containing solid components ( Figure 3G ). Solid components are a key indicator of poor prognosis in ACC patients and are closely associated with distant metastasis (10, 41). Further investigation is needed to determine whether the malignant progression of solid-type tumors is related to tumor-infiltrating macrophages.

Tumor-associated macrophages exist in various subtypes, some of which have been shown to correlate with cancer prognosis and resistance to immunotherapy (42). To further investigate the subtypes and functions of macrophages in ACC, we utilized previously published single-cell sequencing data from our team (43). We performed single-cell sequencing on a primary tumor (A) and matched adjacent non-tumor tissue (AP), as well as on a set of lung metastases of different sizes (A1, B1, C1) and matched adjacent lung tissue (F), identifying macrophage populations ( Figure 4A ). Further sub-clustering revealed six macrophage types, with cluster 0 identified as tumor-infiltrating macrophages, and clusters 1 to 4 primarily representing macrophages from adjacent tissues ( Figure 4B ). Analysis of specific macrophage markers showed that cluster 5 highly expressed the M1 marker CD80, while clusters 0 to 4 predominantly expressed M2 markers CD163, CD206, CSF1R, PTGS2 ( Figure 4C ). We further analyzed the top five highly expressed genes in each macrophage cluster ( Figure 4D ), most of which were related to immunoregulatory functions and tissue repair of macrophages. Genes such as MARCO, FN1, ACE, and CD163L1 are likely associated with the immunoregulatory functions of M2 macrophages, which typically promote tissue repair, anti-inflammatory responses, and support tumor growth (44, 45). Cluster 5 highly expressed DNASE1L3, a gene involved in DNA degradation and potentially in the clearance of self-antigens during cell death processes (46). While cluster 5 may act as anti-tumor macrophages, they represent only a small fraction.

Overall, our IHC and ScRNA-seq analyses confirmed that the predominant cell population in the TIME of ACC consists of M2 macrophages, especially in lung metastases, where they exhibit significant immunosuppressive functions.

We performed a relative quantification analysis of immune cells using the CIBERSORT deconvolution algorithm (36). Consistent with previous results, macrophages were the major component of the myeloid immune cell subpopulation, with M1 macrophages representing a small fraction and M2 macrophages predominating in both primary and lung metastatic lesions ( Figure 5A ). Among the three pairs of matched primary and lung metastatic tumors, two cases showed significant amplification of M2 macrophages in lung metastases, while one case showed a significant increase in plasma cells ( Supplementary Figure S3A ).

Based on the relative proportions of M1 and M2 macrophages, we categorized the population into M1_high, M1_low, M2_high, and M2_low groups according to their respective medians. The M1_high group had higher CD8 T-cell levels compared to the M1_low group (P<0.05), but lower levels of activated dendritic cells (DCs) (P<0.0001) ( Figure 5B ). The proportion of M1 macrophages showed a positive correlation with CD8 T-cell proportion (Spearman r=0.4358, P<0.001) and a negative correlation with activated DCs (Spearman r=-0.6436, P<0.0001) ( Supplementary Figure S3B ). The M2_low group exhibited significantly higher levels of plasma cells and activated NK cells compared to the M2_high group (P<0.05), though activated DCs were lower in the M2_low group, albeit not significantly (P=0.1502) ( Figure 5C ). These findings suggest that the reduction of activated DCs may represent another immune escape mechanism in patients with an anti-tumor TAM phenotype. In line with the anti-tumor function of M1 macrophages, the M1_high group had elevated IIS, TIS, APM scores, and CYT scores compared to the M1_low group (all P<0.001) ( Figure 5D ), indicating a stronger immune response. However, the expression levels of immune checkpoint molecules CD274, HAVCR2, CTLA4, and LAG3 were also higher in the M1_high group (P<0.01) ( Figure 5F ). The M2_low group showed higher TIS and APM scores compared to the M2_high group (all P<0.05), but there was no significant difference in CYT scores (P=0.1127) ( Figure 5E ). CTLA4 expression was higher in the M2_low group (P<0.001) ( Figure 5G ).

In conclusion, categorizing ACC patients based on macrophage transcript levels revealed differences in their immune landscape, including variations in anti-tumor immune cell infiltration, antigen presentation capability, cytolytic activity, and immune checkpoint molecule expression. These results suggest that TAMs may serve as potential biomarkers for predicting immune responses in ACC patients.

Due to the unique microenvironment of the lungs, circulating tumor cells (CTCs) of ACC interact with their surroundings, forming a metastatic niche distinct from the primary site. We next conducted a detailed analysis of ACC lung metastases. Unsupervised clustering of lung metastasis samples based on immune cell infiltration scores divided them into two clusters: Cluster I (immune cell-enriched) and Cluster II (immune cell-poor) ( Figure 6A ). Cluster I had significantly higher IIS and TIS compared to Cluster II (P<0.0001 and P<0.05, respectively) ( Figures 6B, C ). Differentially expressed genes between the two clusters were enriched in pathways related to immune response, chemokine signaling, T-cell activation, antigen processing and presentation, and primary immunodeficiency ( Supplementary Figures S4A-C ), implying that Cluster II may have defects in initiating immune responses and recruiting immune cells.

When comparing the lung metastasis-free survival between Cluster I and Cluster II, there was no significant difference (P=0.1168, Log-rank test) ( Figure 6D ). This indicates that even though Cluster I patients seemed to have a more favorable immune microenvironment, it did not translate into better metastasis-free survival. Due to the short follow-up period, no patients had reached the endpoint event, so overall survival differences between the two groups require longer follow-up to observe. Based solely on immune cell abundance, we did not observe prognostic differences in lung metastasis outcomes among ACC patients. Nevertheless, when stratifying patients based on the proportion of M2 macrophages into M2_high and M2_low groups, the M2_low group had significantly longer metastasis-free survival (P<0.05, Log-rank test) ( Figure 6E ). This suggests that high infiltration of M2 macrophages in tumors may be an indicator of early lung metastasis.

In lung metastases, activated CD8+ T cells, NK cells, and activated DCs were all positively correlated with macrophage scores (all P<0.0001) ( Figure 6F ), indicating that the infiltration of anti-tumor immune cells is often accompanied by a higher proportion of macrophages, thereby limiting their immune effects. Additionally, IIS scores were significantly positively correlated with T-cell exclusion scores (P<0.01) ( Figure 6G ). The expression of immune checkpoint molecules VSIR (P<0.0001), HAVCR2 (P<0.0001), LAG3 (P<0.0001), and the immunosuppressive cytokine TGFB1 (P<0.01) were also significantly positively correlated with IIS scores ( Figures 6H, I ). These results indicate that in lung metastasis patients with immune-enriched tumors, multiple immunosuppressive mechanisms may be present, including the enrichment of TAMs, T-cell exclusion, high expression of immune checkpoint molecules, and the production of immunosuppressive cytokines.