Association study of candidate DNA-repair gene variants and acute graft versus host disease in pediatric patients receiving allogeneic hematopoietic stem-cell transplantation

Acute Graft versus Host Disease (aGvHD) grades 2–4 occurs in 15–60% of pediatric patients undergoing allogeneic haematopoietic stem-cell transplantation (allo-HSCT). The collateral damage to normal tissue by conditioning regimens administered prior to allo-HSCT serve as an initial trigger for aGvHD. DNA-repair mechanisms may play an important role in mitigating this initial damage, and so the variants in corresponding DNA-repair protein-coding genes via affecting their quantity and/or function. We explored 51 variants within 17 DNA-repair genes for their association with aGvHD grades 2–4 in 60 pediatric patients. The cumulative incidence of aGvHD 2–4 was 12% (n = 7) in the exploratory cohort. MGMT rs10764881 (G>A) and EXO rs9350 (c.2270C>T) variants were associated with aGvHD 2–4 [Odds ratios = 14.8 (0 events out of 40 in rs10764881 GG group) and 11.5 (95% CI: 2.3–191.8), respectively, multiple testing corrected p ≤ 0.001]. Upon evaluation in an extended cohort (n = 182) with an incidence of aGvHD 2–4 of 22% (n = 40), only MGMT rs10764881 (G>A) remained significant (adjusted HR = 2.05 [95% CI: 1.06–3.94]; p = 0.03) in the presence of other clinical risk factors. Higher MGMT expression was seen in GG carriers for rs10764881 and was associated with higher IC50 of Busulfan in lymphoblastoid cells. MGMT rs10764881 carrier status could predict aGvHD occurrence in pediatric patients undergoing allo-HSCT.


INTRODUCTION
The most frequent immunological complication after allogeneic haematopoietic stem-cell transplantation (allo-HSCT) is acute Graft versus Host Disease (aGvHD), occurring in 15-60% of transplants in children [1]. In patients receiving HSCT from HLA-identical siblings, the long-term survival rates with aGvHD grades III-IV are below 30% [2]. aGvHD begins with host normal tissue damage by the conditioning regimens that causes pro-and anti-inflammatory cytokine secretion which subsequently activate the host antigen presenting cells (this being phase 1 of the pathobiology of aGvHD) [3]. Thus, the intensity or type of conditioning regimen is determined as one of the donor-independent risk factors for aGvHD [4,5].
Busulfan (BU) is frequently used for conditioning children prior to allo-HSCT [6]. BU is a bifunctional alkylating agent (AG) commonly administered with other alkylating agents like cyclophosphamide (CY) or the purine analog Fludarabine (FLU) [7]. BU and other AGs mediate their cytotoxicity by damaging the DNA through formation of covalent linkages between the alkyl groups, mainly the N 7 position of guanine, while the N 3 position of cytidine and O 6 of guanine also serve as nucleophiles [8]. These covalent modifications lead to inter-or intra-strand DNA crosslink formation, which affects the genomic integrity and causes deleterious consequences during DNA replication. That effect is observed in tumor cells but also in normal cells, the latter being linked to the treatment-related toxicities (TRTs) such as aGvHD.
Hence, variants related to BU metabolism such as GSTA1*B [9] and GSTM1-null [10] were described as risk factors for aGvHD. Other genetic polymorphisms within immunological pathways were also described as risk factors [11], some of them in a pediatric population [12].
The DNA damage caused by AG is repaired by various DNArepair pathways [8,13] including base excision repair (BER), mismatch repair (MMR) and homologous recombination or by nonhomologous end joining. Other mechanisms include demethylation of guanine residues by O6-methylguanine-DNA methyltransferase (MGMT) and MMR of small insertions and modifications by identifying the damaged base with the help of Mut L homologue-1 protein [14]. Genomic predictors of interindividual differences in response to DNA damaging agents have previously been demonstrated [15]. Studies have implicated the role of genetic variants and altered expression of genes in the DNA-repair pathways particularly BER and MMR in determining treatment outcomes of AGs [16]. However, less is known about their role in determining clinical outcomes of BU-based conditioning in a pediatric allo-HSCT.
We hypothesized that children receiving allo-HSCT with efficient DNA-repair ability are at reduced risk of developing aGvHD by diminishing the activation of Phase 1 of the pathophysiology of aGvHD. As the DNA damage caused by cross-linking agents like BU is complex and may involve one or more of the abovementioned pathways; in this study candidate genes from key pathways were therefore investigated as possible biomarkers for aGvHD. Selected variants (list of the candidates selected and criteria for the selection of variants are provided in the methods section) among the genes coding for key proteins of demethylating repair pathways, BER pathway genes or doublestrand break repair pathways were selected for the investigation [8,13,14,16].

RESULTS
The characteristics of the study subjects in the exploratory cohort (n = 60) and extended cohort (n = 187) are given in Table 1. The incidence of aGvHD 2-4 was 12% (n = 7) in the exploratory cohort and 22% (n = 40) in the extended cohort.

DNA-repair genetic variants and aGvHD
Hardy-Weinberg equilibrium (HWE) p values and minor allele frequency (MAF) data for each SNP is presented in Table 2. Five SNPs were found to be non-polymorphic in our sample set (ALKBH1 rs17825440; BRCA1 rs28897687; FANCD2 rs9845756; NBN rs1805794; RAD52 rs7310449). Another three SNPs did not pass at the genotyping stage due to unreliable amplification of product (APEX1 rs4585; BRCA1 rs28897687; NBN rs1805800). Forty-three SNPs were carried forward for statistical analyses. Association analyses between genotype (additive or dominant model) against aGvHD 2-4 are illustrated in Fig. 1A. EXO rs9350 and MGMT rs10764881 showed significance (multiple testing p value ≤0.001, with an odds ratio of 11.5 (95% CI: 2.3-191.8) and 14.8 (0 events out of 40). From the extended analysis of these SNPs with aGvHD, only MGMT rs10764881 (p = 0.03) remained significant (Fig. 1B).
Multivariable analysis ( MGMT mRNA expression pre-and post-busulfan exposure MGMT mRNA expression demonstrated no significant change from basal levels (data not shown). Nevertheless, MGMT rs10764881 showed a change in expression levels irrespective of the BU treatment (p = 0.01 pretreatment and 0.03 posttreatment; Supplementary Fig. 1).
Dual luciferase reporter-gene assays. There was no significant differences observed in expression levels between the short plasmid construct containing alleles A and C of variant rs1625649 that is in strong LD (r 2 = 0.85-1.0 in Europeans), with rs10764881, (p = 0.53) (Fig. 2). However, the longer plasmid construct containing variant rs10764881 allele G differed significantly from the plasmid construct containing variant rs10764881 allele A (p = 0.000005) as well as from the shorter fragment without rs10764881 (p = 0.000001) suggesting the presence of an enhancer element near to this SNP, which shows increased dependency with the presence of allele G.

Dexamethasone-mediated activation of hMGMT promoter
Exposing the cells to dexamethasone increased protein expression of luciferase in all the plasmids compared to their non-treated plasmid construct (p < 0.0005). The highest response was seen from the plasmid construct containing variant rs10764881 allele G compared to the other treated constructs (p < 0.0007) (Fig. 2).
Electrophoretic mobility shift assays Electrophoretic mobility shift assays (EMSA) results showed that a shift is established when nuclear protein is added at varying concentrations to the predicted GRE probe (Supplementary Figs. 4,5). The intensity of the shift diminished when unlabeled competitive probes were added at 100 times higher concentration, indicating that the DNA-protein binding occurs with the predicted probe region. Furthermore, interaction appears stronger when allele G is present in the predicted probe.

DISCUSSION
We showed an association of a variant in MGMT (rs10764881, G>A) with aGvHD 2-4 incidence, even after adjusting for other known risk factors (serotherapy, BU AUC and HLA incompatibility). A trend in association of AA and AG genotype carriers at rs10764881 with severity of aGvHD was also seen (no aGvHD, 35% vs. aGvHD grade 1, 55% vs. aGvHD 2-4, 66%). These observations suggest efficient DNA damage repair due to increased MGMT expression and activity in GG carriers. Minor allele ("A") frequencies of rs10764881 among different ethnicities are similar (∼25-30%) except for the African population (<10%) indicating the utility of this genetic marker among non-African ethnicities (Supplementary Table 2).The only other clinical association study that has investigated DNA-repair genes in relation to aGvHD, assessed BER pathway genes and reported a significant association with a variant in RFC1 (rs6844176) in a mixed cohort of adult and pediatric populations [16]. This variant was not significantly associated with aGvHD in our exploratory cohort. Thus, to the best of knowledge, this is the first report that has evaluated candidate DNA-repair gene variants and their association with aGvHD after allo-HSCT in children. This report unfolds the role of DNA-repair pathway gene candidate genes as a biomarker for stratification of patients at a higher risk of developing aGvHD post HSCT. Oligogenic risk score development may include DNA-repair gene variants along with other reported genetic risk factors in immunological mediators for e.g., interleukin 1 [11], interferon-gamma, interleukin 10, and TGF-β [17] or busulfan metabolic pathway [9] to evaluate performance for the prediction of aGvHD risk post HSCT. Understanding the effect of DNA-repair gene variation on normal tissues could also be beneficial for assessing the risk of TRTs [13]. MGMT encodes the DNA-repair protein O 6 -alkylguanine DNA alkyl transferase. It has been studied extensively in association with methylating-agent resistance [18]. Earlier research has mainly focused on evaluating the inhibition of MGMT to augment therapy with alkylating agents. However, Phase I trials showed that inhibitors that inactivate MGMT, improved the efficacy of BCNU (1,3-bis(2-chloroethyl)-1-nitrosourea) on tumor cells but were associated with more systemic toxicity [19], indicating MGMT expression is also vital for protecting from DNA damage in normal cells. BU preferentially induces N 7 , N 3 guanine, and N 3 adenine lesions and MGMT is known to repair lesions of alkylation reactions at oxygen sites, such as O 6 of guanine. No significant change in the IC50 values of BU in the absence and presence of MGMT in HAP1 cells, indicating its limited role in determining the cytotoxicity caused by BU. However, significant differences in IC50s of BU was observed between LCLs carrying different genotypes for rs10764881 in MGMT. GG carriers exhibited higher BU IC50 values compared to AA and AG carriers ( Supplementary Fig. 3). However, this cannot entirely be attributed to rs10764881 genotypes, as several other gene variants at the same time were also associated with changes in BU IC50 in the LCLs. Multiple testing correction was not applied due to the limited number of the samples (n = 58). However, LCLs serve as a good model for investigating the association of potential pharmacogenetic markers. They represent unrelated individuals with the marker of interest. Thus, irrespective of the causative effect, MGMT rs10764881 genotypes could identify cells that are sensitive to BU or resistant to BU defined based on IC50 values. Interestingly, this SNP is not in linkage disequilibrium (LD) with any SNP in the SNPs in the MGMT exonic region, however it is in strong LD with rs1625649 (R2 > 0.85-1.0) in Europeans (except for Finnish population) and Americans, (specifically Mexican ancestry and Peruvian in Lima) and South East Asians (only in Punjabi in Lahore, Pakistan). A better progression free survival was seen in glioblastoma patients carrying rs1625649 "AA" genotypes which was in turn associated with lower MGMT expression [20] and was also shown to reduce the expression of MGMT as a part of investigated promoter haplotypes [21]. Thus, the observed decreased MGMT expression in rs10764881 AG, AA carriers can be partially explained by AA and AC genotypes at rs16265649 locus and as a result of the interaction between these two loci as demonstrated in the luciferase assays in this report (Fig. 2). The allele frequencies of rs1625649 are given in Supplementary  Table 2. In addition to BU induced DNA damage, coadministering agents with BU such as CY and its metabolites also contribute to the DNA damage and thus play a role in the occurrence of aGvHD. MGMT was shown to be involved in repairing damage caused by acrolein, a cyclophosphamide (CY) metabolite [22]. Thus, the association observed in this study might be due to an indirect effect through the interaction of BU with combination chemotherapy. For e.g., depletion of glutathione (GSH) reserves by BU conjugation may reduce the elimination of CY's metabolite acrolein, etoposide's metabolite quinone or melphalan in turn increasing tissue damage [23], which MGMT could repair. Our results demonstrate that higher 1st day BU exposure increases the risk of aGvHD incidence. Though, direct proportional relationship between MGMT and GSH is well known [24], this needs to be further explored the within an HSCT setting in relation to the conditioning regimen. O6-methyl guanine adducts in the absence of MGMT activity could generate point mutations, mismatching base pairs and lead to the formation of DNA double-strand breaks (DSB) [25]. Thus, the inter-play between *GSTA1 genotyping was either performed according to the previously described procedures [9] or with sanger sequencing of the promoter region. GSTA1 metabolic status was based on reporter-gene assays and PK data as described in Ansari et al. [9]. *BU/CY/VP16 was included in this group due to its reported higher toxicity equal to three alkylating agents.     MGMT repair and other DNA-repair pathways in elucidating the cell death and toxicity of alkylating agents used in HSCT conditioning also need to be addressed in future.
Genome wide e-QTL analysis, showed difference in mRNA expression in relation to rs10764881 across diverse human tissues, where allele G demonstrated higher expression profiles [26]. These findings were depicted in LCLs mRNA expression experiments in this report. Our reporter-gene assays indicate that the region further upstream (−700 to −1873) might contain an enhancer element dependent on allele G and to a lesser extent on allele A. MatInspector tool predicted GREs in the human MGMT promoter (Supplementary Fig. 4). Longer plasmid constructs (1873bp) that we used in dual luciferase reporter assays were predicted to have five GRE's, while two GREs were predicted for the shorter plasmid constructs (785 bp). One of these GRE lies within the region encompassing variant rs10764881 and indicates that perhaps this variant could disrupt binding. There are no other SNPs that fall near to any GRE in this region except for rs2782888 (nonpolymorphic in Caucasians but has a MAF of 8% in an African population, Supplementary Fig. 4). Normally, after glucocorticoid receptor-steroid binding occurs, this complex is transported to the nucleus where it can act as a transcription factor enhancing MGMT expression. One could hypothesize that rs10764881 allele A disrupts this binding and results in less efficient MGMT transcription. By treating HEK cells containing the transfected plasmids with dexamethasone we were able to confirm that the plasmid containing SNP rs10764881 A allele demonstrated significantly lower expression levels (p = 0.0007) compared to the rs10764881 G allele construct, in spite of the number of GREs present. Previous studies confirmed the inducible effects of dexamethasone on MGMT mRNA and protein levels through glucocorticoid binding sites [27]. These results were further confirmed by performing an EMSA ( Supplementary Fig. 5). However, in EMSA the protein-DNA band did not disappear completely when the competitive control probe was added suggesting that there could be non-specific binding or other transcription factors that may compete for this site.
One of the limitations of our study is the size of the exploratory cohort that was small and therefore lacked the statistical power to truly reject false negatives. Nevertheless, the significance of MGMT rs10764881 association in an extended cohort indicates true nature of this association. In in vitro MGMT knockout studies neither promoter methylation pattern nor the activity of MGMT in both parent and knockout HAP1cells was measured, Further, multiple clones with and without MGMT expression were not tested that might explain no significant difference in BU IC50 values between cells with and without MGMT. Super shift assays with GRE specific antibodies was also not incorporated in EMSA experiments. The incidences of aGvHD between exploratory and extended cohort were also different, however, the observed association remained significant in the extended cohort. The difference in the observed aGvHD incidences between the exploratory and extended cohorts could be attributed to the differences in the distribution of the stem-cell source, conditioning Association analyses between genotype (additive or dominant models) against aGvHD 2-4 were tested using an X 2 test (Fishers exact test, two sided). Acute GVHD incidences in this cohort were 12%. In the x-axis, the gene and SNP identifications are given and in the y-axis their significance for association with aGvHD 2-4 are shown as p values. B Clinical association analysis between MGMT rs10764881 and aGvHD 2-4 in the extended cohort. Cumulative incidence of acute Graft versus Host Disease (aGvHD 2-4) in the extended sample (N = 182) using competing risk analysis and Cox-regression analysis to calculate the Hazard ratio (HR). Results plotted for MGMT rs10764881 genotype group AA and AG versus GG. The number of patients with aGvHD 2-4 /total number of patients in each group is provided on the plot along with p value and HR for this analysis. The numbers at risk for developing aGvHD 2-4 at each time interval on the x-axis is mentioned below the plot. Variables included in the analysis (backward stepwise conditional coxregression analysis, removing variables with p > 0.2); GSTA1 (rapid and normal metabolizer groups vs. slow metabolizer group); MGMT rs10764881 (GG vs. AA/AG); diagnosis (malignant vs. non-malignant); HLA matching (MRD vs. MUD, MMRD, MMUD); stem-cell source (bone marrow vs. peripheral blood stem cells vs. cord blood); chemotherapy (one alkylating vs. two alkylating agents or three or with VP16); serotherapy (ATG vs. no serotherapy); Day 1 BU AUC (mg × H/L)as a continuous variable harmonized for the dosing schedule (1 × daily or 4 × daily).
regimen used and aGvHD prophylaxis among the patients between the cohorts. These differences could also be attributed to center specific practices. It is known that aGvHD grades 2-4 reduces the risk of relapse in pediatric HSCT, especially in acute lymphoblastic leukemias indicating significant graft versus leukemia effect in children in aGvHD. Chronic GvHD was also shown to reduce the risk of relapse mostly in acute myeloid leukemias (AML), with no survival advantage in both scenarios [28]. We did not observe an association of the MGMT rs10764881 with cumulative incidence of relapse. This may possibly be explained by higher percentage of AML cases in the cohort, varying disease status (heterogeneity) at the time of HSCT and partly explained by the altered expression profiles of the MGMT and its regulation by methylation status of its promoter in cancer cells. In general, cancer tissues exhibit increased MGMT activity compared to that of normal tissues [29], however, concordance between MGMT activity and clinical outcomes in HSCT setting remains to be determined. Several reports in pediatric brain tumors have shown, that there is higher MGMT activity compared to that seen in adults, responded poorly to the alkylating agent temozolamide, and silencing MGMT in such tumors resulted in better response [20] MGMT activity is also determined by hypermethylation of its promoter that results in lower expression of MGMT [30]. We could not investigate the methylation status of the MGMT promoter in this study as varying exposures to various chemotherapy drugs prior to the HSCT in children might have had an impact on the methylation status of the MGMT promoter, that may be more apparent in malignant cells (modulating the clinical outcomes such as relapse) than in the normal cells. The relevance of this genetic variant association in relation to the promoter methylation needs to be investigated in future.
To conclude, children receiving BU-based myeloablative conditioning prior to allo-HSCT and carrying MGMT rs10764881 variant are at increased risk of developing aGvHD 2-4. We hypothesize that children with efficient MGMT function are at lower risk of aaGvHD2-4 possibly by reducing the activation of Phase 1 of the aGvHD cascade. MGMT rs10764881 should be validated in an independent cohort, as a predictive marker of aGvHD, in combination with other associated cytokine polymorphisms and non-genetic factors of the host to perform a pretransplant acute aGvHD risk assessment.

METHODS Patient sample
Exploratory cohort. Sixty children who underwent an allo-HSCT after myeloablative conditioning with BU/CY from 2001 to 2010 at CHU Saint-Justine, Montreal, Canada. This sample was used to genotype fifty-one chosen SNPs within seventeen DNA-repair genes.
Extended cohort. The study was extended to n = 187 by including 122 children who had undergone allo-HSCT at five different centers (Supplementary Table 3). They were genotyped for those SNPs that showed a significant association within the exploratory sample after accounting for multiple testing correction (see Table 1 for the patient's characteristics). The Institutional Review Board at each center approved the study and all patients and/or parents provided informed consent. Details of inclusion criteria are available at Clinicaltrials.gov site (NCT01257854) and the Australian New Zealand Clinical Trials registry (ACTRN12612000544875). Treatment: intravenous. BU (Busulfex ® , Otsuka Pharmaceuticals, Saint-Laurent, Montreal, QC, Canada or Busilvex ® , Pierre Fabre Laboratory, Paris, France) administration was given as a 2 h or 3 h infusion depending on whether the patients were given four times daily or once daily dose, respectively. BU first dose was either age or weight-based and pharmacokinetic (PK) guided dose adjustment was performed in order to obtain a cumulative AUC of 57.6-86.4 mg*h/L (Supplementary Table 1). Co-medication to BU, aGvHD prophylaxis and serotherapy are summarized in Table 1. All patients received non-manipulated grafts.
Clinical outcomes: aGvHD was graded according to established grading criteria [31] and considered up to day 180 post HSCT [32]. Fig. 2 Reporter-gene assay on primary keratinocyte cells treated with and without dexamethasone. Site-specific mutations for MGMT were designed dependent on the SNP of interest in the promoter region. Two fragments of 785 bp (short constructs), both excluding variant rs10764881, but including another SNP, known as rs1625649 which is in high LD with rs10764881 (D′ = 0.99 and R2 = 0.84) were designed. The two other fragments are longer, 1873bp, one that includes rs10764881 allele G and the other with rs10764881 allele A. All four fragments were cloned into pGl4.10 in front of the firefly luciferase gene. Human epidermal primary keratinocytes (HEK cells) were co-transfected with each of the pGL4.10 MGMT constructs and the pRL-SV40 vector that codes for Renilla luciferase for transfection control and normalization. Promoterless pGL4.10 plasmid was used to determine baseline expression. Measurement of Luciferase and Renilla activity was determined by Dual luciferase assays. With the use of the MatInspector tool, a Glucocorticoid responsive element (GRE) was located near to variant rs10764881 and other areas within the plasmid construct region. Thus, to understand whether the enhancing effects are related to corticosteroids, HEK cells transfected with the gene reporter plasmids were stimulated with 0.1 µM dexamethasone (Sigma, D8893) for 15 h and Luciferase expression examined as previously and compared to the non-treated. Difference in promoter activity between the plasmid constructs was assessed by t-test.
Genotyping: Peripheral blood was collected prior to myeloablation and the DNA extracted using a DNA extraction kit (FlexiGene DNA kit, Qiagen GmbH, Hilden). Seventeen genes in total were chosen for investigation ( Univariate analysis: The SNPs were tested for association with aGVHD 2-4 in the exploratory cohort; using a X 2 test in a model (additive or dominant) that best fit the data according to genotype frequency among cases and controls. SNPs with a p value below a 0.001 cut off were retained for further analysis according to the Bonferroni correction. Incidence of aGvHD 2-4 was estimated using cumulative incidence function within competing risk package (cmprsk) using R with death occurring before aGvHD as a competing risk and compared using Grays' test [33] for significant associations. Clinical characteristics tested in the univariate analysis in the extended cohort were HLA compatibility (matched related donor vs. other donors); stem-cell source (bone marrow; peripheral blood stem cells; cord blood); conditioning regimen (based on the number of alkylating agents, one versus 2 or more,); underlying disease (malignant; non-malignant); GSTA1 metabolic capacity based on diplotypes (classified as slow vs. rapid and normal metabolizers);); serotherapy (not received versus received). aGvHD prophylaxis was not assessed separately in the statistical model as it was highly associated with the stem-cell source. Variables with p value <0.1 were retained for multivariable analysis. Clinical characteristics between the exploratory and extended cohorts were compared using X 2 test (categorical) or Mann-Whitney U test (continuous variables).
Multivariable analysis: Significant SNPs were subsequently re-analysed by competing risk analysis to compare the cumulative incidence of aGvHD 2-4. If still significant, they were retained for estimating the Hazard Ratios with a 95% confidence interval (CI) in Cox-regression multivariable analyses with additional risk factors (using a backward stepwise conditional method). Risk factors in the multivariable analysis included: HLA compatibility, stem-cell source; conditioning regimen, use of serotherapy; baseline disease; GSTA1 metabolic capacity based on diplotypes [9]. Two additional PK measures were included as continuous variables: first day BU AUC and cumulative BU AUC. Both have previously been reported in relation to toxicity [5,9]. Additionally, the latter was included in order to control for the variation in the target AUC across the conditioning regimens. Multinomial logistic regression was used to assess the role of the associated SNPs on aGvHD severity, by considering aGvHD grade 1 as reference. For cellular assays, inhibitory concentration 50 (IC50) was determined by nonlinear curve fitting of percent cell survival against concentrations of BU for each cell line in GraphPad Prism, version 7.02.

Investigations examining the functional role of the associated gene and variant(s)
MGMT mRNA expression pre and post Busulfan exposure. 22 CEPH lymphoblastoid cell lines (Coriell Institute, New Jersey, USA) were obtained with known genotypes extracted from the 1000 genome project (https:// www.ncbi.nlm.nih.gov/variation/tools/1000genomes) and used for mRNA expression differences pre and post BU exposure (see Supplementary Material Methods for details).
Cell viability studies on HAP1 MGMT knockouts. To understand the cellular sensitivity of a cell when MGMT is knocked out and after exposure to BU, the near-haploid human cell line HAP1 (Horizon discovery, Cambridge, UK) was used. One Human HAP1 parental control cell line (C631) and one human MGMT knock out, edited by CRISPR/Cas9 to contain a 20 bp deletion in a coding exon of MGMT (HZGHC000430c006). These cell lines were treated (at passage 2-4) with BU concentrations (25, 50, 100, 250, and 500 µM) for 48 h.
Experiments were performed in triplicate on three occasions. Real-Time Cell growth inhibition was evaluated using the CellTiter 2.0 assay (Promega Corporation, 2800 Woods Hollow Road, Madison, USA).

Dual luciferase reporter-gene assays
Functionality of the rs10764881 or for other SNPs that are in LD with this SNP were assessed by site-specific mutations in MGMT promoter pGL4.10 luciferase reporter plasmid constructs (1.8 kb). Luciferase reporter activity was measured by transfecting primary keratinocytes (see Supplementary Material Methods).