Signatures of GVHD and relapse after posttransplant cyclophosphamide revealed by immune profiling and machine learning

Delayed immune reconstitution is associated with enhanced infectious and relapse risk after allogeneic blood or marrow transplantation (alloBMT)^1-3 and was an early barrier to successful alternative donor transplantation.^4-7 Posttransplantation cyclophosphamide (PTCy) enabled HLA-haploidentical (haplo) alloBMT by limiting excess nonrelapse mortality (NRM) and has been increasingly used in the HLA-matched setting.^8-13 Mechanisms responsible for the reduction in acute and chronic graft-versus-host disease (aGVHD and cGVHD) after PTCy and characterization of the unique immune milieu evident in the different patterns of GVHD remain poorly understood.¹⁴ 

PTCy has been postulated to act as an alloreactive T cell-depleting agent. However, recent studies have revealed relative sparing of alloreactive T-regulatory cells (Tregs) and association with hyporesponsive effector T cell phenotypes, suggesting GVHD-prevention mechanisms distinct from that of conventional prophylaxis.^15-19 Furthermore, data showing that PTCy may promote NK cell depletion followed by recovery to a functionally immature state underscores the importance of understanding not only the kinetics of immune recovery, but also markers of functionality and their influence on clinical outcomes.^20,21 Our previous work in PTCy-treated patients showed a correlation between day 30 levels of STimulation-2 (ST-2), CD25, tumor necrosis factor receptor 1 (TNFR1), and regenerating islet-derived 3-α with NRM and CD25 with aGVHD development.²² However, the biologic underpinnings of these findings, along with the interactions among proteomic, immunologic, and transcriptional signatures and their impact on survival after PTCy-based allografting have not been defined.

One of the difficulties in addressing complex posttransplant interactions is the abundance of contributing factors and the many time points at which they can be examined. Machine learning employs data-driven statistical modeling approaches that integrate multiple variables in a predictive model to identify underlying patterns and cut-points without predefined assumptions. This approach eliminates the effect of selection bias in choosing and categorizing variables and identifies associations that may not emerge in univariate analyses. Machine learning tools have been used to profile the association of immunological variables with clinical outcomes in several disease states.^23-25 Here, we use them to identify the link between biologic variables and clinical outcomes, guide correlative and preclinical studies, and generate data supporting future translational interventions.

This unplanned analysis was performed on prospectively collected samples on 75 patients who received PTCy-based alloBMT using an HLA-matched sibling (MSD) or unrelated (MUD) donor (NCT00809276²⁶) and 70 patients who received haplo alloBMT (NCT00796562²⁷). All patients were treated with myeloablative conditioning and received T cell replete bone marrow (BM) grafts. Clinical data reflect the patient population as a whole with detailed clinical outcomes previously published.^26,27 A total of 17/92 patients were excluded from the HLA-matched clinical trials and 26/96 patients from the HLA-haploidentical clinical trials due to graft failure, early death, or too few samples.

Classification and regression tree (CART)²⁸ machine-based learning was applied for overall survival (OS) treated as a time-to-event outcome, and for aGVHD dichotomized as a binary outcome. The candidate factors for subgroup identification included 10 baseline patient and donor characteristics, 18 immune subsets, and 7 plasma markers at day 28 after alloBMT for OS and aGVHD. All the variables are widely used in the field and have been established to correlate with alloBMT outcomes.^29-38 Day 28 was chosen because earlier predictive tools may allow for the development of impactful adaptive interventions. The splitting criteria in our tree models were maximum depth of 2 levels with at least 30 observations in a node necessary for a split and at least 15 observations in any terminal node for predicting OS and aGVHD after day 28. For the GVHD prediction model, we excluded patients who experienced grade 2 to 4 aGVHD within 28 days after alloBMT as patients already developed the queried outcome and to limit possible overlap with engraftment syndrome. All engrafting patients were included for OS analyses, and a conditional inference tree by binary recursive partitioning³⁹ was applied for detecting and classifying risk groups. The stopping rule for the OS tree model was measured by Bonferroni-adjusted P values <.05 for the partial null hypothesis of a single input marker.

For random forest predicting models,^40,41 1000 trees were built, and the parameters for the optimal number of variables to try for splitting at each tree node and node size were tuned using out-of-bag (OOB) errors. Estimated variable importance (VIMP)⁴² of each predictor in the corresponding random forest model was based on the differences of prediction error between predictor when randomly permuted vs the observed values. We used the top-ranked predictors from our tree/random forest models to guide selection of subsequent RNA sequencing (RNA-Seq) transcriptional studies. Due to a limited number of events, our random forest models were evaluated via OOB estimators.

Detailed descriptions of the platforms, outcome definitions, and group comparisons, as well as correlative analysis methods (flow cytometry, serum biomarker, and RNA-Seq) are outlined in the Supplemental Material available on the Blood Web site. Murine alloBMT studies were pursued to validate select RNA-Seq findings with regard to GVHD and establish the platform for future preclinical investigations. Analyses were carried out with the statistical software R version 3.6.0 (R Foundation for Statistical Computing, Vienna, Austria). All reported P values were 2-sided with significance level <.05 for hypothesis generating.

We compared lymphocyte subset recovery among patients who received MSD MUD, and haplo alloBMT with PTCy (description of the groups provided in Table 1). First, we focused on characterizing the dynamics of the global lymphocyte reconstitution (Figure 1). Since there were no significant differences in analyzed parameters between MSD and MUD alloBMT recipients, MSD and MUD were combined into a single group matched for all subsequent analyses. Mean absolute lymphocyte counts (ALCs) on day 56 were 872 cells/µL after matched alloBMT and 1058 cells/µL after haplo alloBMT. Median ALCs were in the intraquartile range for normal donors (1500-2262 cells/µL) by 1 year after both haplo (median 1580 cells/µL) and matched (median 1320 cells/µL) alloBMT (Figure 1A). CD3⁺ T cells followed a similar recovery pattern and reached pre-BMT levels by 6 months to 1 year (Figure 1A). B cell counts exceeded pre-BMT levels and were comparable to normal donors by 6 months after haplo and HLA-matched alloBMT (Figure 1A). NK cells recovered to normal donor levels by day 56 in all cohorts. There were no significant differences for reconstitution in ALC, CD3⁺ T cells, B cells, and NK cells across the 3 PTCy cohorts (Figure 1A). However, analysis of T cell subsets revealed substantial differences at day 28 for both CD4⁺ and CD8⁺ T cells between HLA-matched and haplo alloBMT (P values < .001). By 2 months after alloBMT, CD4⁺ T cells were similar after HLA-matched and haplo alloBMT. However, by 1 year, CD4⁺ T cells were significantly lower after HLA-matched when compared with haplo alloBMT (218 cells/µL and 358 cells/µL, P < .05), and remained lower than the median donor counts (632, both P values < .0001). After day 28, CD8⁺ T cells were comparable between haplo and HLA-matched alloBMT, exceeded pre-BMT levels, and reached normal donor levels by 6 months.

Within CD4⁺ and CD8⁺ T cell subsets, naïve and central memory (CM) T cells were significantly lower than normal donors and pre-BMT levels for all time points. CM CD8⁺ T cells were significantly lower after haplo alloBMT for all time points examined (P < .0001), but naïve and CM CD4⁺ subsets were similar after HLA-matched and haplo alloBMT (Figure 1B-C). Effector memory (EM) and effector memory expressing RA (TEMRA) CD4⁺ and CD8⁺ were significantly lower at day 28 after haplo compared with HLA-matched alloBMT, but were not significantly different between donor type after that time point and reached the lower interquartile range of normal donor levels by 2 months. Overall, the dynamics of global lymphocyte reconstitution were similar in all cohorts with EM being the predominant phenotype of recovering T cells consistent with previous reports in HLA-matched alloBMT with conventional prophylaxis⁴³ and HLA-matched alloBMT using PTCy.¹⁷ 

Total numbers of Tregs were not significantly different by donor type, while conventional T cells (Tconv) were significantly lower at day 28 and significantly higher at 1 year (P < .001 and P < .05, respectively) after haplo alloBMT compared with HLA-matched recipients (Figure 1A-C). Development of aGVHD after day 28 was associated with decreased B cell counts but did not significantly affect overall immune reconstitution of other subsets (supplemental Figure 1). As such, subsequent analyses with additional markers sought to illuminate the interplay of molecular and cellular factors in the development of aGVHD.

Machine-learning approaches were then applied to investigate the relationship between cellular immune subsets and proteomic markers with clinical outcomes. For the whole cohort, incidences of grade 2 to 4 and grade 3 to 4 aGVHD were 38% (55/145) and 11% (16/145), respectively. We first used CART analyses to determine which cellular immune subsets and proteomic markers were associated with aGVHD development. We found that at day 28, higher CXCL9 levels accompanied by elevated Tconv levels were present in 28% of recipients and strongly associated with aGVHD development, with 60% of those patients experiencing aGVHD after day 28 (Figure 2A). In contrast, 52% of patients had low CXCL9 levels in whom there was an 11% incidence of aGVHD. Finally, 20% of patients fell into an intermediate risk category, having elevated CXCL9 but low Tconv, which was associated with a 24% incidence of aGVHD. Tconv and CXCL9 levels selected by CART analysis were used to perform cumulative incidence curves and significantly correlated with both grade 2 to 4 aGVHD and grade 3 to 4 aGVHD (P < .00001 and P = .0002, respectively) (Figure 2B).

Next, we performed VIMP conducted via random forest analysis to allow the identification of additional pertinent factors (Figure 2C). Similar to the CART analysis for aGVHD, CXCL9 was the top-ranked variable, followed by Tconv. The additional top-ranked variables included ST2, memory B cells, CD4⁺ T cells, BMT year, Reg3α, CD38⁺ CD4⁺ T cells, pre-BMT disease status, and age at BMT.

Evaluation of Treg homeostasis in both haplo and HLA-matched cohorts suggests that the development of aGVHD was not associated with insufficient Treg numbers or impaired Treg/Tconv ratios (Figure 3A). To evaluate for molecular signatures of post-PTCy aGVHD, we pursued RNA-Seq of recovering lymphocyte subsets on peripheral blood mononuclear cells (PBMCs) collected on day 28 from the HLA-matched alloBMT recipients. To assess the longitudinal changes following post-PTCy aGVHD, available samples collected on day +180 were also sorted and sequenced.

Differential gene expression analysis comparing Tregs from a subset of patients with or without aGVHD revealed 56 differentially expressed genes (DEGs) (30 downregulated) (Figure 3B; supplemental Table 3). In-depth analysis of the identified genes revealed that those necessary for bonafide Treg function and aGVHD protection, such as IFNG and GZMA,^44,45 were expressed at a significantly higher level in patients who developed aGVHD (Figure 3C). To define additional Treg transcriptional signatures and the relation of DEGs to aGVHD, we conducted gene set enrichment analysis (GSEA) using published gene sets. We found that the signatures of active and highly suppressive Treg profiles^46,47 were robustly enriched in Tregs from the cohort who developed aGVHD (NES 1.79, 2.16, and 3.38). Significant enrichment for these signatures was also present at later points after successful aGVHD treatment, but the difference was less pronounced (NES 1.50, 1.78, and 1.98) (Figure 3D).

Next, we examined and compared the transcriptional profiles of CD4⁺ Tconv and CD8⁺ effector T cells from patients according to the development of aGVHD. Applying the same significance criteria as above, we identified a small number of DEGs in the CD8⁺ T cells at day +28 (19 upregulated, 40 downregulated) (supplemental Figure 2; supplemental Table 4). However, we identified a total of 277 DEGs in the effector Tconv subset (Figure 4A; supplemental Table 5), of which 67 were upregulated in patients subsequently developing aGVHD. DEG analyses highlighted Th1-skewing and upregulation of interferon-regulated genes previously described as influencing GVHD development (CXCR3, TBX21).^48-50 Several inflammatory genes (F2R, LGLS1) were upregulated while those involved in apoptosis (IER3, BCL2) and the differentiation of pathogenic follicular helper cells (C5AR1)⁵¹ were downregulated in GVHD cases. To evaluate long noncoding RNAs, which are influential in T cell immunity, we identified Malat1⁵² and Linc00402⁵³ as downregulated in the aGVHD cohort (supplemental Figure 3).

The observation of CXCR3 upregulation reinforced our finding of increased serum levels of CXCL9, a CXCR3 ligand, as a predictor of aGVHD in our machine-learning analysis. To corroborate observed molecular changes linking Tconv and CXCL9 response, we performed flow cytometric analysis of patient PBMCs. Patients who developed grade 2 to 4 aGVHD showed significant enrichment for both CXCR3-expressing and Th17 (CD146⁺CCR5⁺)-prone Tconv when compared with patients without aGVHD and healthy donors (P < .05) (Figure 4B). Patients with GVHD also had a higher percentage of CD146⁺CCR5⁺ Tconv cells that coexpressed CXCR3. To extend the biological relevance of this finding, we evaluated changes in CXCR3 expression in a CD4-driven murine MHC-matched alloBMT model (B10.D2→BALB/c) (supplemental Methods). We confirmed CXCR3 upregulation during post-PTCy GVHD, observing increased numbers of CXCR3⁺ Tconv in the spleens and livers of mice developing breakthrough GVHD (supplemental Figure 4).

We next compared observed molecular changes with datasets curated in MSigDB (Figure 4C). GSEA revealed that in the aGVHD cohort Kyoto Encyclopedia of Genes and Genomes (KEGG) GVHD signature enrichment was present in both Tconv and CD8⁺ T cells, while Tconv were strongly enriched (NES >2) for the signatures of effector T cells (Figure 4C). We did not identify enrichment in the aGVHD cohort for pathways entertained as critical drivers of alloimmunity in aGVHD (Notch and IL6/JAK/STAT) (Figure 4C). Since the prior studies of PTCy activity in murine models suggested the creation of long-lasting T cell hyporesponsiveness¹⁸ reminiscent of exhaustion, we indexed our data to landmark datasets that define exhaustion profiles occurring during the course of chronic infection and tumor-induced T cell dysfunction.^46,47,54-56 Lack of GVHD was not associated with an enrichment for gene expression signatures of exhausted CD8⁺ T cells or Tconv (Figure 4F), suggesting that exhaustion is an unlikely mechanism behind PTCy protection against aGVHD. However, effector T cells sorted from patients who developed aGVHD exhibited the hallmarks of activation, including Th1/17 skewing in the Tconv compartment and a loss of a naïve CD4⁺ Tconv signature. To evaluate if establishment and persistence of operational tolerance after PTCy may be related to delayed development of an exhaustion-like phenotype, we pursued RNA-Seq analyses of CD8⁺ T cells and Tconv subsets from the aGVHD and the aGVHD-free cohort on day 180 after alloBMT. At that time, lack of enrichment for exhaustion in both CD8⁺ T cells and Tconv persisted, but in contrast to earlier time points, there was enrichment for a naïve T cell signature (Figure 4E-F).

T cell activation during GVHD is a biologically robust process that requires extensive adaptation to augment cellular energy generation capabilities. Both glycolysis⁵⁷ and oxidative phosphorylation (OXPHOS)⁵⁸ are suggested as dominant processes in preclinical aGVHD models, and metabolic reprogramming is increasingly considered a putative therapeutic target in aGVHD.⁵⁹ We identified a robust enrichment for OXPHOS, fatty acid metabolism, and less so glycolysis in Tconv subset isolated from patients with aGVHD after HLA-matched alloBMT (Figure 5A). As selective enrichment was dominant in Tconv over the CD8⁺ T cells subset (Figure 5A), we hypothesized that Tconv cells acquire unique functionality during post-PTCy aGVHD.

Accordingly, we pursued metabolic flow cytometric studies⁶⁰ to evaluate key glycolysis-related transporters (GLUT1 and HK2), mitochondrial transport proteins (TOMM20, necessary for import and assembly of respiratory chain complexes⁶¹), and mitochondrial porins (VDAC, linking glycolysis and OXPHOS⁶²). Consistent with transcriptional studies, we observed an enhanced OXPHOS signature with TOMM20 and VDAC expression in Tconv (less so in CD8⁺ T cells and Treg subsets) in the aGVHD cohort and an overall augmented metabolic profile when compared with patients without GVHD and donor controls (Figure 5B). These same patients demonstrated a trend to enhanced glycolysis on day 28, although with heterogenous effects across the lymphocyte subsets (Figure 5C).

Having probed the predictors of aGVHD after PTCy use, we next focused on detecting those that impact OS and identified NK cell recovery as the main determinant. To give context, median numbers of NK cells in HLA-matched and haplo alloBMT cohorts at day 28 were similar at 54 cells/µL and 66 cells/µL, respectively. Using CART, NK cell count ≤50.5 cells/µL at day 28 emerged as the only significant risk marker for OS (Figure 6A). Two-year OS and progression-free survival (PFS) for recipients with NK cells >50.5 cells/µL were 81% (95% CI: 72% to 90%) and 76% (95% CI: 67% to 86%) compared with 50% (95% CI: 39% to 64%) and 44% (95% CI: 34% to 59%) for those with NK cells ≤50.5 cells/µL (P < .00001 and P ≤ .0001, respectively) (Figure 6B). We then performed cumulative incidence curves by NK cell counts >50.5 cells/µL (Figure 6B) and found that both relapse (21% compared with 35% at 2 years, P = .04) and NRM (4% compared with 21%, P = .004) were significantly lower in those with higher NK cell counts. Next, we performed VIMP via random forest analysis to identify additional predictors of OS (Figure 6C). NK cell count was again the top selected factor. Other selected factors in order of the strength of association included pre-BMT disease status, ST2, Naïve CD4⁺ T cells, TNFR1, graft CD3⁺ cell dose, CXCL9, Reg3α, alloBMT year, and plasmablast counts.

We next sought to evaluate the phenotype and functional hallmarks of reconstituting NKs after PTCy-based alloBMT. Initial characterization of the recovering NK cells in haplo and HLA-matched alloBMT recipients confirmed that the most immature population of NK cells (CD56^brightCD16^neg) predominated at day 28 in all cohorts (supplemental Figure 5A), which was consistent with our previous findings and those of Rambaldi et al in different haplo cohorts using PTCy.^20,21 CD56^brightCD16^int NK cells comprised a higher percentage of the NK cells at all time points after alloBMT when compared with healthy donors. Recovery of CD56^dimCD16^pos NK cells in both donor cohorts was slow and only reached the normal distribution at year 1 post-BMT (supplemental Figure 5B). The percentage of CD56^negCD16^pos was similar to that of normal donors throughout all time points.

Prior studies have shown that loss of NK functionality (decreased production of TNF-α and IFN-γ) early after alloBMT predicts enhanced relapse.⁶³ To test if a similar phenotype occurred with post-PTCy relapse, we focused on the haplo alloBMT cohort and sort-purified NK cells from PBMCs obtained on day 28 after alloBMT (Figure 7A,B). Relapse was associated with profound alterations in multiple inflammatory signaling pathways, with key impacts on TNF-α signaling (Figure 7C). NK cells retrieved from relapsing patients exhibited transcriptional aberrancies, distinguishing them from those identified in naïve and functional NK cells⁶⁴ (Supplemental Figure 6). In addition, we profiled CD8⁺ T cells isolated from the same patients and identified enrichment for a molecular signature of exhaustion, with the loss of naïve CD8⁺ T cell hallmarks (Figure 7). Taken together, these data suggest that disease recurrence may be heralded by both the numerical disadvantage of NK cell reconstitution and transcriptional hallmarks of NK dysfunction and CD8⁺ T cell exhaustion.

As the first study to directly compare immune reconstitution after PTCy-based MSD, MUD, and haplo alloBMT, we provide critical insights into the patterns of immune recovery and the landscape of PTCy-induced tolerance. First, we found comparable quantitative reconstitution of NK and B cells after all platforms, but a transient early delay in CD4⁺ and CD8⁺ T cell recovery in the haplo alloBMT cohort that may be attributable to their additional immunosuppression and abrogated by discontinuation of MMF at day 35. Rapid EM and TEMRA CD4⁺ and CD8⁺ T cell recovery supports the low rates of PTLD and infectious deaths after PTCy.^13,65

Our work also highlights a distinct landscape of tolerance that arises with PTCy use. Enrichment of bonafide effector Treg signatures and preserved Treg numbers argue against the loss of Treg-driven tolerance as the key mechanism of post-PTCy breakthrough aGVHD. We hypothesize that the expansion and functionality of Tregs during an aGVHD event after PTCy may be compensatory, contribute to steroid responsiveness, and prevent the development of severe aGVHD or cGVHD, although this requires testing in future studies. Tconv enrichment for activation pathways, accompanied by lack of such enrichment in CD8⁺ T cells, adds to the growing evidence of the key role of CD4⁺ T cells in alloimmunity.^66,67 Further, we identify Tconv recovery, when coupled with elevated CXCL9, a CXCR3 ligand, and type 1 interferon-regulated chemokine that promotes leukocyte tissue egress, as a critical feature in aGVHD development after PTCy. Transcriptional and proteomic upregulation of CXCR3 on Tconv in both patients and mice further underscores the importance of this partnership in post-PTCy aGVHD. Whether the CXCL9-CXCR3-Tconv axis can be exploited as a prognostic or predictive biomarker or as a therapeutic target remains to be tested in future studies.

While PTCy-driven aGVHD protection was attributed to CD8⁺ T cell hyporesponsiveness in murine models,¹⁸ our data suggest that the role for CD8⁺ T cells may be more limited in post-PTCy aGVHD. In our cohort, the classical hallmarks of T cell exhaustion were not associated with protection from GVHD, but may be relevant to the risk of disease recurrence, a finding that requires further testing. Since PTCy has been associated with improved survival and less aGVHD in patients exposed to immune checkpoint blockade pre-BMT,^68,69 our data may support the use of checkpoint blockade after PTCy-based alloBMT to augment graft-versus-leukemia responses without an augmented risk of aGVHD.

Our immunophenotypic studies also yield insight into relapse risk after PTCy allografting, uncovering decreased early NK counts as a key predictor of reduced OS due to increases in both relapse and NRM. These data are augmented by evidence of impaired inflammatory functioning in NK cells and exhausted CD8⁺ T cells in patients who ultimately relapsed. The observed global dysfunction included a multitude of cytokine-signaling driven pathways and acquisition of dysfunctional NK phenotype. We also demonstrated that gene expression in NK cells on day 28 in patients who did not relapse was consistent with the functional capacity reported by Pical-Izard et al.⁶³ These observations lend further credence to the pursuit of relapse prevention with the transfer of expanded mature NK cells⁷⁰ and provide directions for future biomarker-driven adaptive therapy instituted early after BMT. Shortcomings of our study include the relatively small numbers of patients available for bulk immune subset sequencing, emphasizing the need for future transcriptome analyses in aGVHD, including single-cell sequencing of effector T cells in target tissues. Furthermore, future studies should also seek to independently validate our machine-learning findings.

In summary, we identify the limited impact of HLA-matching on the robust and rapid immune reconstitution patterns following PTCy-based myeloablative conditioning alloBMT. Slight differences in the timing of T cell recovery are best explained by the use of additional immunosuppression and younger median recipient age in the haplo cohort. While our data does not allow for a direct comparison with allografting using conventional GVHD prophylaxis, a recent study lends further credence to our emphasis on CD4⁺ T cell dysregulation during post-PTCy GVHD.⁷¹ We also provide molecular insights into the putative operational mechanisms of PTCy by coupling transcriptional observations with the results of machine-learning studies to identify Tconv as key mediators of aGVHD and NK numbers and functional signatures as critical to underlying disease control. Our work highlights the unique immune landscape after PTCy and hints at GVHD mechanisms distinct from those reported with calcineurin-based allografting. Whether the same paradigms apply to the nonmyeloablative or PBMC-based allografting strategies, where G-CSF mobilization-related imprinting creates distinct immune circuitry, remains to be explored. Given that immune reconstitution is a dynamic process, selective use of day 28 studies to feed the machine-learning algorithms may be another limitation. However, this served well our goal of identifying early posttransplant changes that may be exploited for development of adaptive strategies to prevent or treat GVHD and to optimize graft-versus-tumor effects. For instance, targeting the CXCR3-CXCL9 axis may prove beneficial in aGVHD control or prevention. In addition, harnessing NK cells to minimize relapse will be explored in our ensuing work.

The authors thank Edus H. Warren, Borje S. Andersson, and Paul V. O’Donnell for contribution to the original NCT00809276 studies and Sofia Berglund for NK phenotypic analysis.

This work was supported by grants from the NIH/NHLBI (R01HL110907) and Otsuka Pharmaceutical to L.L., V.R. (K08HL145116), J.D.P. (R01HL110907), and NIH/NCI to S.P. (R01CA168814) and to L.L., R.J.J., and K.R.C. (NIH/NCI P01CA225618, P30CA006973). Career Development Award from the American Society for Transplantation and Cellular Therapy (V.R.) and Huntsman Cancer Foundation (V.R.). Work presented in this publication was made possible through the use of the NIH/NCI P01CA225618 Immune Monitoring Core (L.L.); JH-TIE Shared resource (NIH/NIBIB P41EB028239); University of Utah and the Huntsman Cancer Institute Shared Resource Services (NIH/NCI P30CA042014; NIH/NCI 1S10RR026802-01).

Contribution: S.R.M., V.R., and L.L. were involved in the study’s conception and design, data collection, analysis, laboratory analyses and interpretation, primary writing of the manuscript, and manuscript revisions; H.L.T. and R.V. statistically analyzed and interpreted the data and revised the manuscript; C.G.K., A.V., S.I., E.T., B.L., and D.A. were involved in data collection and laboratory analyses; and C.G.K., S.P., K.R.C., R.J.J., J.D.P., and H.J.S. were also involved in analysis, conception and design, and revised the manuscript; and all authors approved the manuscript.

Conflict-of-interest disclosure: S.R.M. Grant/research/clinical trial support: Gilead Sciences, Aprea Therapeutics, McCabe Fund. V.R. Grant/research/clinical trial support: Syndax Pharmaceuticals. Consultant/advisory board: Regeneron Pharmaceuticals. S.P. Patent holder on “Biomarkers and assays to detect chronic graft versus host disease” (U.S. patent 10 571,478 B2). L.L. Patent holder: WindMiL Therapeutics. Grant/research/clinical trial support: Genentech. Consultant/advisory boards: Gilead Sciences, Rubius Therapeutics, PrecisionBiosciences, Talaris Therapeutics. The remaining authors declare no competing financial interests.

Correspondence: Leo Luznik, Johns Hopkins University, The Bunting-Blaustein Cancer Research Building, 1650 Orleans St Room 2M08, Baltimore, MD 21287; e-mail: luznile@jhmi.edu; and Vedran Radojcic, University of Utah Huntsman Cancer Institute, 2000 Circle of Hope Dr, Rm 5263, Salt Lake City, UT 84112; e-mail: vedran.radojcic@hsc.utah.edu.