Gut microbiome modulates cytokine release syndrome and therapeutic response to CAR-T therapy in hematologic malignancies

Chimeric antigen receptor (CAR)-T cell therapy has emerged as a promising immunotherapeutic treatment for hematologic malignancies. By comparing the diversity and composition of the gut microbiome during different stages of CAR-T therapy, signicant changes were detected, not only in patients with relapsed/refractory multiple myeloma (MM; n = 43), but also in those with acute lymphocytic leukemia (ALL; n = 23) and non-Hodgkin lymphoma (NHL; n = 12). Analysis of treatment responses revealed signicant temporal differences in diversity and abundance of Bidobacterium, Prevotella, Sutterella, and Collinsella between MM patients in complete remission (n = 24) and those in partial remission (n = 11). Furthermore, we found that patients with severe cytokine release syndrome (CRS) exhibited higher abundance of Bidobacterium, Leuconostoc, Stenotrophomonas, and Staphylococcus. This study has important implications for understanding the biological role of the microbiome in the CAR-T treatment of patients with hematologic malignancies (ChiCTR1800017404).


Introduction
B-cell-derived hematologic malignancies, including acute lymphoblastic leukemia (B-ALL), non-Hodgkin lymphoma (B-NHL), and multiple myeloma (MM), carry a high probability of relapse after conventional chemotherapy 1 . With novel therapeutic strategies incorporating monoclonal antibodies, bispeci c T-cell engager (BiTE) antibodies, and hematopoietic stem cell transplantation (HSCT), treatment outcomes have greatly improved 2,3,4 . However, some patients progress to relapsed/refractory (r/r) status, with a poor prognosis 5 . The 5-year overall survival (OS) rate generally is < 10% with a median OS of 3-6 months for patients with r/r B-ALL 6, 7 . The complete response (CR) rate is 7% with a median OS of 6.2 months for r/r diffuse large B-cell lymphoma (DLBCL) 8 . For r/r MM patients, the 1-year OS is about 40% 9 . There is an urgent need to explore novel treatment strategies for these malignancies.
Chimeric antigen receptor (CAR) T-cell therapy (approved by the U.S. Food and Drug Administration) recently emerged as promising for r/r B-ALL, DLBCL, and mantle cell lymphoma (MCL) 10,11,12 . In multiple myeloma, investigations targeting the B-cell maturation antigen (BCMA) yielded encouraging outcomes with reversible toxic effects such as cytokine release syndrome (CRS) and pancytopenia 13,14,15,16,17 . However, the e cacy and toxicity have been inconsistent. No biomarker has been identi ed that predicts outcome and associated toxicities after CAR-T in patients.
Several studies have reported that the differences in diversity and composition of the gut microbiome might in uence cancer immunotherapy response 18,19,20,21 . After analyzing fecal samples from 43 melanoma patients treated with anti-programmed cell death 1 protein (PD-1) immunotherapy, signi cantly higher alpha diversity and abundance of Clostridiales/Ruminococcaceae were found in responders, whereas Bacteroidales were signi cantly enriched in non-responders 19 . In hematologic malignancies, intestinal bacteria also modulate the risk of graft-versus-host disease (GVHD) and infection after allogeneic hematopoietic stem cell transplantation (allo-HSCT). Greater bacterial diversity and abundance of the genus Blautia were associated with reduced GVHD-related death and improved OS 22,23 . However, no study has shown a potential role for the intestinal microbiota in the e cacy and toxicity of CAR-T therapy for B-cell malignancies.
The primary aims of this study were to understand the intestinal microbiome changes in patients with r/r B-cell-derived hematologic malignancies undergoing CAR-T cell treatment and to investigate associations of the microbiota with clinical responses and CRS severity. Finally, the potential of the gut microbiome to predict treatment outcomes and CRS severity was explored.

Patient cohorts
A total of 92 patients with r/r B-cell-derived hematologic malignancies were screened. Ten patients were not eligible for inclusion. Another four patients were excluded because of lack of su cient 16S sequencing depth. Thus, MM (n = 43), B-ALL (n = 23), and B-NHL (n = 12) patients were included (Fig. 1A).
The median age of the MM patients was 59 (range 39-75) years, and 55.8% were male ( Table 1). The median number of prior lines of therapy was 4 (range 2-8), with all receiving proteasome inhibitor therapy and 95.3% immunomodulatory agents. At enrollment, 39.5% had received autologous stem cell transplantation, and 55.8% had extramedullary disease(s).  1D). Two patients died: one from sepsis caused by Pseudomonas aeruginosa and the other from intracranial hemorrhage (Fig. 1D). Both the BCMA CAR-T/CD3 + T-cell percentages in peripheral blood (PB) and serum concentrations of interleukin (IL)-10 and interferon (IFN)-γ increased during CRS and differed signi cantly in the CR and PR groups (Fig. 1E). Patients' temperature and C-reactive protein (CRP), ferritin, and lactic dehydrogenase (LDH) concentrations were elevated, and IL-6 and IFN-γ concentrations were signi cantly different in grade 3 vs grade 1 CRS ( Fig. 1F and Supplementary Fig. 1A-C). The serum immunoglobulins (IgG, IgA) and immunoglobulin κ and λ light chain concentrations decreased dramatically after CAR-T (Supplementary Fig. 1D-F). Figure 1G shows the differences of positron emission tomography-computed tomography (PET-CT) scans and plasma cells detected by Wright's stain of a bone marrow smear (43.5% vs. 0), as well as ow cytometry (68.9% vs. 0) of bone marrow cells before and after CAR-T infusion for a representative subject.

Changes in the intestinal microbiome during CAR-T cell therapy
To detect changes in the gut microbiota during CAR-T, we collected fecal samples from each patient at ve times (FCa, FCb, CRSa, CRSb, and CRSc; Fig. 1C), where FCa denotes the baseline when patients were rst enrolled; FCb after chemotherapy; CRSa after CAR T-cell infusion but before the onset of CRS; and CRSb and CRSc the peak and during the recovery phase of CRS, respectively.
We rst evaluated the diversity of the gut microbiota in all subjects during CAR-T cell therapy. There was a signi cant decrease in diversity (measured by the Simpson index) during and after CRS (at CRSb and CRSc) compared with baseline ( Fig. 2A). This decrease was observed in the microbiome of patients receiving CAR-T therapy for r/r ALL ( Supplementary Fig. 4A) or r/r NHL ( Supplementary Fig. 4B). Refer to Supplementary Table 1 for details on the characteristics of r/r B-ALL and B-NHL patients. To further assess the similarity of composition between different therapy stages, we performed pairwise Spearman correlation analysis of operational taxonomic unit (OTU) level bacterial abundance (Fig. 2B) and found that stronger correlations emerged during the early stages with a ρ value of 0.71, 0.73, and 0.68, respectively, at FCa, FCb, and CRSa. Correlations between late stages (CRSb and CRSc) and early stages were weaker, suggesting that changes in microbiome composition might be related to CRS.
We next explored community structure and temporal shift of bacterial abundance at multiple taxonomic levels during CAR-T therapy. In myeloma, bacterial communities were dominated by Firmicutes and Bacteroidetes at the phylum level (Fig. 2C) and characterized by signi cant enrichment of Firmicutes and depletion of Bacteroidetes at the last two timepoints (Fig. 2D, E and Supplementary Fig. 4C). By applying the longitudinal analysis in the Qiime2 microbiome analysis platform, we detected changes in the gut microbial communities at taxonomic levels from phylum to genus ( Fig. 2F and Supplementary Table 2).
We further employed a negative binominal (NB) regression model-based time-course analysis to identify genera with signi cant temporal changes (Supplementary Table 3). Five genera were detected by both Qiime2 and maSigPro procedures, which included increases in Enterococcus, Lactobacillus, and Actinomyces and decreases in Bi dobacterium and Lachnospira ( Supplementary Fig. 4D). Most changes were aggravated during the late stages. Moreover, by checking changes in the ve genera in ALL and NHL patients, we observed consistent shift trends in NHL (four genera; Supplementary Fig. 4E) and ALL (two genera; Supplementary Fig. 4F), respectively.
Association between microbial communities and clinical response to CAR-T therapy We next determined whether microbial compositions or changes were associated with the response to CAR-T. Because we wanted to identify maximum differences and only six subjects presented in the VGPR group, we performed comparisons only between the CR and PR groups.
Notable differences in microbial alpha and within-sample diversity were observed in patients with CR and PR (Fig. 3A, B). Although no differences were detected at baseline, PR patients descended more dramatically in alpha diversity and had signi cantly lower Shannon indices than CR patients after CAR-T infusion (Fig. 3A). As the degree of differences between CR and PR groups changed across therapeutic stages, we characterized the periods with greater differences by summarizing the amount of CR/PRenriched OTU at each timepoint. The most pronounced differences occurred at CRSb (Fig. 3C).
To explore longitudinal differences between CR and PR across all therapeutic stages, we identi ed OTU features with differential dynamic pro les by applying negative binominal regression-based time-course differential analysis with the maSIgPro package. In total, 125 OTUs were found to have differential timecourse patterns between CR and PR patients ( Fig. 3D and Supplementary Table 4). The signi cant OTUs were further grouped into three clusters according to pro les of their abundance. Most of these OTUs were in clusters 1 and 2 ( Fig. 3E). Cluster 1, characterized by enrichment in the CR group, was comprised mainly of OTUs, which belong to the phyla Firmicutes and Bacteroidetes and the orders Clostridiales and Bacteroidales. Cluster 2 was comprised of OTUs from a broader taxonomy, which included the orders Clostridiales, Bacteroidales, Lactobacillales, and Actinomycetales (Fig. 3F).
We identi ed 30 genera with differential time-course patterns in patients with CR and PR after CAR-T (Supplementary Table 5). To explore these differences further, we divided the therapeutic period into before and after CAR-T infusion and performed genus-level class comparisons using linear discriminant analysis (LDA) of effect size (LEfSe) 24 . We detected 34 genera with differences in abundance in the CR and PR groups (Fig. 4A). Eighteen genera were detected by both procedures (Supplementary Fig. 5A). Consistent with the results from OTU-level pattern analysis, most of the signi cant genera such as Faecalibacterium, Roseburia, and Ruminococcus were enriched in CR patients after CAR-T. The genera Bi dobacterium, Prevotella, Sutterella, Oscillospira, Paraprevotella, and Collinsella had a higher abundance in CR versus PR patients both before and after CAR-T ( Fig. 4A and Supplementary Fig. 5B).
We also took patients with VGPR into consideration and analyzed the above-mentioned genera before and after CAR-T infusion. The bacterial abundance in VGPR patients fell somewhere between CR and PR patients, but no statistical signi cance was evident for most of genera ( Fig. 4B and Supplementary  Fig. 5D).
To explore whether early bacterial abundance was indicative of therapeutic response, we used RF feature selection to identify key discriminatory genera for responses 25 . By de ning the stages before CAR-T infusion as early, we applied feature selection procedures individually at both baseline (FCa) and postchemotherapy (FCb) and identi ed gut microbiome signatures comprising 8 and 14 discriminatory genera separately for baseline and post-chemotherapy (Fig. 4C, D and Supplementary Fig. 5C). The area under the receiver operating characteristic curve (ROC) of the two RF models using these discriminatory features was 0.73 and 0.85, respectively (Fig. 4E, F). Prevotella, Collinsella, Bi dobacterium, and Sutterella were enriched in CR versus PR both before and after CAR-T infusion and were identi ed by RF analysis as signi cant at baseline and post-chemotherapy. This indicates potential associations between these genera and the response to CAR-T.
We also checked the abundance of these genera in r/r NHL and ALL patients. In NHL, Faecalibacterium, Bi dobacterium, and Ruminococcus were signi cantly (or almost signi cantly) enriched in CR versus PR and in patients not having a remission (NR), consistent with our results in myeloma ( Supplementary   Fig. 5E). However, for ALL, we observed enrichment of Bi dobacterium, Roseburia, and Collinsella in NR ( Supplementary Fig. 5F), which differed from the results for MM and NHL but might be determined by the small NR sample.
To further demonstrate the association between these taxa and outcome, we assessed progression-free survival (PFS) following CAR-T therapy. By stratifying patients by tertile of bacterial abundance, we observed that for Sutterella, patients in the highest-abundance tertile had signi cantly prolonged PFS (Fig. 4G). Even after strati cation by timepoints, this association remained signi cant ( Supplementary  Fig. 6A). However, for genus Faecalibacterium, which was reported to be signi cantly associated with PFS and anti-PD-1 therapy 19 , we did not observe an association ( Supplementary Fig. 6B, C).
We performed pathway analysis using Phylogenetic Investigation of Communities by Reconstruction of Unobserved State (PICRUSt) and identi ed signi cant changes in amino acid metabolism (Fig. 4H Fig. 7B and Supplementary Fig. 8B). By analyzing associations between CRS grade and taxa at the genus level, we identi ed signatures discriminating severe from mild CRS, including decreases in amount of Bi dobacterium and Leuconostoc in patients with severe CRS (Fig. 5A and Supplementary Table 8). Bi dobacterium was increased in patients with worse CRS, not only during the window of CRS, but also at early stages (Fig. 5A, B). Leuconostoc was signi cantly enriched during the window in patients with high CRS grade (Fig. 5A, B). In addition, the abundance of Stenotrophomonas and Staphylococcus differed severe vs moderate CRS during the window (Supplementary Fig. 8D and Supplementary Table 9).
Comparisons of KEGG pathways across CRS groups showed that the gut microbiome of patients with severe CRS had high metabolism or biosynthesis related to in ammatory compounds, including several pathways associated with amino acid synthesis and metabolism, purine metabolism, lipoic acid metabolism, and biosynthesis of lipopolysaccharide and peptidoglycan (Supplementary Fig. 9 and Supplementary Fig. 10).
Primary in ammatory markers of CRS are cytokines, such as IL-6, IL-2, IL-10, interferon gamma (IFN-γ), and tumor necrosis factor-α (TNF-α). Various cytokines are elevated in the serum of patients experiencing CRS after CAR-T cell infusion 31 . By assessing serum cytokine concentrations and immune cell numbers during CAR-T, we observed signi cantly increased amounts of serum in ammatory cytokines (IL-6, CRP, IFN-γ, D-dimer, ferritin) but low numbers of immune cells (monocytes, lymphocytes, neutrophils, leukocytes) in severe CRS (Fig. 5C). We also compared serum cytokine concentrations and immune cell numbers in CR and PR, observing signi cant differences for many of them (see Supplementary Fig. 11A).
To explore further associations between the gut microbiome and CRS during CAR-T therapy, we determined whether serum cytokine concentrations and numbers of PB immune cells correlated with the abundance of gut microorganisms (Fig. 5D). The abundance of the genus Leuconostoc, previously linked to CRS grade, correlated positively with ferritin and D-dimer concentrations. The abundance Bi dobacterium correlated signi cantly negatively with PB monocytes (Fig. 5E). We also found a correlation between in ammatory markers and bacteria associated with the clinical response and PFS. For example, Sutterella correlated negatively with serum concentrations of CRP and D-dimer ( Supplementary Fig. 11B). Prevotella correlated negatively with the number of multiple PB immune cells but positively with the serum D-dimer concentration ( Supplementary Fig. 11B). Faecalibacterium correlated negatively with the serum concentrations of D-dimer and IFN-γ ( Supplementary Fig. 11B).

Discussion
Although several studies have revealed the critical role of the gut microbiome in treatment responses and survival after administration of another important immunotherapy -immune checkpoint inhibitor (e.g., PD-1, PD-L1) therapy 20 , no study has reported on the association between the gut microbiome and CAR-T therapy. In this study, we describe the changes of the gut microbiome during CAR-T therapy and associations with treatment responses and CRS severity in CAR-T-treated patients with B-cell malignancies.
Some of the bacterial genera with differences in abundance in CR versus PR patients have been reported to be involved in the regulation of the immune response, including to immunotherapy. Faecalibacterium, reported to enhance antitumor immune responses and survival after anti-PD-1 therapy in melanoma 19,32 , was in this study associated with CR. Multiple species within the genera Bi dobacterium and Collinsella increased in responders to anti-PD-1 therapy for melanoma 33 , resulting in depleted peripherally derived colonic regulatory T cells, increased Batf3-lineage dendritic cells (DCs), and augmented T-helper 1 cell (Th1) responses and thus better immune-mediated tumor control 34 . Here, we observed an increased abundance of these two bacteria in CR patients, suggesting a similar response-associated effect of these taxa on the immune system across cancer types and therapeutic strategies.
Nevertheless, some taxa might have effects that are speci c for cancer or therapy types. For example, high abundance of genus Sutterella was associated with both CR and prolonged survival after CAR-T therapy. However, previous studies reported higher numbers of Sutterella in non-responders versus responders in non-small-cell lung cancer (NSCLC) treated with nivolumab 35 . Besides, in this study, we observed contradictory results for the genus Bi dobacterium, Roseburia, and Collinsella in three types of hematologic malignancy (Supplementary Fig. 2F). This indicates a potentially distinct involvement or function of some bacteria in different cancer types and treatments. But these ndings require con rmation in studies with larger cohorts.
Gut microbial communities contribute to inter-individual variation in cytokine responses 36 . We propose that gut microbes are related to the intensity of CRS during CAR-T therapy. Bi dobacterium, Leuconostoc, Stenotrophomonas and Staphylococcus were enriched in myeloma patients with severe CRS. Additional studies also demonstrated an association between these microbes and cytokine production. Previous research showed that Bi dobacterium correlated with the production of multiple cytokines (e.g., IFN-γ) in The mechanisms through which gut microbes modulate host immunity are largely unknown. Gut microbial communities modulate host defenses mainly through the release of intermediary metabolites rather than by direct interaction between speci c microorganisms and immune cells 36 . Multiple bioactive gastrointestinal metabolites produced by gut microbes, such as amino acids, short-chain fatty acids (SCFAs; e.g., butyrate), and bile acids, exert immunomodulatory functions through immune cell metabolic reprogramming or transcriptional and epigenetic modulation of immune-related genes 26 . Lipopolysaccharide (LPS) from some pathogens is a well-known endotoxin that can stimulate the release of a variety of cytokines/chemokines 40,41 . Peptidoglycans in bacterial cell walls are a conserved PAMP that trigger innate in ammatory responses throughout the body 42 .
In addition to myeloma, CAR-T therapy has been applied other blood cancers and solid tumors. The link between the gut microbiome and different cancer types needs to be studied systematically. Our research describes associations between changes in the gut microbiome of CAR-T patients and clinical responses and survival. This will open an avenue for investigating the interaction of the gut microbiome and CAR-T cells and lead to novel ways to improve the therapeutic e cacy of CAR-T therapy by targeting the gut microbiome.
As one of the most prominent treatment strategies for hematologic malignancies, CAR-T cell therapy has recent received great attention. Here for the rst time, we found that the dynamic changes in the gut microbiome correlated signi cantly with therapeutic response and CRS during CAR-T treatment of hematologic malignancies (B-ALL, B-NHL, and MM). These ndings will aid the development of novel biomarkers for predicting treatment outcome and CRS severity, thereby optimizing the management of these patients while reducing potential toxicities. Peripheral blood mononuclear cells (PBMCs) were obtained from each patient by leukapheresis for CAR-T cell preparation. The puri ed CD3 + T cells were transduced with lentiviral vector to express BCMA CAR (Fig. 1B). Then the engineered T cells were expanded ex vivo under interleukin-2 stimulation. All patients received lymphodepletion with udarabine (30 mg/m 2 of body surface area daily on days − 4, -3, and − 2) and cyclophosphamide (500 mg/m 2 daily on days − 3 and − 2) followed by an infusion of BCMA CAR-T cells on day 0. The primary response outcome, de ned by the guidelines from the International Myeloma Working Group (IMWG) as a complete response (CR), very good partial response (VGPR), or partial response (PR) in the third month after CAR-T treatment 43,44 . CRS was graded by the Lee criteria 30 .

Microbiome sample collection and restoration
Gut microbiome samples were collected at ve timepoints (Fig. 1C)

Assessment of serum cytokine concentrations
All blood samples were stored at 4°C until centrifugation at 5000 rpm for 6 min. The supernatant liquids were quanti ed with the BD Cytometric Bead Array Human Th1/Th2/Th17 Cytokine Kit and its corresponding software (BD Biosciences) according to the manufacturer's instructions.

Assessment of CAR-T cell expansion and persistence
Serial PB samples were collected in BD Vacutainer K 2 EDTA tubes (BD Biosciences) before and after CAR-T cell infusion. The expansion of CAR-T cells in vivo was determined by detecting the CAR-T ratio continuously in PB as described 45,46 . BCMA CAR-T expression was assessed using biotin-SP-conjugated F(ab')2 fragment goat anti-mouse IgG, F(ab')2 fragment-speci c antibody, and the secondary staining Amplicon data processing Sequenced reads were demultiplexed according to barcodes. Paired-end reads were merged with the fastq_mergepairs command from VSEARCH v. 2.4.4 47 .The minimum length of overlap between pairedend reads was set to 5. Merged reads were then imported into Qiime2 (v. 2020.2) 48 . Jointed reads were processed by the qiime quality-lter q-score-joined command to lter sequences with low-quality scores.
Sequences were denoised with the Deblur work ow 49 . Amplicon sequence variants (ASVs) were summarized with the feature-table summarize command. To calculate phylogenetic diversity, a rooted phylogenetic tree was constructed using the align-to-tree-mafft-fasttree pipeline from the q2-phylogeny plugin within Qiime2. The pipeline performed a multiple sequence alignment of the ASV sequences and then masked the alignment to remove positions that are highly variable. The masked alignment was used to generate a phylogenetic tree by FastTree program 50 . Alpha and beta diversity matrices were generated through the q2-diversity plugin using the above-mentioned ASV feature

Bioinformatics and statistical analysis
Comparisons of alpha diversity and taxonomic abundances between two groups were conducted with the Wilcoxon rank-sum test, while comparisons among three or more groups were conducted using the Kruskal-Wallis rank-sum test. For beta diversity analysis, a PCoA plot was generated with weighted Unifrac distances. To test the signi cance of between-sample diversity alternation, permutational analysis of variance (PERMANOVA) was performed with the adonis function within the R package vegan.
The feature-volatility plugin 53 within Qiime2 was applied to implement longitudinal analysis to identify features that are associated with therapy stages. In this pipeline, supervised learning regressor was used to identify important features and assess their ability to predict therapy states. Unclassi ed taxonomic features, features absent in more than 90% of all samples, and features with low abundance (< 0.01%) were all excluded from the analysis. Net average change scores and importance scores, which denote the correlation between input features and therapy stages, were exported and visualized in a volcano plot.
Only features with net average change scores more than 0.2% and importance scores within the rst tertile of distribution were retained.
For time-course differential analysis, the R package maSigPro 54, 55 was used to nd taxonomic features with signi cant temporal changes and signi cant differences between experimental groups (e.g., clinical response and CRS grade groups). Speci cally, the maSigPro algorithm de ned a generalized regressive model by dummy variables followed by two regression steps: the rst one selects features with non-at pro les by the least-squared technique and the second step creates best regression models for each feature by using stepwise regression to identify features with different pro les between experimental groups. We used as input, the normalized relative abundance (scaled to 100 million) and excluded features that did not occur in more than 90% of all samples. We employed a negative binominal regressive model for the microbial counts data and ran maSigPro on therapy stages with a degree of 4. All features with a signi cant group difference were exported. The signi cant features were further clustered together using the hclust function method according to the patterns of their relative abundance. For each cluster, a median pro le and tted curve of all included features were summarized to visualize the pro le pattern.
The LAD effect size (LEfSe) algorithm 24 was employed to identify differentially abundant features between groups (e.g., between clinical response and CRS grade). The method rst detected features with signi cant differential abundance using the non-parametric factorial Kruskal-Wallis rank-sum test with pre-de ned α of 0.05. Signi cant features were then used to build a Linear Discriminant Analysis (LDA) model for estimating the effect size of each differentially abundant feature. The LDA score threshold for discriminative features was set to 2.0.
To identify early predictive biomarkers with respect to clinical response (PR vs. CR), we implemented a random forest (RF) feature selection procedure within the R package caret. The recursive feature elimination (RFE) algorithm with 5-fold cross validation was applied for feature selection. An optimized number of feature sets was determined by performance of 5-fold cross validation. To depict the receiver operating characteristic (ROC) curve and calculate the area under the curve (AUC), the pROC package was utilized.
For progression-free survival (PFS) analysis, subjects were classi ed as high, medium, or low based on tertiles of the distribution of speci c taxa abundance (e.g., genus Sutterella). Time to progression was de ned as the interval (in days) from the date of CAR T-cell infusion to the date of disease progression.