Distinct Tumor-Associated Macrophage Signatures Shape the Immune Microenvironment and Patient Prognosis in Renal Cell Carcinoma

1. Introduction

Tumor-associated macrophages (TAMs) are a critical component of the tumor immune microenvironment (TIME), playing diverse and often conflicting roles in cancer development and progression. TAMs originate from resident tissue-specific macrophages or monocytes recruited to the tumor site, where they adapt to the TIME and exhibit either pro-inflammatory or anti-inflammatory phenotypes [1,2,3]. These phenotypes, commonly referred to as M1 and M2 macrophages, respectively, influence tumor progression by modulating immune responses, promoting angiogenesis, and facilitating immune evasion [4,5,6]. Given their central role in cancer biology, TAMs have emerged as a key focus for understanding tumor behavior and developing novel therapeutic strategies.

In renal cell carcinoma (RCC), a prevalent form of kidney cancer, TAMs significantly contribute to tumor progression and patient prognosis. TAMs in RCC promote tumor growth and metastasis by enhancing epithelial-mesenchymal transition (EMT) and cell invasion through pathways such as IL-6/STAT3 signaling [7,8]. Additionally, they suppress immune responses by producing cytokines such as IL-10 and CCL2, impairing natural killer (NK) cell function and fostering an immunosuppressive microenvironment [9,10]. High infiltration of certain TAM subsets correlates with poor patient outcomes, while others may associate with improved survival, highlighting the dual roles TAMs play in RCC progression [11].

Beyond their roles in tumor biology, TAMs offer significant potential as therapeutic targets and biomarkers. Strategies aimed at modulating TAM activity include depleting TAMs, inhibiting their recruitment, or reprogramming them from a pro-tumor M2 phenotype to an anti-tumor M1 phenotype, with the goal of reducing tumor growth and enhancing the efficacy of immunotherapies [5,8,12]. Emerging approaches, such as macrophage-derived nanosponges, have shown promise in modulating TAM behavior and enhancing anti-tumor immunity in RCC [13]. Furthermore, the expression profiles of TAM-related genes and their association with immune checkpoints provide opportunities to refine patient stratification and predict responses to immune checkpoint inhibitors (ICIs) [11,14].

Despite these advances, the mechanisms by which specific TAM subtypes contribute to RCC progression remain incompletely understood. This study focuses on elucidating the roles of TAM signatures in RCC, leveraging their expression profiles to stratify patients based on risk scores, and investigating their interactions with immune components within the TME. By integrating computational deconvolution methods and survival analyses, this work aims to provide insights into the prognostic and therapeutic significance of TAMs in RCC. These findings could pave the way for the development of TAM-targeted therapies, offering new avenues for improving outcomes in RCC patients. A 27-gene risk score derived from TAM signatures across independent datasets can effectively classify patients into low- and high-risk categories, distinguish overall survival, and investigate the potential mechanisms underlying the role of TAM in the prognosis of RCC patients. The workflow of the present study is illustrated in Figure 1.

2. Materials and Methods

2.1. Data Utilization

Level 3 TCGA RNA-seq data and clinical information for RCC were obtained from TCGA on FireBrowse (gdac.broadinstitute.org/). TCGA MAF files for gene mutation analyses were obtained from https://gdc.cancer.gov/about-data/publications/pancanatlas. All genes in which non-silent mutations occurred were considered to be mutated. Macrophage regulation scores, leukocyte and lymphocyte infiltration scores, and IFN-γ response and TGFβ response scores for TCGA-KIRC samples were downloaded as a supplementary file (Supplementary Table S1) from prior work [15].

2.2. Curation of Immune-Related Genes (IRGs)

Immune-related genes (IRGs) were obtained from Charoentong et al. [16], Bindea et al. [17], and Xu et al. [18]. All genes from immune cells were collected (marker genes attributed to cancer cells were excluded) and combined into a single list of 831 IRGs genes.

2.3. Immune Cell Inference

Immune cell fractions in tumor samples were quantified using the quanTIseq method [19], an RNA-seq deconvolution algorithm implemented in R. quanTIseq applies constrained least squares regression to estimate proportions of immune cell populations, including macrophages, neutrophils, dendritic cells, CD8+ T cells, regulatory T cells (Tregs), and other immune subsets [19]. Using patient bulk gene expression data, quanTIseq enabled identification and correlation analysis between TAM signatures and specific immune cell types, providing insights into immune cell interactions within the tumor microenvironment. Immune infiltration scores of six immune cells were calculated using Binding Association with Sorted Expression (BASE) [20,21,22,23,24], a rank-based gene set enrichment method. Previous publications have detailed and validated immune cell infiltration using this method [6,7,8,9,10,20]. BASE uses immune cell-specific weight profiles and patient gene expression data to infer immune cell infiltration for each patient and immune cell type. Briefly, BASE orders genes for a patient's gene expression profile from high to low expression and then uses weights from each immune cell weight profile to weigh the patient's gene expression values. BASE calculates two running sums, one representing the cumulative distribution of the patient's weighted gene expression values (foreground function) and another representing the cumulative distribution of the patient's complementary weighted (1-weight) gene expression values (background function). In the presence of a high amount of infiltrate from a specific immune cell type, the foreground function increases quickly, as the highly expressed genes in a patient's profile tend to be the ones with high immune cell weights, while the background function increases slowly. The maximal absolute difference between the foreground and background functions represents the immune infiltration level and, after a normalization procedure, results in the final immune infiltration score. Full details on the calculation and validation of the immune infiltration scores can be found in [7,25]. Similarly, BASE was used to calculate single cell-based TAM scores using TAM signatures (see next section).

2.4. Generation of TAM Signatures

RCC single-cell RNA-seq datasets from human RCC were obtained from previous publications [26,27]. Cluster annotations were also obtained from these publications. For each RCC cluster, a list of marker genes was provided by identifying genes that are over-expressed in the corresponding cluster than all the other clusters. These cluster-specific marker gene sets were used as TAM signatures. In total, eight human TAM signatures were defined (Supplementary Table S2). Given an RCC gene expression dataset, the BASE algorithm was used to calculate sample-specific TAM scores for each signature. The TAM signatures were represented as gene sets without assigning weights to the member genes. In this case, the BASE algorithm degenerated into a method like the single-sample GSEA analysis [28]. A high TAM score indicates that the corresponding TAM cells are strongly infiltrated in the tumor.

2.5. Lasso Cox-Regression

The TCGA-KIRC dataset was randomly divided into a training and testing set with a 1:1 ratio. The training set was analyzed to identify potential prognostic genes, and both the testing set, and the entire set were used for validation. First, univariate Cox-proportional hazards regression analysis was used to evaluate the association between the expression of IRGs and overall survival. Genes with a P value of < 0.05 based on the log-rank test were selected as candidate genes. Second, Least Absolute Shrinkage and Selection Operator (Lasso) Cox-regression analysis from the R glmnet package was employed to screen the IRGs most associated with overall survival in a multivariate model with 10-fold cross-validation, which resulted in Supplementary Table S3. These 27 genes composed the final risk score, which is described as follows:

$$\mathrm{Riskscore}=\sum _{\mathrm{i}=0}^{\mathrm{n}}{\mathrm{\beta }}_{\mathrm{i}}{\mathrm{x}}_{\mathrm{i}}$$

where ${\beta }_i$ refers to the coefficients of each gene and $x_i$ represents the expression value of the gene.

2.6. Survival Analysis

For univariate and multivariate survival analyses, Cox proportional hazards models were calculated using the “coxph” function from the R “survival” package. Survival curves were visualized using Kaplan-Meier curves using the “survfit” function from the R “survival” package. Median immune cell infiltration scores were used to stratify patients into “high” and “low” groups for univariate analyses. For multivariate analyses, an infiltration score of 0 was used as a separator to stratify patients into “high” and “low” groups. Differences in survival distributions in each Kaplan-Meier plot were calculated using a log-rank test using the “survdiff” function from the R “survival” package.

2.7. Statistical Analyses

The Spearman correlation coefficient (SCC) was reported for all correlation analyses as the assumptions underlying the Pearson correlation (i.e., normal distribution, homoscedasticity, or linearity) were not met. SCC was calculated using the R function cor, and significance was assessed using the cor.test. Principal component analysis (PCA) was performed using the prcomp R function. Principal component coordinates for each sample were extracted using the factoextra R package (https://github.com/kassambara/factoextra). Principal component 1 (PC1) was used to represent TAM infiltration. The sensitivity and specificity of the diagnostic and prognostic prediction models were analyzed by the ROC curve and quantified based on the area under the ROC curve (AUC). All statistical tests were two-sided; P values < 0.05 were considered statistically significant. All statistical analyses were performed using R software (version 4.2.0).

3. Results

3.1. Characterization of Tumor-Associated Macrophage Subtypes and Their Immune Microenvironment in Renal Cell Carcinoma

The TAM signatures in RCC samples exhibited substantial variability, reflecting the molecular heterogeneity within the tumor microenvironment (Figure 2A). Hierarchical clustering analysis revealed distinct expression patterns across samples, suggesting that TAM profiles may contribute to the stratification of tumor subtypes and provide insights into their functional roles.

The relationship among the eight TAM subtypes (TAM 1, 2, 3, 4, 5, 6, 7, and 8) was evaluated using pairwise Spearman correlation analysis. Strong inter-TAM correlations were observed, except for TAM 5 and 6, which exhibited weaker associations with the other subtypes (Supplemental Figure S1A). To further investigate this distinct behavior, correlation analysis was conducted between TAM signatures and immune cell fractions using the quanTIseq package in R (Figure 2B). This analysis revealed that while most TAM signatures exhibited strong correlations with macrophage populations, TAM 5 and 6 displayed lower correlation coefficients with classical macrophage markers and showed closer associations with neutrophils. These findings raised questions regarding their cellular identity and functional role within the tumor microenvironment. Except for TAM 5 and 6, positive correlations were identified, particularly with macrophages, NK cells, and CD8+ T cells, suggesting cooperative interactions among these subtypes in shaping the immune landscape of RCC.

Given the distinct behavior of TAM 5 and 6, further analysis was conducted by excluding these subtypes to focus on the primary TAM populations. The correlation analysis between the selected TAMs (TAM 1, 2, 3, 4, 7, and 8) was conducted, and a strong correlation was shown (Figure 2C). PCA was performed to capture overarching variation among the selected TAMs with PC1 explaining the majority of variance (90.8%), indicating a dominant pattern in TAM distribution (Figure 2D). Correlation analysis between PC1 and selected TAM subtypes confirmed that the selected TAMs were strongly associated with PC1 (Figure 2E). Further correlation between PC1 and immune cell populations derived from quanTIseq analysis indicated a strong association between the primary TAM axis and macrophage populations, as well as CD8+ T cells, highlighting the immunomodulatory role of TAMs in RCC (Figure 2F). The identity of these TAM signatures was further validated using an independent deconvolution method, BASE [20], as shown in Supplemental Figure S1B-C. We also performed the survival analyses for these TAM signatures (Supplemental Figure S1D-I). The P values are not significant in RCC cancers, which are consistent with the previous studies [29,30,31,32]. TAMs play a complex role in cancer progression, with their impact varying across different tumor types. While high TAM density is generally associated with poor prognosis in many cancers, including gastric, breast, bladder, ovarian, oral, and thyroid cancers [31,32] some studies report conflicting results. For instance, in colorectal cancer, high TAM density correlates with better overall survival [31]. The location of TAMs within the tumor microenvironment is crucial, as nest TAMs in gastric cancer are linked to increased tumor cell apoptosis and improved patient survival [29]. The prognostic value of TAMs may also depend on their phenotype (M1 or M2) and the biomarkers used for detection [30]. Overall, the impact of TAMs on patient survival is complex and tumor-specific, highlighting the need for further research to elucidate their role in different cancer types [30].

3.2. Principal Component Analysis of Tumor-Associated Macrophage Signatures and Their Immunological Associations in RCC

PC1 showed positive correlations with both M1 and M2 macrophage marker genes, indicating that the selected TAM signatures encompass features of both pro-inflammatory and immunosuppressive macrophage phenotypes (Figure 3A). Additionally, PC1 exhibited strong associations with immune checkpoint gene expression, suggesting a link between TAM abundance and regulatory immune signaling mechanisms (Figure 3B).

To assess the relationship between TAM-related variance and broader immune responses, PC1 was correlated with key immune-related pathways and scores. Significant correlations were observed with macrophage regulation, leukocyte fraction, stromal fraction, lymphocyte infiltration signature score, IFN-γ response, TCR richness, Transforming growth factor beta (TGF-β) response, and gamma delta T cells (γδ T cells). (Figure 3C–J). These findings suggest that TAM-associated variability is strongly linked to immune infiltration, immune regulation, and broader tumor microenvironment dynamics.

3.3. A 27-Gene Risk Score for Prognostic Prediction in RCC

The risk model analysis identified 27 risk-associated genes through LASSO Cox regression, utilizing 10-fold cross-validation to select the lambda value that minimized partial likelihood deviance. These genes exhibited non-zero coefficients, defining a gene set relevant for prognostic prediction (Figure 4A). Kaplan-Meier survival analysis stratified patients into low-risk and high-risk groups based on the risk scores. The survival curves displayed distinct differences in overall survival between the groups, with the high-risk group exhibiting shorter survival durations (Hazard Ratio = 4.64, p < 0.001) (Figure 4B). Time-dependent ROC curves were generated to assess the model's predictive accuracy for 1-year, 3-year, and 8-year survival, and AUC values of 0.82, 0.79, and 0.82 were observed, respectively, demonstrating consistent predictive performance across time points (Figure 4C). A forest plot illustrated univariate Cox regression results, including the risk score and clinical parameters such as age, gender, tumor stage, and pathologic T stage (Figure 4D). The correlation analysis between the risk scores and immune cell populations, derived from quanTIseq deconvolution, revealed notable patterns (Figure 4E). Higher risk scores were positively associated with M1 macrophages, CD8+ T cells, and Tregs indicating their enhanced prevalence in patients with elevated risk scores. In contrast, neutrophils and CD4+ T cells showed negative correlations, but these associations were not statistically significant, indicating that immune activation rather than suppression predominantly characterizes the high-risk group's tumor microenvironment. The relationship between the risk scores and immune pathways was also explored, highlighting key associations with pathways such as Cytotoxic Score, Proliferation, Stromal Fraction, Lymphocyte Infiltration Signature, and Leukocyte Fraction. These pathways showed strong positive correlations with the risk scores, suggesting an enrichment of immune infiltration and stromal activity in patients with higher risk scores (Figure 4F).

3.4. The TAM Risk Model Can Evaluate RCC Patients Across Different Clinicopathological Factors

Stratified survival analyses of the risk score model were conducted to evaluate its prognostic value across various clinicopathological factors. In male patients, the Kaplan-Meier survival analysis demonstrated a clear stratification of survival outcomes between the low-risk and high-risk groups, with the high-risk group exhibiting significantly shorter survival times (Figure 5A). Similarly, in female patients, the risk model effectively stratified survival outcomes, further validating its predictive capability across genders (Figure 5B). Age-stratified analyses revealed that the risk model consistently distinguished survival outcomes in both the elderly population (age > 50; Figure 5C) and younger patients (age ≤ 50; Figure 5D). High-risk patients in both age groups displayed poorer survival compared to their low-risk counterparts. The prognostic utility of the risk score model was also examined in relation to tumor stage. Patients with high tumor stages (III and IV) were distinctly stratified into high- and low-risk groups, with significant survival differences (Figure 5E). Similarly, for patients with low tumor stages (I and II), the risk model effectively categorized survival outcomes (Figure 5F). To further dissect the risk model's predictive ability, the TNM staging system was analyzed. In high T stage (T3–T4) patients, the risk model demonstrated significant stratification of survival outcomes (Figure 5G), while low T stage (T1–T2) patients exhibited comparable results (Figure 5H). For lymph node involvement, the model did not significantly stratify survival in patients with high N stage (N1–N2), potentially due to the limited sample size (n = 16) in this subgroup (Figure 5I). However, effective survival stratification was observed in patients with low N stage (N0) (Figure 5J). Finally, the risk model demonstrated similar predictive power for distant metastasis, with high M stage (M1) patients showing marked differences in survival outcomes between the risk groups (Figure 5K), while low M stage (M0) patients followed the same trend (Figure 5L).

3.5. High-Risk Patients with Significantly Down-Regulated TAM and Poor Prognosis in RCC

The distribution of risk scores among RCC patients demonstrated a clear stratification between low-risk and high-risk groups (Figure 6A). Patients were ranked by increasing risk scores, with a noticeable separation at the median threshold. Those classified as high-risk exhibited markedly elevated scores compared to their low-risk counterparts. Analysis of survival times revealed distinct differences between patients based on their risk classification (Figure 6B). A scatter plot of survival time against risk scores highlighted a concentration of shorter survival durations among high-risk individuals. Patients classified as "dead" predominantly clustered in the high-risk group, while those still alive were more frequent in the low-risk category. The vertical median line further emphasized the survival disparities between these groups. To explore the relationship between risk scores, PC1 groupings, and survival outcomes, an alluvial diagram was constructed (Figure 6C). This visualization illustrated the transitions between low- and high-risk groups, their associated PC1 classifications, and survival statuses. The diagram demonstrates that both high and low PC1 groups contribute similarly to poor survival outcomes ('dead' status). However, the high-risk score category, encompassing patients from both PC1 groups, distinctly shows a greater contribution to unfavorable clinical outcomes, showing the risk score's stronger association with survival compared to PC1 alone. The difference between high-risk and low-risk patients in their TIME, was further examined using the ESTIMATE algorithm, in Immune, Stromal, and ESTIMATE scores, respectively (Figures 6D–6F). High-risk patients displayed significantly higher immune scores, indicative of increased immune cell infiltration within the tumor microenvironment (Figure 6D). In contrast, the analysis of Stromal scores did not reveal a statistically significant difference between the high-risk and low-risk groups, suggesting comparable stromal content across both groups (Figure 6E). The ESTIMATE scores, which combine Immune and Stromal scores, were significantly elevated in the high-risk group, reflecting an overall enriched immune and stromal tumor microenvironment in these patients (Figure 6F). These findings align with the notion that high-risk patients exhibit a more active and complex tumor microenvironment.

4. Discussion

This study provides a refined characterization of TAMs in RCC and their relationship with the tumor immune microenvironment. A key finding of this study is the heterogeneity among the eight TAM subtypes, particularly the distinction of TAM 5 and 6, which exhibited weaker correlations with the remaining TAM signatures (Supplemental Figure S1A). A number of recent immunoprofiling studies corroborate that not all myeloid cells classified as “TAMs” share lineage markers, implying that certain subsets may indeed be more closely aligned with neutrophils or transitional myeloid phenotypes [33,34]. Further analysis revealed that these subtypes were less associated with macrophage markers and more closely related to neutrophils (Figure 2B), suggesting that they may represent a different myeloid lineage or a transitional phenotype rather than classical TAMs. These data reinforce emerging evidence that tumor-infiltrating myeloid cells are far from uniform, complicating attempts at a simple M1/M2 dichotomy [35,36]. This distinction among highly heterogeneous TAM subtypes challenges the conventional classification of TAMs as a relatively uniform population and indicates the need for a more nuanced categorization when assessing their functional roles in tumor progression.

To avoid confounding effects from these distinct TAM subtypes, we performed PCA using selected TAMs and identified PC1 as a dominant feature explaining the majority of variance (Figure 2D). Interestingly, PC1 correlated with both M1 and M2 macrophage markers (Figure 3A), reinforcing the concept that TAMs in RCC do not conform to a strict M1/M2 dichotomy but rather exist along a spectrum influenced by the tumor microenvironment. Numerous studies have shown that the microenvironmental cues, such as cytokines and metabolic byproducts, drive a continuum of macrophage phenotypes with overlapping M1- and M2-like features [37]. Additionally, PC1 showed strong associations with immune checkpoint expression (Figure 3B) and key immune-related scores and pathways (Figure 3C–J), further emphasizing that TAMs play a crucial role in shaping immune regulation within the tumor microenvironment. This aligns with emerging data demonstrating that macrophage-driven checkpoint ligand expression can potentiate T-cell exhaustion and immune evasion in several cancers, including RCC [38,39].

One of the most significant findings of this study is the correlative relationship between TAM activity, immune infiltration, and patient prognosis. The risk model developed in this study stratified patients into high-risk and low-risk groups based on a 27-gene signature, demonstrating a clear survival difference (Figure 4B). Notably, higher risk scores were positively correlated with CD8+ T cells, Tregs, and M1 macrophages (Figure 4E), suggesting that an immune-inflamed tumor environment does not necessarily confer a survival advantage [40,41,42,43]. Recent analyses of T-cell–inflamed tumors show that infiltration alone may be insufficient for clinical benefit, especially where the immunosuppressive axes are dominant [40,41]. This contradicts the conventional assumption that greater immune infiltration predicts better treatment outcomes. Indeed, studies have demonstrated that even large populations of CD8+ T cells or M1-like macrophages can fail to control tumor progression if these populations are functionally exhausted or co-opted by suppressive signals, resulting in a net immunosuppressive microenvironment despite high immune infiltration [41,44]. The strong correlation between high-risk scores and immune-related pathways, including cytotoxic activity and lymphocyte infiltration (Figure 4F), reinforces this notion and highlights the complexity of immune interactions in RCC.

The alignment between risk scores, TAM-associated PC1, and survival outcomes further shows the prognostic significance of TAMs in RCC (Figure 6C). High-risk patients with lower PC1 scores were predominantly associated with poor survival, indicating that TAM activity, rather than overall immune infiltration, is a key determinant of tumor progression. Studies using multiplex immunohistochemistry and single-cell transcriptomics have similarly observed that increased macrophage functional states (e.g., high antigen-presenting but also immunosuppressive signatures) portend worse outcomes [45,46]. Moreover, immune score analysis revealed that high-risk patients exhibited significantly higher immune infiltration (Figure 6D), yet stromal scores did not differ significantly (Figure 6E), suggesting that immune cell composition, rather than stromal content, plays a more decisive role in shaping patient outcomes.

These findings contribute to the growing recognition that TAMs are central players in the immune landscape of RCC. Unlike previous studies that broadly categorized TAMs, our approach provides a more refined classification, emphasizing the need to account for their heterogeneity in both research and clinical applications. Furthermore, our findings suggest that targeting TAMs may be a viable strategy for modulating the immune microenvironment in RCC. Future studies incorporating single-cell RNA sequencing could further delineate TAM plasticity and functional states, while experimental validation could clarify their precise roles in tumor progression. Prospective clinical validation of the risk model in independent cohorts would also strengthen its utility for guiding patient stratification and treatment decisions.

5. Conclusions

This study provides a comprehensive characterization of tumor-associated macrophages in RCC, revealing their heterogeneity, immune associations, and prognostic significance. By integrating transcriptomic deconvolution, PCA, and risk modeling, we identified distinct TAM subtypes and demonstrated that their activity, rather than overall immune infiltration, plays a dominant role in shaping patient outcomes. These findings emphasize the need to consider TAM diversity when evaluating immune regulation in RCC and highlight the potential of TAMs as therapeutic targets.

Despite these insights, several limitations should be acknowledged. First, the transcriptomic deconvolution methods used to infer TAM abundance rely on bulk RNA sequencing data, which may obscure intratumoral heterogeneity. Single-cell RNA sequencing or spatial transcriptomics could provide a more granular understanding of TAM subsets and their functional states. Second, while the risk model demonstrated strong prognostic value, further validation in independent RCC cohorts is necessary to establish its clinical applicability. Additionally, experimental studies are needed to confirm the biological mechanisms underlying TAM-driven immune modulation and their influence on tumor progression.

Future research should focus on refining TAM classifications through single-cell approaches, exploring functional validation using in vitro and in vivo models, and assessing the therapeutic potential of targeting TAMs to improve RCC treatment strategies. Integrating these findings into personalized medicine approaches may enable more precise risk stratification and immunotherapeutic interventions tailored to TAM-mediated immune dynamics in RCC.