Recipient and donor genetic variants associated with mortality after allogeneic hematopoietic cell transplantation

Many studies have suggested that genetic variants in donors and recipients are associated with the overall survival (OS), progression-free survival (PFS), nonrelapse mortality (NRM), or relapse mortality (RM) of patients after allogeneic hematopoietic cell transplantation (HCT). The associations identified in these studies, however, were not confirmed in a large independent HCT cohort of unrelated donors and recipients from the Center for International Blood and Marrow Transplant Research (CIBMTR).^1,2  The CIBMTR studies of candidate single-nucleotide polymorphisms (SNPs)¹  and genome-wide exome polymorphisms²  were limited to patients with unrelated donors, and the genome-wide testing did not include nonexonic polymorphisms. Moreover, our review of the updated and previously published literature identified 73 variants that were not included in the CIBMTR studies. Consequently, the utility of testing genetic variation in donors and recipients for risk stratification or understanding mechanisms leading to mortality after HCT has not been established.

We took a twofold approach to address gaps in the CIBMTR studies. First, from a review of the literature through December 2019, we compiled a list of 130 variants, associated with 79 genes, reported as individual SNPs, SNP haplotypes, SNP combinations, or indels having statistically significant associations with mortality after HCT. Forty-seven of these variants were included in the CIBMTR analysis,¹  which encompassed reports published through 2016. Sixty-eight other variants were identified in reports published before 2017 but were not included in the CIBMTR analysis, and 15 were identified in reports published after 2016. We tested 122 of these candidate recipient and donor variants for associations with NRM (ie, death before recurrent or progressive malignancy) and RM (ie, death after recurrent or progressive malignancy) in a cohort of 2560 HCT recipients of European ancestry with related or unrelated donors at our center. Associations discovered in this cohort were tested for replication in a separate cohort of 1710 HCT recipients. Second, we conducted a genome-wide association study (GWAS) to identify variants associated with NRM and RM using the same discovery and replication approach.

All recipient and donor blood samples were collected before HCT according to research protocols approved by the Fred Hutchinson Cancer Research Center (FHCRC) Institutional Review Board. Project-specific institutional review board approval was obtained for the use of clinical data and research biospecimens.

Table 1 summarizes demographic, clinical, and transplantation characteristics of patients in the study cohort. The overall study cohort included 4270 recipients of European ancestry who had a first allogeneic HCT at the FHCRC or Seattle Cancer Care Alliance from 1990 through 2011. Syngeneic or cord blood donors were excluded. A single prior autologous HCT was allowed. The study was limited to recipients of European ancestry because the number of available non-European recipients was too small for a meaningful analysis in this study. European ancestry was defined by using the minimum covariant determinant method as implemented by Conomos et al.³  HLA-A, -B, -C, -DRB1, or -DQB1 mismatching was present in 7% of the related donor-recipient pairs and in 37% of the unrelated pairs. Indications for HCT included hematological malignancy or myelodysplasia. Conditioning regimens were categorized as myeloablative or nonmyeloablative according to the intensity of chemotherapy and total-body irradiation. Grafts used for HCT included bone marrow cells or growth factor–mobilized blood cells. Supplemental Figure 1 shows the cumulative incidence frequencies of NRM and RM in the study cohort.

We performed a PubMed search to identify studies that reported a recipient or donor SNP, SNP haplotype, or SNP combination associated with OS, NRM, or PFS after allogeneic HCT at a statistical significance level of P ≤ .05. Other genetic variants such as tandem repeats associated with these outcomes were excluded because the genotyping arrays used for our study are not informative for such variants. Supplemental Table 1 shows the candidate variants and reported phenotypes, with supporting references.

Details regarding sample preparation have been described previously.⁴  Donor samples were not available for 214 pairs (2%), and recipient samples were not available for 370 pairs (4%). Genomic DNA was extracted from blood mononuclear cells or Epstein-Barr virus–transformed B-lymphocyte cell lines using a Puregene kit (Qiagen, Valencia, CA). Genotyping assays used 3 different platforms: the Affymetrix 5.0 Human GeneChip (cohort 1; 1048 recipients and 1029 donors), the Illumina 1M Quad (cohort 2ab; 1904 recipients and 1828 donors), and the Illumina 2.5M BeadArray (cohort 2c; 688 recipients and 926 donors). Amplification and hybridization for the Affymetrix 5.0 array were performed at the Affymetrix Service Laboratory (Santa Clara, CA), and amplification and hybridization for the Illumina BeadArrays were performed by the FHCRC Genomics Shared Resource Laboratory. The genotypes of the candidate variants were determined separately for each platform using the BRLMM algorithm for the Affymetrix array⁵  and the GeneCall algorithm for the Illumina arrays.⁶ 

Quality assurance and quality control (QC) were primarily performed separately for variants and samples from each cohort following standard methods.⁷  Briefly, for each cohort, we evaluated batch quality (groups of DNA samples processed together), individual sample quality, ancestry differences, and missing data. We identified large chromosomal anomalies (duplications, deletions, and acquired uniparental disomy) that may cause errors affecting genotype calls relative to the germ line. For the variants in each cohort, we evaluated the missing call rate, Hardy-Weinberg equilibrium, Mendelian errors, heterozygosity, and allele frequency differences between males and females and discordant genotypes between duplicate sample pairs. We used a common subset of variants genotyped on each platform to estimate the pairwise relatedness of all samples and identify samples with questionable identity, expected and unexpected duplicate samples, and incorrectly annotated relationships. Samples with questionable identity were removed, and a fully consistent pedigree was constructed for all related samples. This effort was further facilitated by using the SNP genotypes to impute the HLA type of each sample for comparison with the clinical typing records.

After QC, cohort 1 had 2813 samples and 390 412 variants, cohort 2ab had 4072 samples and 1 003 751 variants, and cohort 2c had 2215 samples and 2 327 361 variants. The resulting post-QC sample set contained a total of 3850 recipients and 3977 donors.

For each sample, we calculated a uniform set of 4-digit allele assignments for HLA-A, -B, -C, -DRB1, and -DQB1 using HIBAG software⁸  to impute the HLA alleles from the SNP genotypes. In general, the imputed alleles had high concordance with the prior genotyping methods used for the clinical typing. A human expert reviewed and resolved all discrepancies. This approach also resolved ambiguous HLA codes, typed previously unidentified alleles, and provided a uniform set of 4-digit allele assignments.

We used SHAPEIT2⁹  and IMPUTE2¹⁰  software to separately prephase and impute, respectively, the samples in each cohort to a common set of variants from the 1000 Genomes Project phase 3.¹¹  We also imputed a small number of sporadic missing genotype calls. We limited the list of imputed variants to those with at least 4 minor allele copies in any of the 5 1000 Genomes Project superpopulations (AFR, AMR, EAS, EUR, SAS). The posterior probability of the most likely genotype was calculated as the probability of observing an unobserved genotype at the imputed locus, given all observed genotypes in the flanking region. The imputed SNP genotype was retained only if the posterior probability of the most likely genotype exceeded 0.9, and it was otherwise set to missing.

To harmonize the genotyped and imputed variants and samples from all 3 genotyping platforms into a single integrated data set, we ensured, for each variant, that the genomic location (chromosome and position) and plus-strand alleles¹²  matched across platforms and aligned to genome build GRCh37/hg19. We used IMPUTE2 metrics (information score and average posterior probability) to retain only well-imputed variants and our own measures (Mendelian errors between sibling pairs, missing call rate, and duplicate genotype discordance across cohorts) to further remove low-quality variants and/or those exhibiting batch effects between cohorts. The result was a data set with 9100 samples and 29 826 485 post-QC variants: 28 889 566 autosomal and 936 919 on the X chromosome. Of these, 15 282 884 passed genotyping and imputation QC in all 3 cohorts, an additional 2 279 374 passed QC in the 2 Omni cohorts but not in the Affymetrix cohort, and another 12 264 227 passed QC in only 1 of the 3 cohorts. Among the 29 826 485 autosomal and X chromosomal post-QC variants, 8 832 649 had minor allele frequencies (MAFs) >0.01. Of these, 8 010 049 were SNPs and 822 600 were indels; 180 017 were in an exon, and 28 535 were annotated as encoding a change in amino acid sequence.

Outcomes tested in this study were NRM and RM. Recurrent or progressive malignancy after HCT was a competing risk for NRM, and NRM was a competing risk for RM. Follow-up of recipients without NRM or RM was censored at last contact and ranged from 1.3 to 25.6 years after HCT (median, 12.9 years). Evaluation of association of genotype with outcomes was based on cause-specific HR analysis using Cox regression. All candidate variants in the recipient and donor genomes were evaluated for allelic and genotypic (recessive and dominant) associations with NRM and RM. For a SNP with a major allele a and a minor allele b, the recessive model tested the hypothesis that the genotype bb would be associated with a higher or lower risk compared with the collective genotypes ab and aa used as the reference. The dominant model tested the hypothesis that the collective genotypes bb and ab would be associated with a higher or lower risk compared with the genotype aa used as the reference. The allelic model tested the hypothesis that the minor allele b would be associated with a higher or lower risk compared with the major allele a, and the number of copies of the minor allele was modeled as an additive effect. The GWAS included both the recipient and donor genomes but was limited to the allelic model. For all SNPs, the minor allele was defined as the less frequent allele in the combined recipient and donor, discovery, and replication data set.

All analyses of candidate variants were conducted in 2 phases. The discovery phase included 60% of the subjects in the data set, and the replication phase included the remaining 40% of the subjects. Candidate variants with P ≤ 5 × 10⁻³ for association with NRM or RM were tested in the replication phase with Bonferroni adjustments for multiple comparisons. In the GWAS discovery phase, variants with MAF >0.01 and P ≤ 1.0 × 10⁻⁶ for association with NRM or RM were tested in the replication phase with Bonferroni adjustments for multiple comparisons. Bonferroni adjustment was applied separately for the 4 genome–end point combinations. Post hoc power estimates to detect an HR ≥1.5 or ≤0.67 at a 2-sided significance level of .05 were based on the estimated standard error of the log HR in the replication cohort.

Supplemental Table 1 summarizes published reports describing 130 variants having statistically significant associations with OS, PFS, or NRM after allogeneic HCT, and supplemental Table 2 summarizes the QC assessment of the candidate variants genotyped or imputed on the 3 platforms used for our study. Six SNPs were not analyzed: 3 did not pass QC (rs333, rs3969913, and rs4715333), 2 were monomorphic (rs4986893 and rs671), and 1 was not genotyped or imputed (rs396991). We were also unable to analyze the H1-homozygous CCR5 haplotype reported by McDermott et al¹³  or the CYP2C19 combination reported by Elmaagacli et al.¹⁴ 

Supplemental Table 3 shows results of testing candidate variants for association with NRM. In the discovery cohort, 2 recipient variants met criteria for replication testing (Table 2), but results were not statistically significant in the replication cohort. Post hoc power was 100% in testing rs10975123 in the allelic model. The HR point estimates for rs10975123 in the recessive model and rs429916 were notably on opposite sides of 1.0 in the discovery and replication cohorts. Taken together, these results suggest that the lack of replication did not result from insufficient power.

In the discovery cohort, 2 donor variants met criteria for replication testing (Table 2). Results for rs3212227 were not statistically significant in the replication cohort, and here again, the HR point estimates were on opposite sides of 1.0. In contrast, results for rs1051792 (MICA) met the Bonferroni-adjusted 0.025 threshold of statistical significance in the replication cohort. The respective HR point estimates in the discovery and replication cohorts were 0.84 and 0.85, respectively, for association of the donor rs1051792 A allele with NRM (Table 2).

Supplemental Table 4 shows results of testing candidate variants for association with RM. In the discovery cohort, only 1 donor SNP met criteria for replication testing (Table 2), but results were not statistically significant in the replication cohort. The HR point estimates were similar, and post hoc power was 100%, indicating that the lack of replication did not result from insufficient power.

In the discovery phase, the GWAS of NRM yielded 6 candidate variants that met criteria for replication testing (Figure 1). In people of European ancestry, rs9350084 and rs9350085 are in close linkage disequilibrium (LD) on chromosome 6 (r² = 1.0), and rs11559982, rs10876550, and rs11451044 are in close LD on chromosome 9 (r² ≥ 0.98; supplemental Tables 5-6). Within each LD group, we selected the variant with the lowest P value for replication testing (Table 3). Replication testing thereby encompassed 3 variants, 1 from each of the 2 LD groups (recipient rs9350085 and donor rs11451044, respectively) plus donor rs9492413. Power for testing in the replication phase was 100% for all 3 variants, but none of the associations with NRM was statistically significant (P = .65, 1.0, and .86, respectively).

In the discovery phase, the GWAS of RM yielded 6 candidate variants that met criteria for replication testing, but rs148290359 and rs74782705 are in LD on chromosome 11 (r² = 0.93), and rs141251357 and rs60125430 are in LD on chromosome 9 (I² = 0.75; Figure 1). In replication testing, none of the 4 variants that remained after pruning for LD had a statistically significant association with RM, although power was limited for all but 1 of these variants. Supplemental Tables 7 and 8 and supplemental Figure 2 summarize the combined discovery and replication GWAS analysis of variants with MAF >1% and P ≤ 1 × 10⁻⁶ for association with NRM or RM with the allelic model for others to test for replication.

The association of the donor rs1051792 A allele with NRM remained unchanged after adjustment for clinical covariates, including donor relation, HLA mismatch, female-to-male sex mismatch, recipient and donor age, source of stem cells, and conditioning intensity (HR, 0.85; 95% CI, 0.74-0.98; P = .02; Table 4). Isernhagen et al¹⁵  previously reported that the rs1051792 A allele in the recipient was associated with lower overall mortality after allogeneic HCT, whereas our discovery and replication identified an association of NRM with the rs1051792 A allele in the donor. Because of HLA matching of donors and recipients and the close LD between HLA-B and MICA, the donor and recipient genotypes for rs1051792 are highly correlated. For pairs with both genotypes, the concordance was 96%, making it statistically impossible to distinguish donor vs recipient associations. The recipient rs1051792 A allele showed an association with NRM similar to that of the donor allele in the replication cohort, although this did not reach statistical significance in the somewhat different set of patients with genotype data (HR, 0.88; 95% CI, 0.8-0.1.0; P = .07; Figure 2A-B).

Previous reports have indicated that MICA mismatching can contribute to NRM.^16,17  Two observations suggest that mismatching cannot explain the association of donor rs1051792 genotypes with NRM in our study. First, the association of the donor rs1051792 genotypes with NRM was similar in HLA genotypically identical sibling recipients who had no MICA mismatching and in unrelated recipients who could have had MICA mismatching. The respective HR estimates were 0.82 and 0.87. Second, 27 of the 1203 HLA-A, -B, -C, -DRB1, or -DQB1 matched unrelated donor-recipient pairs had graft-versus-host or host-versus-graft rs1051792 allele mismatching, which was not associated with NRM (HR, 0.91; 95% CI, 0.43-1.93; P = .81). Similarly, 27 of the 67 unrelated pairs with a single recipient HLA-B mismatch had rs1051792 mismatching, which was also not associated with NRM (HR, 0.51; 95% CI, 0.23-1.13; P = .10). Taken together, these results support the hypothesis that the association of rs1051792 with NRM is attributable to an intrinsic protective effect of the minor allele.

Because the association with NRM might be explained by differences in the risk of acute or chronic graft-versus-host disease (GVHD), we evaluated the incidence of grade 2 to 4 acute GVHD, grade 3 to 4 acute GVHD, NRM among patients with grade 3 to 4 GVHD, and chronic GVHD according to the donor and recipient rs1051792 genotypes in the combined discovery and replication cohorts (Table 5). The results showed that the donor rs1051792 A allele was not associated with the incidence of grade 2 to 4 acute GVHD but was marginally associated with a decreased risk of grade 3 to 4 acute GVHD (HR, 0.90; 95% CI, 0.81-1.01; P = .05; Figure 2C) and a decreased risk of NRM in patients with grade 3 to 4 GVHD (HR, 0.88; 95% CI 0.76-1.01; P = .06; Figure 2E). The donor rs1051792 A allele was not associated with the risk of chronic GVHD or recurrent malignancy. The recipient rs1051792 A allele showed no statistically significant association with any of these end points (Table 5; Figure 2D,F). Taken together, these results suggest that the lower risk of NRM associated with the donor rs1051792 A allele may be explained in part both by a lower risk of grade 3 to 4 GVHD and a lower risk of mortality among patients with grade 3 to 4 acute GVHD.

With 1 exception, broadly applicable associations with NRM or RM could not be demonstrated for previously studied donor and recipient genetic variants. The 1 exception was rs1051792, the minor allele A of which encodes methionine and major allele G of which encodes valine at position 129 in the MICA protein. No previous study has demonstrated an association of donor MICA-129Met with NRM after allogeneic HCT.^16,18  Fuerst et al¹⁶  showed an association of donor-recipient MICA-129 mismatching with NRM after unrelated HCT, but this result was not replicated in other studies^17,18  or in our study, although MICA mismatching at large has been associated with increased risks of NRM and acute and chronic GVHD and a decreased risk of recurrent malignancy after HLA-A, -B, -C, -DRB1, or -DQB1 matched unrelated HCT.¹⁷  Isernhagen et al¹⁵  reported that the presence of MICA-129Met in the recipient was associated with improved survival after allogeneic HCT, particularly in patients who had not been treated with antithymocyte globulin as part of the pretransplantation conditioning regimen. The hazard of overall mortality associated with each additional rs1051792 A allele in the recipient was 0.77 (95% CI, 0.6-1.0; P = .04), but the association with the hazard of NRM was not statistically significant (HR, 0.51; 95% CI, 0.2-1.1; P = .09). In addition, the recipient MICA-129Met variant was associated with an increased frequency of grade 2 to 4 acute GVHD but lower risk of fatality in patients with acute GVHD. MICA-129Met has high avidity for binding to its ligand NKG2D, with more efficient initial NKG2D signaling, but the high-avidity binding reduces NKG2D expression and function in NK cells and CD8⁺ T cells over time, whereas the lower-avidity binding of MICA-129Val does not have this effect.^15,19  Isernhagen et al¹⁵  suggested that the lower NRM associated with MICA-129Met in the recipient could be explained by the attenuation of NKG2D expression and function over time after HCT.

Our results showing decreased risks of grade 3 to 4 GVHD and NRM in patients with grade 3 to 4 GVHD associated with donor MICA-129Met but not with MICA-129Val implicate MICA-129Met in the donor as an explanation for the decreased risk of NRM after HCT. Although MICA expression has been widely recognized as a marker of epithelial stress, it is also expressed by activated T cells and monocytes. Cerboni et al²⁰  showed that antigen-activated human T cells express MICA and other NKG2D ligands, thereby becoming susceptible to fratricidal killing by autologous NK cells. Subsequent studies showed that donor NK cell–mediated killing of alloactivated donor T cells attenuates the severity of acute²¹  and chronic²²  GVHD in mice without impairing graft-versus-tumor effects,²¹  although in some circumstances, exogenous hyperactivation can cause NK cells to produce proinflammatory cytokines that sustain induction of GVHD by T cells.²³  Based on these results, we hypothesize that alloactivated donor T cells with MICA-129Met have increased susceptibility to killing by donor NK cells or CD8⁺ T cells²⁴  as a result of stronger binding to NKG2D when compared with donor T cells with MICA-129Val. Our results do not refute effects mediated by recipient MICA-129Met. Instead, they highlight a previously unrecognized mechanism that could explain effects mediated by donor MICA-129Met. Although linkage of MICA with the major histocompatibility complex would not allow clinical application of our results through donor selection, it is possible that addition of donor NK cells to the graft could attenuate the severity of acute GVHD.^21-23 

Our study and the CIBMTR studies^1,2  had similarities and many differences but nonetheless reached the same general conclusion regarding the association of genetic variants with mortality-related end points after HCT. Similarities included the selection of candidate variants from a comprehensive screen of the literature reporting statistically significant SNP associations with OS, PFS, NRM, and RM, the focus on recipients of European ancestry, the use of myeloablative and reduced-intensity conditioning regimens, and the use of marrow or growth factor–mobilized blood cells for grafting.

At the same time, our approach differed from the CIBMTR approach in many respects, as summarized in Table 6. The CIBMTR candidate SNP analysis included 47 variants statistically associated with survival outcomes reported in previous publications, whereas our analysis included more than twice that number of variants (n = 122). The Karaesmen et al¹  study drew a distinction between replication and validation.²⁵  Replication testing requires matching the same demographic, disease, and donor type inclusion and exclusion criteria and the end points used for each prior study, whereas validation allows differences to determine whether results apply more broadly. In our study, prior reports were used solely to identify candidate variants, and we made no attempt to match the inclusion and exclusion criteria used in prior studies. The scope of the CIBMTR GWAS was limited to exons, whereas our GWAS encompassed all variants that passed QC and had an MAF >1%. The only statistically significant result replicated or validated in the CIBMTR cohort was the association of the donor IL6 rs1800795 C allele with OS (HR, 1.11; 95% CI, 1.0-1.2; P = .018 in the allelic model).¹  In our cohort, however, rs1800795 donor and recipient genotypes showed no statistically significant association with the risk of NRM or RM in recessive, dominant, or allelic models (supplemental Tables 3 and 4).

The negative results of our study and other studies^1,26-28  give rise to the question of whether NRM or RM can be considered as a true phenotype in HCT recipients. With respect to NRM, HCT recipients have variable vulnerabilities related to age and the residual morbidity caused by the underlying disease and prior treatment, and they are subjected to variable stressors such as pretransplantation conditioning regimen, GVHD, and infections. Vulnerability to RM depends on the nature of the underlying disease, the extent to which prior treatment selected for malignant cells that are resistant to the pretransplantation conditioning regimen, and the burden of disease at the time of HCT, whereas higher-intensity conditioning regimens and graft-versus-leukemia effects decrease the risk of recurrent or progressive malignancy after HCT. As a result, death after HCT has many causes. Accordingly, 1 SNP could affect NRM or RM only if it has wide-ranging effects across many pathophysiological mechanisms or such a large effect on any single pathophysiological mechanism that an association with NRM or RM could be detectable with the size of the available HCT cohorts that have extensive genotype data available. Results showing that previously reported genetic associations do not apply broadly should discourage additional studies treating death as a phenotype. Future studies should focus instead on associations with the individual intermediate phenotypes that precede death.

The authors thank Veronika Groh for assistance with the interpretation of our results.

This work was supported by National Institutes of Health grants from the National Institute of Allergy and Infectious Diseases (AI33484 and AI149213), National Cancer Institute (CA015704 and CA18029), and National Heart, Lung, and Blood Institute (HL087690, HL088201, HL094260, HL105914, and K23HL69860).

Contribution: J.A.H., P.J.M., D.M.L., and B.E.S. were responsible for the study concept and design, interpretation of results, and writing of the manuscript; S.C.N. assisted with imputation and genomic phasing; X.D. assisted in curation of candidate variants and GWAS analyses; D.M.L. and B.E.S. were responsible for data acquisition, imputation, informatics analyses, quality control, and statistical analyses; and P.J.M., D.M.L., B.E.S., S.C.N., and X.D. revised draft manuscripts and approved the final version of the manuscript.

Conflict-of-interest disclosure: The authors declare no competing financial interests.

John A. Hansen died on 31 July 2019.

Correspondence: Paul J. Martin, Fred Hutchinson Cancer Research Center, 1100 Fairview Ave N, Seattle, WA 98109-1024; e-mail: pmartin@fredhutch.org.