Diversity of intratumoral regulatory T cells in B-cell non-Hodgkin lymphoma

Abstract

Tumor-infiltrating regulatory T cells (Tregs) contribute to an immunosuppressive tumor microenvironment. Despite extensive studies, the prognostic impact of tumor-infiltrating Tregs in B-cell non-Hodgkin lymphomas (B-NHLs) remains unclear. Emerging studies suggest substantial heterogeneity in the phenotypes and suppressive capacities of Tregs, emphasizing the importance of understanding Treg diversity and the need for additional markers to identify highly suppressive Tregs. Here, we applied single-cell RNA sequencing and T-cell receptor sequencing combined with high-dimensional cytometry to decipher the heterogeneity of intratumoral Tregs in diffuse large B-cell lymphoma and follicular lymphoma (FL), compared with that in nonmalignant tonsillar tissue. We identified 3 distinct transcriptional states of Tregs: resting, activated, and unconventional LAG3+FOXP3- Tregs. Activated Tregs were enriched in B-NHL tumors, coexpressed several checkpoint receptors, and had stronger immunosuppressive activity compared with resting Tregs. In FL, activated Tregs were found in closer proximity to CD4+ and CD8+ T cells than other cell types. Furthermore, we used a computational approach to develop unique gene signature matrices, which were used to enumerate each Treg subset in cohorts with bulk gene expression data. In 2 independent FL cohorts, activated Tregs was the major subset, and high abundance was associated with adverse outcome. This study demonstrates that Tregs infiltrating B-NHL tumors are transcriptionally and functionally diverse. Highly immunosuppressive activated Tregs were enriched in tumor tissue but absent in the peripheral blood. Our data suggest that a deeper understanding of Treg heterogeneity in B-NHL could open new paths for rational drug design, facilitating selective targeting to improve antitumor immunity.

Graphical abstract

Figure 1.

Transcriptional landscape of CD4⁺T cells in B-NHL and healthy donor tonsils. (A) Schematic representation of the workflow of scRNA-seq, TCR single-cell variable diversity joining sequencing (scVDJ-seq), and CITE-seq of sorted CD4⁺ T cells. (B) Single-cell data derived from sorted CD4⁺ T-cell populations from healthy donor tonsils (n = 3), FL (n = 3), and DLBCL (n = 3) projected onto uniform manifold approximation and projection (UMAP) by combining scRNA-seq and CITE-seq using weighted nearest neighbor (Wnn) method, provided 13 distinct clusters based on gene and protein expression differences for 18 771 cells passing the quality control. The clusters were annotated based on a combination of gene and protein expression as naive (KLF2, CCR7, LEF1, and protein expression of CD45RA), 3 memory clusters (CD69, SELL, and lack of CD45RA protein expression), GZM⁺ CD4⁺ T cells (GZMK/A, NKG7, and CST7), 2 T follicular helper (Tfh) clusters (PDCD1, CXCR5, IL21, TOX2, and protein expression of PD-1, CXCR5), a T helper cell cluster (CXCR4 and KLF6), 2 clusters of Tregs (FOXP3, IL2RA, CTLA4, and protein expression of CD25 but lack of CD127), and 1 cluster of unconventional FOXP3⁻LAG3⁺ Tregs (LAG3, CTLA4, IL10, and lack of CD127 protein). (C) UMAPs of the single-cell data shown in panel B, divided based on tissue origin and color coded based on patient sample. (D) Dot plots showing average expression of top 5 differentially expressed genes (DEGs) for each cluster. (E) Expression of selected genes (purple) and surface protein expression (green) overlaid onto the Wnn UMAP coordinates from panel B. (F) Pie charts showing the cellular abundance of each cell cluster per patient sample.

Figure 2.

Characterization of TCR diversity and clonotype expansion by scVDJ-seq. (A) Clonal expansion mapped onto the UMAP from Figure 1B. (B) Shannon diversity index estimated from TCR sequencing results for 45 random cells per cluster. Values are shown for each cluster identified in Figure 1B. (C) Histograms showing the percentages of Treg subsets with expanded clonotypes within the entire Treg compartment and across samples of tonsils (n = 3), FL (n = 3), and DLBCL (n = 3). Expanded clonotypes were categorized as small expansion (1 < X ≤ 3) and medium expansion (3 < X ≤ 6). (D) Distribution of CDR3 length for each Treg cluster as per the tissue origin. The mean CDR3 lengths for activated Tregs (actTregs), resting Tregs (restTregs), and LAG3⁺ Tregs were 28.1 amino acids (aas), 28.3 aas, and 27.9 aas, respectively in FL; 28.0 aas, 27.7 aas, and 28.3 aas in DLBCL; and 28.3 aas, 28.0 aas, and 28.1 aas in tonsils. (E) Venn diagrams showing the number of clonotypes within and shared between Treg subsets. For scaling purposes, numbers are shown for the Tregs found among 1000 random CD4⁺ T cells with gene expression and single TCRα/TCRβ gene.

Figure 3.

Characterization of Treg subsets by scRNA-seq. Volcano plots based on the scRNA-seq data from healthy donor tonsils (n = 3), FL (n = 3), and DLBCL (n = 3), showing key significantly DEGs between (A) activated and resting Tregs and (B) activated and LAG3⁺ Tregs. Each dot indicates 1 gene; red dots represent a log₂ fold change (log₂FC) of >0.35, and blue dots depict a log₂FC < −0.35. Significant genes were selected with adjusted P value < .05. Full list of DEGs is shown in supplemental Tables 6 and 7. (C) Dot plot of selected genes from scRNA-seq data, showing average gene expression of checkpoint receptors, transcription factors, cytokines, and DEGs across the 3 distinct cellular subsets of Tregs. (D) Heat map of gene set enrichment analysis (GSEA) highlighting highly enriched pathways among the 3 Treg clusters based of Hallmark and KEGG (K) libraries, shown as enrichment score.

Figure 4.

Intratumoral Tregs are phenotypically and functionally heterogeneous. (A) The expression of 18 checkpoint receptors on CD4⁺FOXP3⁺CD25⁺CD127⁻ Tregs from healthy donor PBMCs (n = 3) and tonsils (n = 3), FL (n = 3), and DLBCL (n = 2) were assessed by t-distributed stochastic neighbor embedding (tSNE) analysis and FlowSOM clustering of mass cytometry data (supplemental Table 1). (Top) Population analysis of Tregs colored based on FlowSOM metaclusters. (Bottom) Density plots of checkpoint receptor protein expression on concatenated Tregs. (B) Heat map of checkpoint receptor expression in activated Tregs (metacluster 2) and resting Tregs (metacluster 1) from the same samples as in panel A. The heat map shows arcsinh ratio of median fluorescence intensity (MFI) normalized to column minimum. (C) Spearman correlation of Treg subset distribution assessed by flow cytometry (5-marker tSNE) with the distribution assessed by mass cytometry (18-marker tSNE). (D) Frequency of activated Tregs among CD4⁺ T cells assessed by flow cytometry in PBMCs (n = 5) and tonsils (n = 7) from healthy donors, and in FL (n = 15), DLBCL (n = 16), and MCL (n = 10). Horizontal line represents median frequency. (E) Frequency of the activated Treg subset among CD4⁺ T cells in paired samples, PBMCs, and LNs, from patients with FL (n = 9). (F) Representative flow cytometry histogram of CellTrace Violet dilution in effector CD4⁺ or CD8⁺ T cells when stimulated with T-cell-expansion beads (stim. Teff) and cocultured in presence of either actTregs or restTregs purified from an FL sample (left), or calculated as percent suppression of stim. Teff proliferation when cocultured in the presence vs absence of Tregs, (n = 9). (G) (Left) A representative experiment showing TGF-β expression for activated and resting Tregs in the presence or absence of different stimuli. (Right) Summary of TGF-β median fluorescence intensity for unstimulated, activated, and resting Tregs from FL (n = 4) and tonsillar (n = 5) samples. Statistical differences were calculated using Kruskal-Wallis with Dunn multiple comparison test for panel D and Wilcoxon for panels E-G. ∗P < .05; ∗∗P < .01; ∗∗∗P < .001.

Figure 5.

Spatial single-cell analysis of FL tissues. Imaging mass cytometry was performed on tissue section from formalin-fixed paraffin-embedded FL samples (n = 4). (A) Shown is 1 representative image with FOXP3 (green, Tregs), CD4 (pink), CD20 (blue; tumor cells), and CD21 (yellow; follicular dendritic cells) delimiting the follicles and T-cell zone. Scale bar, 100 μm. (B) Scatter plots depicting the median distance of different cell populations to the nearest actTreg, localized either in the T-cell–rich zone or in the intrafollicular area. Statistical differences were calculated using Kruskal-Wallis test and Dunn multiple comparisons (supplemental Table 8). ∗∗∗P < .001 and ∗∗∗∗P < .0001. (C) actTregs can colocalize with CD8 T cells. Shown is a selected area from the T-cell zone of panel A, with FOXP3 (green), TIGIT (pink), and CD8 (blue). Arrow heads point at FOXP3⁺TIGIT⁺ actTregs in proximity to CD8⁺ T cells. (D) actTregs can colocalize with dendritic cells (DCs) and macrophages. Shown is 1 representative area from the image in panel A, with FOXP3 (green), CD11c (pink, DCs), CD68 (blue, macrophages), and CD21 (white). Arrow heads point at Tregs in proximity to DCs, whereas arrows point at Tregs in proximity to macrophages.

Figure 6.

Validation of Treg clusters in separate scRNA-seq cohorts and prognostic significance by imputation of Treg subset frequencies in bulk RNA-seq cohorts. (A) Heat map of top 30 DEGs across the 3 Tregs subsets; actTregs, restTregs, and LAG3⁺ Tregs were plotted after normalization using scTransform. (B) Unique gene signature matrices were developed for CD4⁺ T-cell clusters including actTregs, restTregs, and LAG3⁺ Tregs, using the computational framework of CIBERSORTx. The matrices were then used to reannotate cell clusters in our scRNA-seq data sets and are shown as histograms of the frequency of Treg subsets in healthy donor tonsils (n = 3), FL (n = 3), and DLBCL (n = 3) (left); and in external scRNA-seq data sets (right), healthy donor tonsils (n = 8), FL (n = 7), and DLBCL (n = 6). (C) Pearson correlation of Treg subsets distribution in the external scRNA-seq data sets with the distribution in the scRNA-seq data generated in this study. (D) Abundance of the Treg subsets (left) and Kaplan-Meier survival curve for progression-free survival of patients with FL in the bulk RNA-seq cohort, stratified above the 85th percentile for the actTreg population (right). (E) Abundance of the Treg subsets (left) and Kaplan-Meier survival curve for failure-free survival of patients with FL in the external bulk RNA-seq cohort (right), in which patients were stratified using the same actTreg abundance threshold as in the discovery cohort (shown in panel D). ∗∗∗P < .001.