Chimeric antigen-receptor (CAR) T cells lead to high response rates in myeloma, but most patients experience recurrent disease. We combined several high-dimensional approaches to study tumor/immune cells in the tumor microenvironment (TME) of myeloma patients pre– and post–B-cell maturation antigen (BCMA)-specific CAR T therapy. Lower diversity of pretherapy T-cell receptor (TCR) repertoire, presence of hyperexpanded clones with exhaustion phenotype, and BAFF+PD-L1+ myeloid cells in the marrow correlated with shorter progression-free survival (PFS) following CAR T therapy. In contrast, longer PFS was associated with an increased proportion of CLEC9A+ dendritic cells (DC), CD27+TCF1+ T cells with diverse T-cell receptors, and emergence of T cells expressing marrow-residence genes. Residual tumor cells at initial response express stemlike genes, and tumor recurrence was associated with the emergence of new dominant clones. These data illustrate a dynamic interplay between endogenous T, CAR T, myeloid/DC, and tumor compartments that affects the durability of response following CAR T therapy in myeloma.

Significance:

There is an unmet need to identify determinants of durable responses following BCMA CAR T therapy of myeloma. High-dimensional analysis of the TME was performed to identify features of immune and tumor cells that correlate with survival and suggest several strategies to improve outcomes following CAR T therapy.

See related commentary by Graham and Maus, p. 478.

This article is highlighted in the In This Issue feature, p. 476

T-cell redirection with chimeric antigen-receptor (CAR) T cells or bispecific antibodies has led to high rates of clinical response in patients with relapsed/refractory multiple myeloma (RRMM; ref. 1). CAR T cells targeting B-cell maturation antigen (BCMA) have achieved regulatory approval by the U.S. Food and Drug Administration as a therapeutic option in multiple myeloma (MM) patients progressing after four or more prior lines of therapy (2). However, in contrast to the experience in patients with acute lymphoblastic leukemia or B-cell non-Hodgkin lymphoma, CAR T therapy in MM has not led to durable unmaintained responses. Instead, MM patients treated with BCMA-specific CAR T cells remain at ongoing risk of disease recurrence even after complete response (2, 3). Therefore, there is an unmet need to identify host, tumor, or therapy-related factors that affect the risk of recurrent disease following current MM CAR T therapies.

Prior studies have shown that the expansion of CAR T cells in vivo following infusion correlates with the initial response to therapy, and the depth of response correlates with the outcome (2, 4). Cellular pharmacokinetic studies have shown variable but limited durability of CAR T persistence in MM (2). Although BCMA downregulation has been observed following BCMA CAR T infusion (4), and biallelic deletion of BCMA was linked to resistant clones in some patients (5), most patients recur with BCMA-expressing tumors. Properties of the tumor microenvironment (TME) have been linked to response and resistance to several immune therapies, particularly immune-checkpoint blockade in solid tumors (6). A recent study has shown that the TME in MM is characterized by the progressive emergence of exhausted T cells and attrition of TCF1+ stemlike memory T cells with proliferative potential (7). The impact of the TME on outcomes following BCMA CAR T therapy in MM is not known. In this study, we combined several high-dimensional single-cell approaches, including cellular indexing of transcriptomes and epitopes by sequencing (CITE-seq), single-cell transcriptomics, mass cytometry, and T-cell-receptor sequencing, to evaluate bone marrow cellular phenotypes and states that correlate with the durability of response following BCMA CAR T infusion in a prior clinical trial (4).

Pre- and posttherapy bone marrow (BM) specimens from MM patients enrolled in a clinical trial of BCMA-specific CAR T therapy (4) exhibiting a clinical response to therapy were analyzed with a combination of high-dimensional single-cell approaches, including CITE-seq, single-cell transcriptomics, mass cytometry, and T-cell-receptor sequencing. A total of 28 specimens were analyzed, of which 23 underwent CITE-seq and single-cell transcriptomics and 24 underwent mass cytometry (Fig. 1A; Supplementary Fig. S1A and S1B). These data were then correlated with the durability of clinical response. Clinical aspects of the trial have been described previously (4). The 23 specimens that underwent CITE-seq and single-cell transcriptomics yielded 151,054 single cells meeting predefined quality control criteria (Supplementary Fig. S2). Unsupervised clustering of single-cell transcriptomic data identified 63 clusters, which were classified into broad categories of tumor, T, natural killer (NK), B, myeloid/dendritic cell (DC) and progenitors based on the expression of lineage markers detected by concurrent antibody staining (Fig. 1B and C; Supplementary Fig. S2A and S2B).

In order to correlate immune data with outcome, patients were split into two broad groups, based on PFS <6 months (median 125 days, range, 57–150 days, n = 7) or >6 months (median 752 days, range, 190–1,053 days, n = 4; Supplementary Table S1). CAR T-cell vector copy number in peripheral blood at peak expansion and at day 28 (D28) in bone marrow did not differ significantly between the two groups (Supplementary Table S1). As expected, tumor cells as analyzed by CITE-seq declined at D28 after CAR T-cell therapy in all patients when compared with pretreatment samples (mean pretherapy 35% vs. D28 0.5% P = 0.03; Supplementary Fig. S3A). Although there were no significant differences in the proportion of tumor cells when comparing samples from patients with long or short PFS at pretreatment or day 28 time points, the proportions of tumor cells at 3 months following infusion were higher as expected in patients with shorter PFS (PFS <180 24% vs. 0.6% for PFS >180 P = 0.03; Supplementary Fig. S3B).

In order to account for the differences in the proportion of tumor cells at different time points, we focused on the nontumor compartment to analyze changes in BM immune composition. Patients with longer PFS had a higher proportion of T cells at D28 (21% vs. 62% P = 0.06), which was not observed in patients with shorter PFS (Supplementary Fig. S3C). Proportions of myeloid cells declined over time in patients with longer PFS, which was not observed in those with short PFS (Supplementary Fig. S3D). Similar differences in the cellular composition of the bone marrow were also observed with mass cytometry (Fig. 1D). For example, at D28, patients with longer PFS had a higher proportion of T cells and a lower proportion of CD14+ myeloid cells, as detected by mass cytometry (Fig. 1D). Data from mass cytometry showed a strong correlation with estimates of cell populations detected by CITE-seq (Supplementary Fig. S4). Transcriptome-based uniform manifold approximation and projection (UMAP) plotted based on both specimen time point (classified as pretherapy, D28 or 3 to 6 months following CAR T) and PFS status (PFS <180 days or >180 days) highlights differences in T cells and myeloid/DC clusters after BCMA-specific CAR T infusion when comparing patients by PFS (Fig. 1E). Together, these data suggest that changes in cellular composition, and particularly T and myeloid compartments, may influence the durability of response following CAR T therapy.

Binding to soluble BCMA was utilized to detect BCMA-specific CAR T cells, by both mass cytometry and CITE-seq, as in other studies (8). A representative gating strategy for detection by mass cytometry is shown in Supplementary Fig. S5. Detection of CAR T cells correlated with previous qPCR data for CAR vector copies (ref. 4; Supplementary Fig. S6A and S6B). ViSNE analysis of the mass cytometry data revealed increased proportions of both CD4+ and CD8+ T cells following CAR T therapy (Fig. 2A). The increased proportion of CD8+ T cells after CAR T infusion was dominated by a distinct population of T cells with CD27+TCF1+ granzyme B–CD57 phenotype (Fig. 2A and B). The increased proportion of this population was detected in patients with both short or longer PFS (Supplementary Fig. S7A and S7B), and this population included CAR T as well as non-CAR T cells (Supplementary Fig. S8). This T-cell population is distinct from granzyme B+KLRG+ terminal effector T cells that are prominent in the pretreatment marrow and previously described as increased in MM marrow (7).

Hierarchical consensus clustering of CD4+ T cells analyzed by mass cytometry using FlowSOM identified 15 clusters, of which the proportion of cluster 4 was increased at D28 (Supplementary Fig. S9A and S9B). T cells in this cluster were of the CD27+CD127+TCF1+ phenotype (Supplementary Fig. S9C). Clustering of CD4+ T cells by transcriptome also identified an increased proportion of a distinct cluster (cluster 0) at D28, which was enriched for expression of CD27, CD127, and TCF1 (Supplementary Fig. S9D–S9F). Clustering of mass cytometry data for CD8+ T cells similarly identified 10 distinct clusters (Supplementary Fig. S10A). Of these, the proportion of cluster1 (expressing CD27 and TCF1) was increased and that of cluster 8 (expressing T-bet and granzyme B) decreased at D28 following CAR T therapy (Supplementary Fig. S10B and S10C). Similarly, clustering of CD8+ T cells by transcriptome identified 10 clusters (Supplementary Fig. S10D). Of these, the proportion of CD8+ T cells in clusters 6 and 7 was increased at D28. T cells in this cluster expressed TCF1/7 and CD27 more highly compared with other clusters. They also lacked expression of granzyme B and TBET/TBX21 as seen with mass cytometry (Supplementary Fig. S10E and S10F). CAR T cells were detected as a component of transcription-based clusters for both CD4 and CD8+ T cells and did not form unique clusters (Supplementary Fig. S11A and S11B). Together, these data provide a detailed analysis of the heterogeneity of T cells in the bone marrow of patients undergoing CAR T therapy and illustrate that CAR T therapy is associated with distinct changes in bone marrow T cells as early as one month following therapy, which can be detected by orthogonal methods such as mass cytometry and CITE-seq.

Analysis of T cells by CITE-seq/transcriptome also demonstrated transcriptional differences in T cells and CAR T cells based on PFS. CD8+ CAR T cells at D28 from patients with short PFS exhibited higher expression of lytic granule genes consistent with more terminal effector differentiation (Fig. 2C; Supplementary Table S2 for full gene list). Similar changes, albeit less prominent than those in CAR T cells, were also observed in non-CAR CD8+ T cells (Fig. 2D; Supplementary Table S3 for the full gene list). As with CD8+ T cells, CD4 CAR T cells at D28 from patients with short PFS also demonstrated a more differentiated/effector phenotype, with a higher expression of genes such as IL32 (Supplementary Fig. S12; Supplementary Table S4 for the full gene list). Transcriptomes of CD8+ T cells from patients with longer PFS also differed from their counterparts with shorter PFS at 3 months after infusion (Fig. 2E; Supplementary Table S5 for the full gene lists). CD8+ T cells from these time points in longer PFS patients had higher expression of genes such as CXCR4 and CD69, previously associated with human bone marrow residence/retention (Fig. 2E; ref. 9). Correlation of increased CXCR4 and CD69 expression on T cells at this time point with longer PFS was also confirmed by concurrent antibody staining, both by mass cytometry (Supplementary Fig. S13A) and by CITE-seq (Supplementary Fig. S13B). Together, these data show that early T-cell expansion/infiltration in the bone marrow following CAR T therapy is dominated by a distinct population of T cells. These cells are characterized by a less exhausted phenotype and expression of tissue-retention genes over time in patients with longer PFS.

We also examined the T-cell receptor (TCR) diversity of marrow T cells by integrating CITE-seq/transcriptome data with TCR sequencing. TCR clonotypes were classified as hyperexpanded if they consisted of >10% of all TCRs sequenced. Mapping the hyperexpanded TCRs to one of the publicly available databases (tcrex.biodatamining.de) did not find these to be specific against common pathogens, such as cytomegalovirus. Hyperexpanded clones exhibited a phenotype of KLRG+HLA-DR+ terminally differentiated T cells expressing coinhibitory markers (Fig. 2F) and were present at baseline only in patients with lower PFS; they were undetectable at D28 but were also detected at later follow-up (Fig. 2G). As with bulk T cells, most of the CAR T cells were detected as single TCRs (Fig. 2H). TCR diversity was increased at D28 (Fig. 2I). Shannon diversity at baseline correlated with PFS and was lower in patients with shorter PFS (Fig. 2J). Together with the immunophenotypic studies, these data show that both baseline and posttreatment features of non-CAR T cells as well as CAR T cells correlate with the durability of the therapeutic response.

TCRs from patients with serially available data, including baseline marrow samples, were utilized to track changes in non-CAR T TCRs over time following CAR T therapy. Circos plots of the top 20 and all shared TCRs identified several TCRs that underwent expansion or contraction (in terms of clonal size) over time (Supplementary Fig. S14). Analysis of differentially expressed genes (DEG) or antibodies revealed that the TCRs that increase over time (e.g., change from single TCRs to small- or medium-size clones) exhibit an activated phenotype (reflected by higher HLA-DR expression) at baseline, consistent with prior immune recognition (Supplementary Fig. S15). Together, these data demonstrate that CAR T therapy can lead to the expansion of specific TCRs within tumor-infiltrating lymphocytes in vivo. When a similar analysis was performed for CAR T TCRs, shared TCRs between D28 and 3-month time point could be identified in only one patient (Supplementary Fig. S16).

CITE-seq/transcriptome analysis also suggested differences in the myeloid/DC populations based on response duration. UMAPs of myeloid/DCs based on transcriptomes identified 11 clusters, one of which was identified as plasmacytoid DCs (pDC; Fig. 3A, cluster 10). The remaining clusters fell into two distinct groups: one (group 1) enriched for monocyte/macrophage markers (such as CD14 and CD11b) was increased in patients with shorter PFS at D28 time point (Fig. 3A and B). In contrast, group 2, composed of clusters expressing markers associated with DCs (such as CLEC9A and CD1c), was enriched in patients with longer PFS, particularly at the 3-month time point (Fig. 3A and B). Group 1 monocyte/macrophage clusters included cells with expression of myeloid-derived suppressor cell-associated genes (such as S100A9, clusters 0 and 4) and a CD16+ cluster (cluster 2) expressing interferon response genes (Fig. 3C). Notably, these clusters were associated with higher expression of BAFF and PD-L1 among the top differentially overexpressed genes/markers (Fig. 3D). DEG analysis revealed higher expression of several genes implicated in immune-suppressive myeloid cells (such as S100A family and versican) in cells belonging to group 1 myeloid cells (Fig. 3E). Together, these data preliminarily suggest that differences in myeloid/DC populations in the bone marrow may correlate with PFS in MM patients treated with BCMA-specific CAR T therapy. Further studies are needed to explore this possibility.

Finally, we explored the changes in the tumor cell transcriptome in patients undergoing therapy. Comparison of transcriptional profiles of pretherapy tumor cells from patients with short versus longer PFS revealed that tumors in patients with longer PFS had higher expression of mature plasma cell genes (such as syndecan-1/SDC1, BCMA/TNFRSF17, and XBP1) and interferon response genes at baseline, whereas those with shorter PFS had higher expression of genes associated with less mature B cells (e.g., PBXIP) consistent with a less-differentiated phenotype (Fig. 4A; Supplementary Table S6 for the full gene list). Pathway analysis of DEGs revealed the enrichment of the epithelial–mesenchymal transition (EMT) pathway in patients with shorter PFS (Fig. 4B). Although there was a marked reduction in tumor cells at D28 post-CAR T infusion, a residual tumor could be detected in most patients at D28. Pretreatment and D28 posttherapy tumor cells were reclustered to analyze the properties of residual disease (Fig. 4C). Whereas pretreatment tumors form patient-specific clusters as in prior studies (10), posttreatment residual tumors from different patients cluster together (cluster 9 in Fig. 4D), suggesting the emergence of a shared transcriptional signature in residual disease. DEGs in this cluster included downregulation of plasma cell genes such as syndecan-1, BCMA, and XBP1, and upregulation of genes such as SPARC and MYL9 was implicated in EMT (Fig. 4E; Supplementary Table S7 for the full gene list). The pathway analysis of DEGs identified the enrichment of EMT and the downregulation of UPR in this cluster (Fig. 4F). In order to further evaluate patient-specific transcriptional changes in tumors over time, tumor cells from different time points in each patient were clustered together to identify transcriptionally distinct subclusters. These analyses revealed that in patients with short PFS, specific subclusters are similar between baseline and 3-month follow-up time points. However, in patients with longer PFS, follow-up samples were associated with the emergence of new dominant subclusters (Fig. 4G). In aggregate, these data demonstrate the plasticity of myeloma gene expression and suggest that myeloma cells may be evolving over time in the response to selective pressure imparted by CAR T cells and contribute to recurrent disease.

These studies utilize complementary single-cell tools to identify biological variables that correlate with the durability of response to BCMA CAR T therapy in MM. Our analysis of responding patients receiving an experimental BCMA-specific CAR T-cell therapy suggests that durable responses may depend not just on properties of CAR T cells but also on endogenous T-cell states, and properties of the TME, including the myeloid/DC compartment and tumor cells themselves. The initial expansion of T cells within the bone marrow following CAR T infusion includes both CAR T cells and non-CAR T cells, with a preferential expansion of a diverse TCR repertoire. The observation that a CD27+ TCF1+granzymeCD57 T-cell subset undergoes preferential in vivo expansion following CAR T therapy is similar to the expansion of TCF1+ T cells that provide the proliferative burst following immune-checkpoint blockade (11). This is also consistent with the finding that a higher diversity of TCRs at baseline and the absence of hyperexpanded TCR clones that bear markers of exhausted T cells correlates with longer PFS. These data are also consistent with prior studies that the phenotype of source T cells affects the potency of the CAR T product in MM (12). We have shown previously that the accumulation of terminally differentiated T cells and depletion of TCF1+ stemlike memory T cells begins early during the evolution of MGUS to clinical MM (7). Therefore, harvesting T cells earlier in the course of the disease evolution (as opposed to during earlier “lines of therapy”) should be considered to improve the quality of CAR T products and outcomes following CAR T therapy.

In addition to T cells, properties of both tumor and myeloid cells were also correlated with outcomes. Tumors with less-differentiated phenotypes were associated with shorter PFS. In addition, the presence of myeloid cells that express genes, such as BAFF, previously implicated in promoting the growth and clonogenicity of MM cells (13), correlated with shorter PFS. Although the impact of the myeloid component of the TME in CAR T-treated MM patients is less studied, a correlation between circulating myeloid cells and outcomes following CAR T therapy has also been reported in lymphoma (14). In contrast to immune-suppressive/growth-promoting myeloid cells, the presence of cells expressing conventional type 1 DC-associated markers, such as CLEC9A, is linked to better outcomes. This finding is consistent with the emerging role for DCs in enhancing the efficacy of adoptive cellular therapies (15). Our data suggest that durable responses depend in part on endogenous immunity mediated by non-CAR T cells; however, further studies are clearly needed to dissect the role of the myeloid/DC compartment in MM and CAR T-cell therapy.

These data also provide insights into properties of the tumor cells that may affect outcomes following CAR T therapy. Residual tumor following early CAR T-mediated cytoreduction bears features of less-differentiated, stemlike cells with expression of genes seen in EMT, similar to other studies in patients treated with immune-modulatory drugs indicating plasticity of myeloma (16). As these residual tumor cells likely contribute to disease relapse, targeting stemlike pathways (such as SOX2) previously implicated in MM pathogenesis should be explored in combination with BCMA CAR T cells to better eradicate residual disease (16–18).

A major strength of this study is the analysis by multiple complementary single-cell methods of tumor microenvironment specimens serially collected at more than one time point following CAR T therapy. This is notable, as most prior correlative studies for CAR T therapy have been based on analyses of circulating immune cells. The use of multiple single-cell analysis approaches also yields corroborative data. Potential weaknesses include the relatively small sample size and the variable nature of CAR T products that are made individually for each patient, which affects CAR T efficacy. Future studies should also address spatial aspects of immune response, which was not examined in the current study. It is appreciated that PFS rates in different BCMA CAR T trials in MM are different, but it is difficult to compare between studies, as these studies also differ in terms of patient selection and cohorts. It is notable, however, that all BCMA CAR T trials in MM, including in this trial, have utilized CAR constructs of similar structure with 4-1BB as the costimulatory domain.

Our results may inform the future design of CAR T approaches to improve outcomes in MM. Baseline TCR diversity may help identify patients most likely to achieve durable remissions; restoring or preserving TCR diversity during current MM therapy may be key to achieving durable immune control following T cell–based immunotherapy. These data also support the evaluation of combination approaches with strategies to either deplete immune-suppressive/growth-promoting S100A9/BAFF-expressing myeloid cells or recruit Clec9a+ type 1 conventional DCs to the tumor bed, to improve outcomes following CAR T therapy. Increased expression of PD-L1 on myeloid cells linked to short PFS also supports future studies exploring the combination of CAR T cells with PD-L1 blockade, especially in light of studies showing the capacity of PD-L1 blockade to reprogram the myeloid compartment and enhance DC function (19). The capacity of residual tumors to evolve following prolonged remissions suggests that early eradication of residual disease may be essential to realize the curative potential of these therapies.

Patients and Specimens

All patients treated in this study were enrolled in a prospective clinical trial (NCT02546167, www.clinicaltrials.gov), as described (4). Bone marrow specimens were collected following informed consent approved by the institutional review board at the University of Pennsylvania. All studies were performed following written informed consent from patients and in accordance with the Declaration of Helsinki. Bone marrow mononuclear cells (BMMNC) were isolated using density gradient centrifugation followed by cryopreservation. The choice of samples utilized for this study was based on the availability of archived specimens, and prioritized patients who responded to BCMA CAR T therapy, as the objective was to explore correlates of durable versus nondurable response (Supplementary Table S1). All specimens were collected at time points defined by the clinical protocol. PBMCs were isolated by density gradient centrifugation and cryopreserved in liquid nitrogen till further use. Aliquots were thawed and utilized for all planned assays together, to help reduce batch effects. All specimens with adequate cell counts after thawing and those that met quality control criteria as defined below under CITE-seq analysis were utilized for analysis. In the setting of limited specimens or yields following thawing, CITE-seq/transcriptome analysis was prioritized over mass cytometry.

Mass Cytometry

Thawed BMMNCs were stained with custom panels of metal-conjugated antibodies at manufacturer-suggested concentrations (Fluidigm; antibodies as noted in Supplementary Table S8). Cells were fixed, permeabilized, and washed in accordance with the manufacturer's cell-surface and nuclear staining protocol as described. After antibody staining, cells with incubated with intercalation solution, mixed with EQ Four Element Calibration Beads (catalog no. 201708), and acquired using a Helios mass cytometer (all from Fluidigm). Gating and data analysis were performed using Cytobank (https://www.cytobank.org). Viable cells and doublets were excluded using a cisplatin intercalator and DNA content with an iridium intercalator. Equal numbers of cells from each donor were utilized when data were concatenated prior to analysis.

CITE-seq and TCR-seq

After thawing BMMNCs, cells were incubated at 37°C for 1 hour to rest, resuspended in FACS buffer (PBS (Corning) with 0.5% FBS (Gibco), 90 μL), and stained with 10 μL of Fc-tagged, biotinylated recombinant human BCMA (Creative BioMart). Following a 30-minute incubation on ice, cells were twice washed and stained using a TotalSeq-C antibody cocktail (Supplementary Table S9) following the 10× Genomics protocol for Chromium Single-Cell Immune Profiling with Feature Barcoding Technology (ver. 1.0). Single cells were isolated using the Chromium Controller (10× Genomics, PN-110203). Gene expression, CITE-seq, and TCR libraries were prepared using the following kits from 10× Genomics: Chromium Single-Cell 5′ Library and Gel Bead Kit (PN-1000006), Chromium Single-Cell 5′ Library Construction Kit (PN-1000020), Chromium Single-Cell V(D)J Enrichment Kit, Human T cell (PN-1000005), Chromium Single-Cell A Chip Kit (PN-120236), Chromium i7 Multiplex Kit (PN-120262), Chromium i7 Multiplex Kit N, Set A (PN-1000084). The quality of the prepared libraries was assessed using the Agilent HS Bioanalyzer 2100, and sequencing was conducted on an Illumina NovaSeq 6000 (paired end 26 × 91 bp) at a targeted depth of 50K reads per cell for gene-expression libraries and 5K reads per cell for CITE-seq and TCR-seq libraries. To minimize variability, all samples were sequenced together.

Data Analysis for CITE-seq

Reads were aligned to the human reference sequence GRCh38/hg38, filtered, deduplicated, and converted into a feature barcode matrix using Cell Ranger 6.0. Samples with low fraction reads (<70%), high reads mapped to antisense strand (>10%), and with problematic barcode rank plot were excluded. Final samples with feature barcode matrix containing data for genes and 55 surface proteins were analyzed using Seurat V3 (20). Cells with fewer than 200 unique sequenced genes, more than 10% mitochondrial gene, and more than 70,000 sequenced features were filtered to exclude dead cells and doublets. The data were batch corrected using the canonical correlation analysis (CCA) method using the Seurat V3 data integration pipeline. A benchmark comparison study of different batch correction methods performed by Tran and colleagues (21) revealed Seurat 3 CCA as one of three preferred batch integration techniques for such data. Gene expression for each cell was normalized using R package SCTransform (22) and fitting the Gamma-Poisson generalized linear model. For each gene, unique molecular identifier (UMI) counts are considered as the response and cellular sequencing depth as the predictor variable to obtain regularized parameter estimates. The estimated regularized regression parameter was used to transform the UMI counts for each gene into Pearson residuals serving as scaled gene-expression data. The final data set used for analyses consisted of 151,054 cells from 23 samples (pretreatment n = 6, day 28 n = 6, subsequent time points n = 11). Protein expression data were treated as compositional data and normalized independently of transcriptomic data using centered log-ratio normalization, where counts were divided by the geometric mean of the corresponding feature across cells, and log-transformed (23). Principal component analysis (PCA) was then performed using the top 3,000 genes ranked by residual variance. The elbow method was utilized to identify the number of principal components that explained maximum variability in data, which resulted in 50 as an optimal number of principal components. Local neighborhood for each cell was defined by taking the 20 nearest neighbors in the kNN-Graph calculated using Euclidean distance in the PCA space. Rare and nonconvex cell populations were further identified using shared nearest neighbor networks by calculating neighborhood overlap using the Jaccard index. Finally, the Louvain community detection algorithm with a resolution of 1.2 was used to identify clusters of cells with similar gene-expression profiles. These clusters were further classified into major categories, B, T/NK, myeloid/DC, precursor, and tumor cell types, based on the protein expression of lineage-associated markers. Cells in T, myeloid, and tumor categories were also analyzed separately, and subclustering within each category was performed based on transcriptomic and protein expression data. The subcluster identity within each category was determined using significantly DEGs and proteins using the Wilcoxon rank-sum test between each cluster and the rest. DEGs were denoted as statistically significant for the FDR adjusted P < 0.05 with a fold change exceeding 1.2× (i.e., ≥1.2 or ≤0.83). The gene set enrichment analysis (GSEA; ref. 24) method was used for pathway analysis based on preranked genes. The Hallmark (25) gene set was utilized to identify pathways, and those with q-value < 0.05 were retained.

Analysis of TCR Sequencing Data

TCR sequencing reads from FASTQ were aligned to vdj_GRCh38_alts_ensembl-5.0.0 using Cell Ranger-6.0.2 to generate single-cell V(D)J (variable, diversity, and joining) segments for alpha and beta chains. The unique combination of alpha- and beta-chain genes is referred to as clonotype. The clonotype data for TCR sequencing were mapped to a Seurat object containing CITE-seq data by matching cell barcodes. Clonotypic frequency obtained from Cell Ranger were categorized as hyperexpanded, large, medium, small, and single as described elsewhere (26). Diversity was assessed by calculating the Shannon equitability index, which is Shannon diversity index divided by the maximum diversity resulting in normalized Shannon diversity index values between 0 and 1. Lower Shannon equitability index values indicate less evenness (i.e., some clonotypes have a higher frequency or are expanded compared with others), whereas 1 indicates complete evenness (i.e., all clonotypes have the same frequency).

Statistical Analysis

Statistical analysis of mass cytometry data was performed using 2D graphing and statistics software GraphPad Prism. Nonparametric Mann–Whitney (for comparing two groups) and Kruskal–Wallis (for comparing three groups) tests with a significance threshold of P < 0.05 were used to compare different cohorts. The Wilcoxon rank-sum test with a significance threshold of P < 0.05 after Benjamini–Hochberg correction for false discovery rate was used to identify DEGs between clusters and disease states in the scRNA-seq data. Data in bar graphs were plotted as mean ± SEM.

Data Availability

The sequencing data included in this article are available on Gene-Expression Omnibus through GSE210079. Other data are available through the corresponding author.

A.D. Cohen reports grants from Novartis during the conduct of the study; grants and personal fees from GlaxoSmithKline, personal fees from Genentech/Roche, Janssen, Takeda, BMS/Celgene, Oncopeptides, Pfizer, AbbVie, and Ichnos outside the submitted work; in addition, A.D. Cohen has a patent for 17/042,129 issued, licensed, and with royalties paid from Novartis and a patent for 16/050,112 issued, licensed, and with royalties paid from Novartis. A.L. Garfall reports grants from Novartis during the conduct of the study; grants and personal fees from Janssen, grants from Tmunity, personal fees from GlaxoSmithKline, Amgen, Legend Biotech, and grants from CRISPR Therapeutics and the Leukemia andLymphoma Society outside the submitted work; in addition, A.L. Garfall has a patent for US15/757,123 pending, licensed, and with royalties paid from Novartis, a patent for US16/764,459 pending, and a patent for US16/768,260 pending. S.F. Lacey reports a patent for Kymriah and related biomarkers licensed to Novartis. J.J. Melenhorst reports grants from Novartis during the conduct of the study; grants and personal fees from IASO Biotherapeutics, personal fees from Poseida Therapeutics, Chroma Medicine, and Gilead outside the submitted work; in addition, J.J. Melenhorst has a patent for methods of making chimeric antigen receptor-expressing cells licensed to Novartis, a patent for methods for improving the efficacy and expansion of immune cells licensed to Novartis, and a patent for biomarkers predictive of therapeutic responsiveness to chimeric antigen receptor therapy and uses thereof licensed to Novartis. C.H. June reports grants and personal fees from Novartis during the conduct of the study; grants from Tmunity Therapeutics outside the submitted work; in addition, C.H. June has a patent for Novartis issued to IPR licensed to Novartis. M.C. Milone reports grants from Novartis during the conduct of the study; other support from Verismo Therapeutics outside the submitted work; in addition, M.C. Milone has a patent for US-2016046724-A1 issued, licensed, and with royalties paid from Novartis. M.V. Dhodapkar reports other support from Sanofi, Lava Therapeutics, and Janssen outside the submitted work. No disclosures were reported by the other authors.

K.M. Dhodapkar: Conceptualization, resources, formal analysis, supervision, funding acquisition, writing–original draft, project administration, data interpretation and analysis. A.D. Cohen: Conceptualization, writing–review and editing, study design, sample collection, data interpretation and analysis. A. Kaushal: Formal analysis, visualization, writing–review and editing, data analysis. A.L. Garfall: Conceptualization, writing–review and editing, sample collection, data interpretation and analysis. R. Manalo: Writing–review and editing, sample processing, single-cell sequencing. A.R. Carr: Writing–review and editing, mass cytometry. S.S. McCachren: Writing–review and editing, data interpretation and analysis. E.A. Stadtmauer: Writing–review and editing, sample collection. S.F. Lacey: Writing–review and editing, sample processing and storage. J.J. Melenhorst: Writing–original draft, sample processing and storage. C.H. June: Conceptualization, writing–review and editing, data interpretation and analysis. M.C. Milone: Conceptualization, writing–review and editing, data interpretation and analysis. M.V. Dhodapkar: Conceptualization, resources, supervision, funding acquisition, writing–original draft.

M.V. Dhodapkar is supported in part by funds from NIH R35CA197603, Specialized Center of Research Award from the Leu­kemia and Lymphoma Society, and from the Paula and Roger Rinney Foundation. K.M. Dhodapkar is supported in part by funds from CA238471 and AR077926. A.L. Garfall is supported in part by LLS Scholar in Clinical Research Award. The authors acknowledge Winship Mass cytometry resource (R. Radzievsky and D. Doxie) and Emory Integrated Genomics Core (supported by P30CA138292). The samples used in this study were collected as part of a clinical trial supported by a sponsored research agreement between Novartis and the University of Pennsylvania, as well as NIH grant 1P01CA214278. We wish to acknowledge the University of Pennsylvania Translational and Correlative Science Laboratory for sample processing and storage, the Clinical Cell and Vaccine Production Facility for the generation of CAR T-cell products, and the Center for Cellular Immunotherapies Clinical Trials Unit for clinical and operational support during the trial.

Note: Supplementary data for this article are available at Blood Cancer Discovery Online (https://bloodcancerdiscov.aacrjournals.org/).

1.
Cohen
AD
,
Raje
N
,
Fowler
JA
,
Mezzi
K
,
Scott
EC
,
Dhodapkar
MV
.
How to train your T cells: overcoming immune dysfunction in multiple myeloma
.
Clin Cancer Res
2020
;
26
:
1541
54
.
2.
Munshi
NC
,
Anderson
LD
Jr
,
Shah
N
,
Madduri
D
,
Berdeja
J
,
Lonial
S
, et al
.
Idecabtagene vicleucel in relapsed and refractory multiple myeloma
.
N Engl J Med
2021
;
384
:
705
16
.
3.
Berdeja
JG
,
Madduri
D
,
Usmani
SZ
,
Jakubowiak
A
,
Agha
M
,
Cohen
AD
, et al
.
Ciltacabtagene autoleucel, a B-cell maturation antigen-directed chimeric antigen receptor T-cell therapy in patients with relapsed or refractory multiple myeloma (CARTITUDE-1): a phase 1b/2 open-label study
.
Lancet
2021
;
398
:
314
24
.
4.
Cohen
AD
,
Garfall
AL
,
Stadtmauer
EA
,
Melenhorst
JJ
,
Lacey
SF
,
Lancaster
E
, et al
.
B cell maturation antigen-specific CAR T cells are clinically active in multiple myeloma
.
J Clin Invest
2019
;
129
:
2210
21
.
5.
Samur
MK
,
Fulciniti
M
,
Aktas Samur
A
,
Bazarbachi
AH
,
Tai
YT
,
Prabhala
R
, et al
.
Biallelic loss of BCMA as a resistance mechanism to CAR T cell therapy in a patient with multiple myeloma
.
Nat Commun
2021
;
12
:
868
.
6.
Sharma
P
,
Hu-Lieskovan
S
,
Wargo
JA
,
Ribas
A
.
Primary, adaptive, and acquired resistance to cancer immunotherapy
.
Cell
2017
;
168
:
707
23
.
7.
Bailur
JK
,
McCachren
SS
,
Doxie
DB
,
Shrestha
M
,
Pendleton
K
,
Nooka
AK
, et al
.
Early alterations in stem-like/resident T cells, innate and myeloid cells in the bone marrow in preneoplastic gammopathy
.
JCI Insight
2019
;
5
:
e127807
.
8.
Reichman
A
,
Kunz
A
,
Joedicke
JJ
,
Hopken
UE
,
Keib
A
,
Neuber
B
, et al
.
Comparison of FACS and PCR for detection of BCMA-CAR-T cells
.
Int J Mol Sci
2022
;
23
:
903
.
9.
Boddupalli
CS
,
Nair
S
,
Gray
SM
,
Nowyhed
HN
,
Verma
R
,
Gibson
JA
, et al
.
ABC transporters and NR4A1 identify a quiescent subset of tissue-resident memory T cells
.
J Clin Invest
2016
;
126
:
3905
16
.
10.
Kaushal
A
,
Nooka
AK
,
Carr
AR
,
Pendleton
KE
,
Barwick
BG
,
Manalo
J
, et al
.
Aberrant extrafollicular B cells, immune dysfunction, myeloid inflammation, and MyD88-mutant progenitors precede Waldenstrom macroglobulinemia
.
Blood Cancer Discov
2021
;
2
:
600
15
.
11.
Im
SJ
,
Hashimoto
M
,
Gerner
MY
,
Lee
J
,
Kissick
HT
,
Burger
MC
, et al
.
Defining CD8+ T cells that provide the proliferative burst after PD-1 therapy
.
Nature
2016
;
537
:
417
21
.
12.
Garfall
AL
,
Dancy
EK
,
Cohen
AD
,
Hwang
WT
,
Fraietta
JA
,
Davis
MM
, et al
.
T-cell phenotypes associated with effective CAR T-cell therapy in postinduction vs relapsed multiple myeloma
.
Blood Adv
2019
;
3
:
2812
5
.
13.
Kukreja
A
,
Hutchinson
A
,
Dhodapkar
K
,
Mazumder
A
,
Vesole
D
,
Angitapalli
R
, et al
.
Enhancement of clonogenicity of human multiple myeloma by dendritic cells
.
J Exp Med
2006
;
203
:
1859
65
.
14.
Jain
MD
,
Zhao
H
,
Wang
X
,
Atkins
R
,
Menges
M
,
Reid
K
, et al
.
Tumor interferon signaling and suppressive myeloid cells are associated with CAR T-cell failure in large B-cell lymphoma
.
Blood
2021
;
137
:
2621
33
.
15.
Spranger
S
,
Dai
D
,
Horton
B
,
Gajewski
TF
.
Tumor-residing Batf3 dendritic cells are required for effector T cell trafficking and adoptive T cell therapy
.
Cancer Cell
2017
;
31
:
711
23
.
16.
Verma
R
,
Branagan
AR
,
Xu
ML
,
Flavell
RA
,
Dhodapkar
KM
,
Dhodapkar
MV
.
Role of MBD3-SOX2 axis in residual myeloma following pomalidomide
.
Leukemia
2021
;
35
:
3319
23
.
17.
Garfall
AL
,
Stadtmauer
EA
,
Hwang
WT
,
Lacey
SF
,
Melenhorst
JJ
,
Krevvata
M
, et al
.
Anti-CD19 CAR T cells with high-dose melphalan and autologous stem cell transplantation for refractory multiple myeloma
.
JCI Insight
2018
;
3
:
e127684
.
18.
Spisek
R
,
Kukreja
A
,
Chen
LC
,
Matthews
P
,
Mazumder
A
,
Vesole
D
, et al
.
Frequent and specific immunity to the embryonal stem cell-associated antigen SOX2 in patients with monoclonal gammopathy
.
J Exp Med
2007
;
204
:
831
40
.
19.
Bar
N
,
Costa
F
,
Das
R
,
Duffy
A
,
Samur
MK
,
McCachren
SS
, et al
.
Differential effects of PD-L1 versus PD-1 blockade on myeloid inflammation in human cancer
.
JCI Insight
2020
;
5
:
e129353
.
20.
Stuart
T
,
Butler
A
,
Hoffman
P
,
Hafemeister
C
,
Papalexi
E
,
Mauck
WM
3rd
, et al
.
Comprehensive integration of single-cell data
.
Cell
2019
;
177
:
1888
902
.
21.
Tran
HTN
,
Ang
KS
,
Chevrier
M
,
Zhang
X
,
Lee
NYS
,
Goh
M
, et al
.
A benchmark of batch-effect correction methods for single-cell RNA sequencing data
.
Genome Biol
2020
;
21
:
12
.
22.
Stoeckius
M
,
Hafemeister
C
,
Stephenson
W
,
Houck-Loomis
B
,
Chattopadhyay
PK
,
Swerdlow
H
, et al
.
Simultaneous epitope and transcriptome measurement in single cells
.
Nat Methods
2017
;
14
:
865
8
.
23.
Hafemeister
C
,
Satija
R
.
Normalization and variance stabilization of single-cell RNA-seq data using regularized negative binomial regression
.
Genome Biol
2019
;
20
:
296
.
24.
Subramanian
A
,
Tamayo
P
,
Mootha
VK
,
Mukherjee
S
,
Ebert
BL
,
Gillette
MA
, et al
.
Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles
.
Proc Natl Acad Sci U S A
2005
;
102
:
15545
50
.
25.
Liberzon
A
,
Birger
C
,
Thorvaldsdottir
H
,
Ghandi
M
,
Mesirov
JP
,
Tamayo
P
.
The molecular signatures database (MSigDB) hallmark gene set collection
.
Cell Syst
2015
;
1
:
417
25
.
26.
Borcherding
N
,
Bormann
NL
,
Kraus
G
.
scRepertoire: an R-based toolkit for single-cell immune receptor analysis
.
F1000Res
2020
;
9
:
47
.

Supplementary data