Gene expression profiles of breast tumors and monocytes were obtained from MIAME-compliant studies (23) that followed Institutional Review Board (IRB)-approved protocols. The assembly, normalization and clinical annotation of the MDACC-701 gene expression dataset has been described previously (15). Briefly, this data set consists of Affymetrix U133A microarray gene expression profiles of 701 breast tumor biopsies (acquired prior to administration of neoadjuvant chemotherapy) and corresponding clinical drug response annotation. Compilation of normalized gene expression data and clinical annotations constituting the breast tumor meta-cohort #1 (MC1) data set has been described previously (8, 20). Briefly, the MC1 data set consists of 1,954 patient samples and corresponding clinical annotations. Expression data were generated using the U133 series GeneChip microarray platform (Affymetrix) containing 22,268 probe sets common to U133A, U133A2 and U133 PLUS 2.0 array platforms. MC1 represents independent cohorts from 16 medical centers in the United States, Europe and Asia. The meta-cohort #2 (MC2) data set includes 616 breast tumor expression profiles from patients in the United States, Asia, Europe and New Zealand. For both cohorts, raw data (CEL files) were downloaded from the NCBI’s Gene Expression Omnibus (GEO) data repository and processed as previously described (8, 20). MC1 and MC2 were analyzed on the same Affymetrix platform, and RMA-normalized and batch corrected using R software package from the Bioconductor project (8, 20). Molecular subtype calls were assigned to MC1 and MC2 cohorts according to Nagalla et al. (8). To identify gene expression changes downstream of TREM-1 activation, the Dower et al. gene expression dataset (24) comparing expression profiles of human peripheral blood monocytes with TREM-1 activation (induced by ligation with agonist antibody, n=11 patients) to those of matched control-treated monocytes (isotype antibody, n=11 patients) was analyzed. Affymetrix U133 PLUS 2.0 CEL files were downloaded from the GEO data repository (accession GSE9988) and RMA-normalized as described above. Differential expression analyses were conducted using Wilcoxon rank-sum test. Genes with Benjamini-Hochberg FDR-adjusted p-values < 0.01 were considered significant.

FID, WID and PID (Favorable, Weak and Poor Immune Dispositions, respectively) are gene expression-based, tumor immune subclasses that reflect the continuum from immunologically “inflamed” (FID) to immune “cold” (PID) tumors (8, 20). The immune subclasses are based on the integration of scores for three distinct effector immune gene signatures, identified in a previous breast tumor study using the MC1 data set (8). These immune signatures reflect the relative abundance of tumor-infiltrating effector cells associated with cytolytic activity (T/NK signature; CD8+ T cells/NK cells), antigen presentation (M/D signature; myeloid/dendritic cells), and humoral immunity (B/P signature; B cells/plasma cells). Each immune signature is scored in tumors using the geometric mean of expression of a set of immune-specialized genes, and each signature was found to provide additive prognostic information in multivariable Cox regression models (8). Patients were assigned to immune subclasses as follows. For each of the three immune signatures, tumors were ranked into population tertiles according to their immune signature scores. Tumors with upper-tertile expression scores for all three signatures, simultaneously, are classified as FID; tumors with lowest-tertile expression scores for any one signature are classified as PID; and tumors with combined intermediate- and upper-tertile signature scores are classified as WID, as established previously (8, 20). Accordingly, patients with FID tumors were observed to experience significantly better distant metastasis-free survival than those with WID or PID tumors, while those with PID tumors experienced significantly worse survival than those with WID or FID tumors. In a similar vein, FID tumors (at biopsy) were found to exhibit a significantly higher rate of responsiveness to neoadjuvant chemotherapy as compared to WID or PID tumors, consistent with reports that breast tumors with elevated TIL at baseline show greater responsiveness to neoadjuvant chemotherapy (25, 26); while PID tumors showed a significantly poorer response to neoadjuvant treatment as compared to WID or FID tumors (15).

Logistic regression was used to identify genes with antagonistic interactions with the effect of FID, WID and PID immune subclasses on tumor response to neoadjuvant chemotherapy. The MDACC-701 data set was utilized, and the analysis was restricted to the most well detected genes defined as being the probe sets whose mean expression levels across tumors were greater than the median of the mean expression among genes (n=11,108 probe sets). Tumor response was defined as favorable (Pathologic Complete Response (pCR) or Residual Cancer Burden (RCB) score = 0 or 1) or unfavorable (RCB score = 2 or 3) as previously described for this data set (15). In the logistic regression model, covariates included tumor response (0, 1), immune subclass (FID, WID, PID), patient age (continuous) and expression profile of each microarray gene (MG, dichotomous (using median as cut point)), while the interaction between MG expression and immune subclass was defined as being the parameter of interest. The magnitude of an interaction was represented by the regression coefficient of the product term. The logistic regression was run in R, and a Wald test was used to evaluate the statistical significance of each coefficient in the model. The Wald test computes a Z statistic as the ratio of the square of the regression coefficient to the square of the standard error of the coefficient and then computes significance. The unadjusted significance of discovered genes ranged from P = 0.05 to P = 0.00005, and gene ontology enrichment analysis was used to identify significantly-associated biological processes in the top genes ranked by Z statistic (Z values ≤ -2.0).

The DAVID Functional Annotation Tool (DAVID Bioinformatics Resources version 6.7, NIAID/NIH) (27) was used to study associations between gene sets and GO terms related to biological processes and pathways.

Molecular and clinical markers were analyzed as continuous or categorical variables for significant associations with DMFS by Cox proportional hazards regression (including multivariable models containing clinical and molecular covariates) and Kaplan-Meier analysis. GraphPad Prism, Sigmaplot 12.0 and the R package “survival” (v3.3-2) were used for these analyses, with DMFS endpoints defined as “no evidence of distant recurrence” (censored) or “clinical evidence of distant recurrence” (event), within 10 years of diagnosis.

To compute tissue-specific gene enrichment scores, we used the procedure and data set described in Benita et al., (28) and the limma package of Bioconductor (29). Briefly, one compares each group to all others and computes the linear model coefficient for each pair, which is a measure of the difference between two groups. For each linear model coefficient, limma also computes an associated P value (Bonferroni-corrected P ≤ 0.05). Genes that are specific to one tissue subset will result in higher enrichment scores. Proportions of immune cell types in tumors of the MC1 data set were estimated from RNA transcripts according to the CIBERSORT method and the LM-22 gene matrix (547 genes that distinguish 22 human hematopoietic cell types) (19, 30). CIBERSORT employs a novel support vector regression (SVR) to deconvolve proportions of distinct cell types found within complex mixtures. The resultant cell proportion estimates were mean-centered and averaged within quartiles of TREM1 expression. The R package gplots was used to visualize the cell proportion averages by heat map (31, 32). Cell proportion estimates were analyzed by Spearman’s rank-ordered correlation for significant associations with TREM1 expression. TREM1 expression in The Cancer Genome Atlas (TCGA)-breast cancer (BRCA) database was analyzed by TIMER (cistrome.org/TIMER/) (33) for significant correlations between TREM1, downstream target genes and markers of myeloid cell populations. The whole dataset (Breast Invasive Carcinoma; n = 1,093) and the Basal-like tumors (BRCA-Basal; n = 139) were utilized.

Human breast tumor tissue microarrays consisting of a variety of invasive ductal and lobular carcinomas were fabricated in-house by the Wake Forest Baptist Comprehensive Cancer Center’s Tumor Tissue and Pathology Shared Resource and used in immunofluorescence staining studies as follows. Formalin-fixed paraffin-embedded (FFPE) sections (5 μm) on slides were deparaffinized, antigen retrieved with sodium citrate method, blocked with 2.5% heat-inactivated horse serum, and incubated with polyclonal goat anti-TREM-1 antibody (ab) (1:50; R&D systems #AF1278) and mouse anti-CD68 (pan-myeloid) monoclonal ab (1:100; DAKO #2015-08, Clone KP1) in a humid chamber overnight at 4°C. After washing with PBS-T (PBS with 0.05% Tween-20), sections were incubated with rabbit anti-goat-Alexa 546 (1:200; Invitrogen) and goat anti-mouse Alexa 488 (1:250; Invitrogen) for 2 hours at room temperature in a humid chamber. After rinsing with PBS-T, sections were stained with DAPI or DRAQ5 (1:1000, 5μM, Cell Signaling) for 30 min, rinsed and mounted in Prolong anti-fade (Invitrogen) or FluoroShield (Sigma). Negative control and isotype matched mouse IgG1 (DAKO # 2015-07) were used to confirm staining specificity. Slides were imaged on an Olympus VS-110 Virtual Microscopy System at the Wake Forest School of Medicine Virtual Microscopy Core Facility.

Remnant surgical tissue from an AJCC stage III triple negative invasive ductal carcinoma resected at Wake Forest Baptist Medical Center (WFBMC) in Winston Salem, NC, was accessed via Wake Forest University Health Sciences Institutional Review Board protocol IRB00048977. Fresh tumor tissue was minced by razor blade into ~1 mm pieces and digested for 1 hr at 37°C with Human Dissociation Kit (Miltenyi Biotec, Auburn, CA) as per manufacturer’s instruction with the aid of gentleMACS Octo Dissociator fitted with heating collar (Miltenyi Biotec). Digested tissue was passed through a 70 µm strainer, washed and resuspended in RPMI plus 10% fetal bovine serum (FBS), then underlaid with Histopaque-1077 (Sigma-Aldrich) and centrifuged at 400 x g for 30 min. The resulting buffy coat was harvested and washed in RPMI plus 10% FBS, resuspended in the same media, and counted. Approximately 1 million cells were frozen in 90% FBS/10% DMSO. A single-cell cDNA library was generated from thawed cells exhibiting >80% viability using the 10X Genomics Chromium Controller (34) maintained by the Wake Forest Baptist Comprehensive Cancer Center’s Cancer Genomics Shared Resource. Indexed libraries were paired-end sequenced using an Illumina NextSeq 500 at an average read depth of 173,842 mapped reads per cell (n = 770 cells). Raw bcl and fastq data were demultiplexed, mapped to the GRCh38 human reference genome, and post-processed using Cell Ranger (10X Genomics) mkfastq pipelines and QC algorithms. Processed data files have been deposited in the Gene Expression Omnibus (GEO; www.ncbi.nlm. nih.gov/geo/) and can be accessed through GEO Series accession number GSE188600. Cells that expressed less than 400 genes or that contained mitochondrial reads comprising >20% of their total count of unique molecular identifiers were excluded. T-distributed stochastic neighbor embedding (t-SNE) (35) was used to reduce data dimensionality and cluster cells based on gene expression patterns. In accordance with the Cell Ranger R software, the statistical significance of genes differentially expressed between cell clusters or identified cell populations was assessed by negative binomial exact test [for genes with low counts, based on sSeq method (36)] and fast asymptotic beta test [for genes with high counts, edgeR method (37)]. Genes with Benjamini-Hochberg FDR-adjusted p-values < 0.01 were considered significant. A linear normalized expression value threshold > 0 was used to identify cells positive for TREM1 expression or expression of other marker genes studied.

In parallel, we analyzed a publicly available single cell RNAseq dataset (the Bassez et al. dataset) that profiled 175,942 cells from 31 treatment-naïve primary breast tumors (European Genome-phenome Archive (EGA) accession no. EGAS00001004809). Consistent with the analysis of the TNBC specimen from WFBMC, this study also utilized the 10X Genomics Chromium system for single cell library production, followed by Illumina sequencing (NextSeq 500 and NovaSeq 6000 platforms), and utilized Cell Ranger data processing workflows and similar quality control filters for cell selection (38). Analysis of the Bassez et al. dataset was facilitated by the BioTuring BBrowser and associated analytical tools (39) and the Seurat v3 R toolkit (40).

To identify genes with potential to antagonize anti-tumor immunity, we investigated breast tumor immune subclasses predictive of neoadjuvant chemotherapy response. Previous studies have shown that the abundance of tumor infiltrating lymphocytes in breast tumors is predictive of chemotherapy response in the neoadjuvant setting (41, 42). In earlier studies we identified three breast tumor immune gene signatures that reflect the relative abundance of tumor-infiltrating effector cell populations associated with cytolytic activity (T/NK signature), antigen presentation (M/D signature) and humoral immunity (B/P signature), each with independent prognostic power (8). These signatures allowed discernment of three tumor immune subclasses, termed FID, WID and PID (i.e., Favorable, Weak and Poor Immune Dispositions, respectively), that reflect the continuum from immunologically “inflamed” (FID) to immune “cold” (PID) tumors (see Materials and Methods). These subclasses were reported previously to differ significantly in terms of patient survival (8, 14, 20) and response to neoadjuvant chemotherapy (15) with FID tumors associating with more favorable clinical outcomes as compared to PID tumors, consistent with protective anti-tumor immunity being associated with FID status. To discover genes expressed in the TME that antagonize immune-mediated protective effects, we utilized a retrospective microarray dataset, MDACC-701 (15), consisting of 701 breast tumor expression profiles from core or needle biopsies taken from patients prior to neoadjuvant chemotherapy. In this data set, the FID immune subclass was found to identify patients with a significantly higher chemotherapy response rate as compared to WID or PID tumors. Using logistic regression (see Materials and Methods), we ranked genes from the MDACC-701 dataset based on the magnitude of antagonistic interaction, where higher gene expression in the FID immune subclass tumors was associated with a decrease in favorable chemotherapy response. Ranking genes by Z statistic, we identified 195 genes (207 probe sets) with interaction Z values ≤ -2.0 (see Supplemental Table 1 ), and analyzed these genes for significant enrichment of Gene Ontology (GO) terms. The one significant biological process identified was Neutrophil Chemotaxis (univariate P = 5.8x10^-5, Benjamini-Hochberg adjusted P = 0.06) which was supported by 7 genes. Among these genes, the most significant by Z statistic was TREM1 [ranked #6 overall for negative Z value (Z = -3.3, P = 0.001, Wald test)]. Moreover, included among the seven genes were the two well-characterized TREM-1 target genes, IL1B and CXCL8 (IL8), both of which showed a highly significant positive correlation with TREM1 expression in the MDACC-701 dataset (Pearson and Spearman correlations >0.6; see Supplemental Figure 1 ).

TREM-1 is a type-1 cell surface receptor expressed on cells of myeloid origin that triggers the transcriptional activation of numerous cytokines and chemokines in response to inflammatory signals, such as pathogen- and damage-associated molecular patterns (43). Moreover, TREM-1 was recently shown to inhibit anti-tumor immunity in a model of hepatocellular carcinoma (44). To further examine the relationship between TREM1 expression and the therapy-predictive power of the immune subclasses, we stratified patients first according to subclass (FID, WID or PID), and then according to low or high TREM1 expression ( Figure 1A ). We observed that the higher frequency of favorable response in FID tumors was dependent on low TREM1 expression. Specifically, the frequency of favorable response was 1.9-fold higher in TREM1-low FID tumors as compared to TREM1-high FID tumors (P = 0.005, Chi-square test), while no difference in response was observed among WID or PID tumors in the context of TREM1 low or high status ( Figure 1B ).

To examine the association between TREM1 expression and patient prognosis, we utilized a multi-institutional data set (Meta-Cohort 1 (MC1)) consisting of 1,954 breast tumor expression profiles characterized on the Affymetrix U133 platform and annotated for distant metastasis-free survival, molecular subtype, and other clinical and demographic features (8, 45). First, we subdivided the MC1 cohort into two evaluation cohorts (termed 977A and 977B, respectively) by randomizing the tumor expression profiles equally into two cohorts. Cox regression analysis was performed for each gene and within each evaluation cohort to identify genes with significant associations with DMFS. Genes significantly associated with DMFS (P < 0.1, Benjamini-Hochberg) in 977A (n = 3,094) and 977B (n = 3,304) were then used as input for hierarchical clustering of samples and genes ( Supplemental Figures 2A, B ). Previously reported prognostic gene clusters were identified, including a cluster of inferior outcome-associated, proliferation-related genes that mirror cancer proliferative capacity (8) and several correlated clusters of immune-specialized genes that reflect infiltrating B cells, T cells and dendritic cells and associate with good survival outcomes (8). TREM1 was found to be significantly inversely associated with DMFS in both 977A and 977B, and reproducibly clustered together with a set of co-expressed inferior outcome-associated genes enriched for myeloid cell functions including SPP1, IL8 (CXCL8), INHBA, GREM1, PLAU, PLAUR, AQP9, ADM, HPSE, and CTSL1. Among the genes of this cluster, TREM1 was found to have the most statistically significant association with inferior DMFS in both the 977A and 977B data sets (P = 5.3x10^-6 and P = 1.6x10^-4, respectively) ( Supplemental Figure 2C ). To further study this association, patients from the full MC1 cohort were stratified into quartiles based on TREM1 expression levels, and the quartiles were examined by Kaplan-Meier analysis. Patients whose tumors ranked in the highest quartile of TREM1 (Q4) experienced poorer DMFS as compared to those ranked in lower quartiles (P = 2.5x10^-11, log-rank test), with progressively worse DMFS observed for each increment of TREM1 quartile ( Figure 2A ). This observation was reproducible in an independent breast cancer data set, Meta-Cohort 2 (MC2), consisting of 616 patients whose tumors were also profiled on the Affymetrix U133 platform (8, 20) (P = 3.1x10^-03; Figure 2B ). Analysis of molecular subtypes showed that the highest quartile of TREM1 expression (Q4) contained proportionally more Basal-like, Claudin-low and HER2E tumors, and fewer Luminal A (LumA) and Normal-like tumors, as compared to the lowest quartile (Q1) ( Figures 2C, D ). However, TREM1 expression distributions only showed moderate variation between different molecular subtypes, with decreased expression of TREM1 observed in LumA and Normal-like tumors, in particular ( Figures 2E, F ). Notably, we observed that TREM1 expression associated with poorer DMFS within certain clinical and molecular subtypes. By univariate Cox regression analysis ( Table 1 ), TREM1 was associated with DMFS in both estrogen receptor positive (ER+) and negative (ER-) breast cancers, as well as lymph node positive (LN+) and negative (LN-) disease. Within molecular subtypes, TREM1 levels showed significant associations with DMFS in Basal-like, Luminal B (LumB) and Luminal A (LumA) tumor subtypes. We further examined the prognostic relevance of TREM1 by adjusting for conventional variables in multivariable Cox models ( Table 2 ). As a continuous variable, TREM1 expression was associated with poorer DMFS, univariately (P = 1.9x10^-13), as well as in the multivariable model (P = 3.2x10^-05), where it remained significant in the presence of molecular subtype, immune subclass, histologic grade, tumor size, lymph node status, ER status, patient age, and treatment status.

We hypothesized that the negative association between TREM1 expression and DMFS may reflect, in part, TREM-1-related immunoregulatory signaling that antagonizes favorable immune-mediated outcomes. In previous work, we showed that the FID, WID and PID immune subclasses exhibited a reproducible significant prognostic stratification in breast tumors classified as immune benefit-enabled (IBE), consisting mostly of highly proliferative Basal-like, HER2-enriched and Luminal B subtypes, but not in immune benefit-disabled (IBD) tumors defined by a reduced proliferative phenotype (20). Using the CIBERSORT algorithm (19) which deconvolves gene expression profiles from complex cellular mixtures and estimates the relative proportions of 22 different leukocyte types, we investigated the prognostic significance of TREM1 in the context of the estimated proportions of CD8+ T cells in IBD and IBE breast tumors ( Figure 3 ). CD8+ T cell scores were used to stratify the MC1 cohort into CD8+ T cell quartiles (Q1=lowest quartile, Q4=highest quartile). While higher quartiles of CD8+ T cells associated with improved DMFS in IBE tumors (P = 0.0002) ( Figures 3A, B ), the prognostic power of TREM1 was most significant in IBE CD8+ T cell-High tumors (defined as quartiles 3 and 4 combined; P = 0.0002) as compared to CD8+ T cell-Low tumors (defined as quartiles 1 and 2 combined; P = 0.58) ( Figures 3C, D ), with TREM1-Low tumors exhibiting a more favorable DMFS. This observation is consistent with a possible role for TREM-1 in attenuating protective T cell-mediated immunity in breast cancer.

To determine the cellular source of TREM1 expression, we first analyzed the expression patterns of TREM1 and the genes comprising the TREM1-associated gene cluster observed in breast tumors using the Gene Enrichment Profiler tool (28), which measures tissue-specific gene enrichment in 126 normal human cell subsets and tissues. This analysis showed that in non-malignant tissues, TREM1 is highly enriched in macrophages and neutrophils ( Supplemental Figure 3 ), consistent with previous reports (46, 47). Similarly, a number of genes comprising the TREM1-associated gene cluster that we identified by hierarchical clustering (ie, SPP1, PLAUR, AQP9, IL8 (CXCL8), ADM, UPP1, and INHBA) were also highly enriched in myeloid cells, and in lipopolysaccharide (LPS)-stimulated macrophages, in particular ( Supplemental Figure 3 ), suggesting a common myeloid origin for TREM1 and this gene expression cassette. In parallel, we investigated TREM-1 expression directly by immunofluorescence staining using a breast tumor tissue microarray consisting of >40 invasive ductal and lobular carcinomas positive or negative for the clinical markers ER, PGR and HER-2. Across tumor types, the predominant TREM-1 staining pattern observed was moderate to strong TREM-1 staining in a fraction of tumor-infiltrating myeloid cells that co-expressed CD68, with negative staining in tumor cells (pattern P-1, Figure 6A ). In some instances we observed only CD68-positive cells negative for TREM-1 staining, which coincided with malignant (and normal) epithelium also negative for TREM-1 (P2). In more rare instances we observed TREM-1 staining in non-myeloid cells of tumors (P3), such as moderate TREM-1 staining in cancer cells (n=2), or moderate to high TREM1 staining in cells morphologically consistent with cancer-associated fibroblasts (n=3). To further investigate the relationship between TREM1 and tumor infiltrating immune cells, we used the CIBERSORT algorithm to study the relationship between immune cell abundance estimates and TREM1 expression quartiles in the MC1 cohort. As shown in Figure 4 , a number of immune cell types showed marked proportional differences across TREM1 expression quartiles. In the highest TREM1 quartile (Q4), we observed significant enrichment of neutrophils, macrophages (M0 and M2, but not M1), mast cells (activated, but not resting) and dendritic cells (activated). By contrast, significant depletion in TREM1 Q4 was observed for CD8+ T cells, monocytes, follicular helper T cells, mast cells (resting), dendritic cells (resting) and plasma and memory B cells. These observations are consistent with a positive correlation between high TREM1 expression, increasing proportions of immunosuppressive myeloid cell types (M2 macrophages, neutrophils and mast cells) and decreasing proportions of anti-tumor effector cells (CD8+ T cells and follicular helper T cells).

Next, we examined the TCGA breast cancer (BRCA) RNAseq data set for TREM1 correlations with known target genes ( Figure 5A ) and myeloid markers ( Figure 5B ). Side-by-side analyses were conducted for the set of all breast tumors, and for the basal-like breast tumors, alone, since they showed the highest average expression of TREM-1 as compared to other subtypes in the MC1 and MC2 cohorts. Significant positive correlations were identified between intratumoral TREM1 expression and downstream TREM-1 target genes that encode cytokines with direct roles in the recruitment and activation of pro-tumorigenic myeloid cells (IL1B, IL8, IL6, MCP-1/CCL2), as well as markers of pro-tumorigenic myeloid cells including CD11B/ITGAM (a general MDSC marker), CD14 [marker for monocytic (M)-MDSCs (48)], OLR1 [marker for PMN-MDSC (49)] and CD206/MRC1 (M2 macrophage marker).

To more precisely characterize TREM1 expression in the breast tumor microenvironment, we profiled 770 cells from a fresh surgical TNBC specimen by single-cell (sc) RNAseq. Of 8 cell types identified by clustering ( Figures 6B, C ), prevalent and high-level TREM1 expression was observed in the myeloid cell population (n=301 cells), while low-level TREM1 was observed sparsely in tumor cells and T cells ( Figure 6D , left panel). Notably, TREM1 expression overlapped sharply with markers of M2 macrophages (MRC1/CD206) and polymorphonuclear (PMN)-MDSCs (OLR1) ( Figure 6D ). To substantiate these findings, we analyzed a large publicly available scRNAseq dataset that profiled 175,942 cells from 31 breast tumors biopsied prior to treatment (38). Consistent with the initial scRNAseq and TREM-1 IF staining analyses, the prevalence and level of TREM1 expression in the Bassez et al. dataset was greatest in the myeloid compartment, with sparse expression observed in other cell types, including cancer cells and fibroblasts, and to a lesser extent, B cells and endothelial cells ( Figures 7A, B ). Similar TREM1 expression patterns were observed among cell types grouped by tumor subtype (TNBC, n=13; ER+, n=15; HER2+ (ER+/-), n=3) ( Supplemental Figure 4 ), and within tumors of individual patients ( Supplemental Figure 5 ), confirming that TREM1 expression in breast cancer is a predominantly myeloid phenomenon. Within the myeloid compartment (n=16,485 cells), TREM1 was significantly overexpressed in cells positive for expression of OLR1 or MRC1 expression, as compared to cells negative for marker expression, respectively ( Figure 7C ; P<0.0001, both comparisons), consistent with our initial observations in the TNBC specimen ( Figure 6 ), as well as the positive correlations observed between TREM1 and OLR1 and MRC1 in the TCGA cohort ( Figure 5 ).

To gain insight into the transcriptional programming downstream of TREM-1, we performed differential expression analysis in the myeloid cell population differentiated by TREM1 expression (positive vs. negative). First, we observed a number of significantly differentially expressed genes that were reproducibly identified between the TNBC specimen dataset and that of Bassez et al. ( Figures 7D, E ). Cross-referencing these genes with genes previously identified by Dower and colleagues (24) as being significantly induced by TREM-1-specific activation in human monocytes isolated from peripheral blood ( Figure 7F ), resulted in a list of confirmed transcriptional targets robustly overexpressed in the TREM1+ cells of the breast cancer TME. Among these genes were known targets of TREM-1 signaling with pro-inflammatory functions (IL1B, CCL7, CXCL3) and known targets with anti-inflammatory or mixed pro- and anti-inflammatory functions (IL1RN, INHBA, IL6, IL8, CCL2). Also identified in the TREM1+ positive myeloid cells of the breast tumors were genes shown by Dower et al. to be uniquely upregulated in monocytes by TREM-1-specific activation (by agonist antibody ligation) but not by LPS stimulation, including MMP19, IL1RN, PLAUR, CCL7 and SPP1. Notably, SLC11A1 and SPP1 were among the top-most differentially expressed genes in all three datasets. Finally, in tumor myeloid cells negative for TREM1 expression, we observed consistently higher expression of MHC class II antigen-presenting genes (HLA-DRA, HLA-DRB1, HLA-DMB, HLA-DPA1, HLA-DQA1) and the Th1-associated chemokine, CXCL9. However, these genes were not observed to be down-regulated by TREM-1 activation in the Dower et al. study, suggesting they represent a phenotype of non-TREM1 expressing cells in the breast TME. Of note, genes encoding canonical immunosuppressive factors such as TGF-β (TGFB1), IL-10 (IL10) and PD-L1 (CD274) were not consistently differentially expressed in our analyses; although IL10 and CD274 were moderately associated with TREM1 expression ( Figure 7D ) or activation ( Figure 7F ), respectively, these observations lacked reproducibility across datasets.

Next, we assembled a TREM-1 target gene panel comprising genes significantly induced by TREM-1 activation (in Figure 7F ) and associated with intratumoral TREM1+ myeloid cells (in Figure 7E ) to serve as an estimator of TREM-1 signaling output in different cell and tumor populations.

Where TREM1 expression was significantly higher in OLR1+ or MRC1+ cells as compared to cells negative for OLR1 or MRC1 expression ( Figure 7C ), the TREM-1 target genes also exhibited a simultaneous increase in frequency and magnitude of expression in the OLR1+ or MRC1+ cell compartments ( Supplemental Figure 6 ), consistent with heightened TREM-1 signaling in these populations. We next used the target gene panel to profile individual tumors and major cell types expressing TREM1, including myeloid cells, cancer cells and fibroblasts ( Figure 7G ). By ranking tumor samples (columns) on TREM1 expression frequency in the myeloid compartment, a clear correlation was observed between the percentage of TREM1-expressing myeloid cells and the frequency and magnitude of expression of the target genes. Moreover, we observed that myeloid TREM1 frequency also associated with clinical breast cancer subtype (TNBC, ER+, HER2+), whereby an overall difference in rank of the percentage of TREM1-expressing cells was significant (P=0.03, three-group Kruskal-Wallis test). This difference in rank was also significant for the TNBC and ER+ subtypes alone (P=0.007, two-group comparison), where higher TREM1 frequency associated with TNBC status while lower TREM1 frequency associated with ER-positive tumor status. Notably, this observation was consistent with our findings in the MC1 and MC2 cohorts where TREM1 expression by bulk tumor analysis was higher in basal-like tumors as compared to luminal tumors ( Figures 2C–F ). Across individual tumors, we examined the frequencies and expression levels of TREM1 and its target genes in cancer cells and fibroblasts ( Figure 7G ) where sparse, yet moderate to high TREM1 expression could be detected in some cells ( Figure 7B ). In contrast to myeloid cells, neither the cancer cell nor fibroblast compartments showed high frequency TREM1 expression in individual tumors or coordinated expression of target genes. Similar results were observed for the B cell and endothelial compartments (data not shown). Together, these findings suggest that TREM-1 signaling is frequently operative in breast cancer, and that myeloid cells are the predominant source of TREM-1 signaling in the TME.

Finally, we investigated the relationships between TREM1 myeloid expression and the expression of other myeloid receptors related to TREM-1 activation and function. Historically, TREM-1 signaling has been widely studied in the context of lipopolysaccharide (LPS)-induced Toll-like receptor (TLR) signaling in innate immunity, where TLR activation via PAMP/DAMP-sensing can potently induce the expression and activation of TREM-1 in a pro-inflammatory amplification loop (21). In the Bassez et al. breast tumor myeloid cells, the expression of 5 TLR genes (TLR1, TLR2, TLR4, TLR7 and TLR8) were detected in >10% of cells and at moderate to high expression levels. Upon comparing TREM1 expression levels between myeloid cells positive or negative for TLR expression ( Supplemental Figure 7A ), no TLRs exhibited a substantial positive association with TREM1 expression, suggesting that TLRs are not major determinants of TREM1 expression in the breast TME. The TREM-2 receptor is a TREM family member with strong homology to TREM-1, known for its various roles in anti-inflammatory signaling (50) which may antagonize TREM-1 activity (51). Like TREM-1, TREM-2 also requires binding to the DAP12 adaptor protein to mediate signal transduction (52), thus potentially competing with TREM-1 in the presence of activating signals. TREM2 expression was observed in approximately one-third of breast tumor myeloid cells, but with no significant association with TREM1 expression ( Supplemental Figure 7B ). Moreover, in myeloid cells either positive or negative for TREM2, TREM1+ cells were consistently associated with greater frequency and magnitude of expression of the TREM-1 target genes ( Supplemental Figure 7C ) consistent with TREM-1 activation occurring independent of TREM-2 status. Together, these findings suggest that neither TLRs nor TREM-2 substantially influence TREM1 expression characteristics in the breast TME.