- 1Department of Nutrition, Division of Molecular Nutrition, Institute of Basic Medical Sciences, University of Oslo, Oslo, Norway
- 2Center for Treatment of Rheumatic and Musculoskeletal Diseases (REMEDY), Diakonhjemmet Hospital, Oslo, Norway
- 3Institute of Clinical Medicine, University of Oslo, Oslo, Norway
- 4Department of Biostatistics, Oslo Centre for Biostatistics and Epidemiology, Institute of Basic Medical Sciences, University of Oslo, Oslo, Norway
- 5Department of Health Management and Health Economics, Institute of Health and Society, University of Oslo, Oslo, Norway
Methotrexate is one of the cornerstones of rheumatoid arthritis (RA) therapy. Genetic factors or single nucleotide polymorphisms (SNPs) are responsible for 15%–30% of the variation in drug response. Identification of clinically effective SNP biomarkers for predicting methotrexate (MTX) sensitivity has been a challenge. The aim of this study was to explore the association between the disease related outcome of MTX treatment and 23 SNPs in 8 genes of the MTX pathway, as well as one pro-inflammatory related gene in RA patients naïve to MTX. Categorical outcomes such as Disease Activity Score (DAS)-based European Alliance of Associations for Rheumatology (EULAR) non-response at 4 months, The American College of Rheumatology and EULAR (ACR/EULAR) non-remission at 6 months, and failure to sustain MTX monotherapy from 12 to 24 months were assessed, together with continuous outcomes of disease activity, joint pain and fatigue. We found that the SNPs rs1801394 in the MTRR gene, rs408626 in DHFR gene, and rs2259571 in AIF-1 gene were significantly associated with disease activity relevant continuous outcomes. Additionally, SNP rs1801133 in the MTHFR gene was identified to be associated with improved fatigue. Moreover, associations with p values at uncorrected significance level were found in SNPs and different categorical outcomes: 1) rs1476413 in the MTHFR gene and rs3784864 in ABCC1 gene are associated with ACR/EULAR non-remission; 2) rs1801133 in the MTHFR gene is associated with EULAR response; 3) rs246240 in the ABCC1 gene, rs2259571 in the AIF-1 gene, rs2274808 in the SLC19A1 gene and rs1476413 in the MTHFR gene are associated with failure to MTX monotherapy after 12–24 months. The results suggest that SNPs in genes associated with MTX activity may be used to predict MTX relevant-clinical outcomes in patients with RA.
1 Introduction
Rheumatoid arthritis (RA) is a chronic autoimmune disease characterized by systemic inflammation manifested in multiple joints (Smolen and Aletaha, 2006). In the absence of early and effective treatment, the condition can lead to permanent damage, disability, reduced quality of life and life expectancy (Aletaha and Smolen, 2018). The management of RA has undergone significant changes leading up to current treat-to-target strategy, which aims for remission or at least low disease activity (Aletaha and Smolen, 2018; Smolen et al., 2020). Current international guidelines endorse methotrexate (MTX) as a first-line therapy because of its efficacy, cost, and availability of long-term safety data (Gispen et al., 1987; Klareskog et al., 2004; Salliot and van der Heijde, 2009). Despite the effectiveness and clinical utility of MTX, approximately 50%–60% of the patients do not show significant benefit from the treatment due to lack of efficacy or adverse effects (Goekoop-Ruiterman et al., 2007; Hørslev-Petersen et al., 2014; Nam et al., 2014; Emery et al., 2015).
Factors shown to be associated with MTX non-response may be patient-related such as gender and history of smoking and disease-related as duration and severity (Saevarsdottir et al., 2011; Gosselt et al., 2021). Furthermore, levels of anti-citrullinated protein (ACPA) and rheumatoid factor (RF) are considered relevant for MTX responses (Martínez et al., 2020; Maciejewski et al., 2021; Wu et al., 2021). Finally, 15%–30% of variation in drug response are thought to be attributable to patient’s genetic factor or single nucleotide polymorphisms (SNPs) (Eichelbaum et al., 2006).
MTX efficacy has been linked to polymorphism in genes encoding proteins regulating MTX influx, efflux, metabolism, and effector pathways (Figure 1) (Qiu et al., 2017). To this end, cellular uptake of MTX is mediated by reduced folate carrier protein 1 (RFC1, also called SLC19A1) while the ATP-binding cassettes (ABCC) family transporters mediate its cellular efflux (Moscow et al., 1995; Ranganathan and McLeod, 2006; Hider et al., 2007). Once internalized, MTX is reversibly converted to active MTX polyglutamates (MTX-PGs), controlled by the polyglutamation-deconjugation cycle mediated by folypolyglutamate synthetase (FPGS) and gamma-glutamyl hydrolase (GGH) respectively (Chabner et al., 1985). Polyglutamated MTX (MTX-PG) competitively inhibits dihydrofolate reductase (DHFR) and other enzymes such as methylenetetrahydrofolate reductase (MTHFR) contributing to the anti-folate effects of MTX, although not by direct inhibition (Kooloos et al., 2010; Malik and Ranganathan, 2013). Finally, MTX-PG inhibits thymidylate synthase (TYMS) and AICAR (5-aminoimidazole-4-carboxamide ribonucleotide) formyltransferase (ATIC), which contribute to the intracellular accumulation of adenosine (Kooloos et al., 2010; Inoue and Yuasa, 2014; Zhu et al., 2014).
FIGURE 1. Metabolic pathways of intracellular metabolism of MTX. MTX, Methotrexate; SLC, solute carrier; ABC, ATP-binding cassette transporters; FPGS, folylpoly-γ-glutaminase synthetase; GGH, γ-glutamyl hydrolase; MTX-PGs, methotrexate polyglutamate; DHFR, dihydrofolate reductase; DHF, dihydrofolate; THF, tetrahydrofolate; 5,10-MTHF, methyltetrahydrofolate cyclohydrolase; 5-MTHF, L-methylfolate; MTR, methionine synthase; MTRR, methionine synthase reductase; MTHFR, methylenetetrahydrofolate reductase; TYMS, thymidylate synthase; dUMP, deoxyuridine monophosphate; dTMP, deoxythymidine monophosphate; ATIC, 5-aminoimidazole-4-carboxamide ribonucleotide formyltransferase; AICAR, 5-aminoimidazole-4-carboxamide ribonucleotide; AMP, adenosine monophosphate; FAICAR, 5-formylaminoimidazole-4-carboxamide ribonucleotide; IMP, inosine monophosphate; AMP, adenosine monophosphate; AIF-1, allograft inflammatory factor 1 and GGH, gamma-glutamyl hydrolase. SNPs included in the current study are marked with blue font and star (*).
The association between polymorphisms in these genes and MTX response have been described (Owen et al., 2013a; Lv et al., 2018; D’Cruz et al., 2020). However, due to heterogeneous study populations and replication failure, different conclusions have been reported. Thus, the identification and development of clinically effective SNP biomarkers for predicting MTX sensitivity remains a challenge. It is therefore a need to replicate and validate candidate SNPs in various and large patient cohorts. In the present study, we aimed to investigate the associations between 23 SNPs and broad clinically relevant outcomes to assess their contributions to MTX efficacy in a cohort of 224 Norwegian patients with RA who were naïve to MTX and included in the ARCTIC trial (Haavardsholm et al., 2016). The outcomes assessed were categorical outcomes, including ACR/EULAR Boolean non-remission, DAS-based EULAR non-response and failure to MTX monotherapy as well as continuous outcomes including disease activity, fatigue and joint pain. We found that variations in genes encoding key components regulating MTX influx, efflux, metabolism pathways may be used to evaluate clinical outcomes in patients with RA.
2 Materials and methods
2.1 Subjects
The current ARCTIC_SNPs study included 224 patients, a sub-cohort from the original ARCTIC clinical trial (Haavardsholm et al., 2016). The ARCTIC trial was a 24-month randomized, open label, prospective, national, multi-centric study (11 centers in Norway) conducted in accordance with the Declaration of Helsinki and the International Conference on Harmonization Guidelines for Good Clinical Practice. Informed consent was taken from all the study subjects. All patients must have satisfied the following criteria: 1. Age between 18 and 75 years; 2. Fulfillment of the 2010 American College of Rheumatology (ACR)/European League against Rheumatism (EULAR) criteria for RA; 3. No prior treatment with DMARD including MTX; 4. Disease duration less than 2 years and 5. Indication for DMARD treatment. Patients with major comorbidities and abnormal liver or kidney function were excluded from the study (Haavardsholm et al., 2016). After diagnosis, all patients received MTX together with tapering dose of prednisolone.
2.2 Outcome assessment
Sensitivity to MTX treatment was evaluated based on both categorical and continuous outcomes. The primary categorical outcome for the study was that the patient did not achieve ACR/EULAR Boolean remission at 6 months of follow-up (‘non-remission’). The patient must satisfy the following criteria in order to achieve ACR/EULAR boolean remission: Ritchie Articular Index ≤1 (RAI≤1); Swollen Joint Count ≤1 (SJC44 ≤ 1); C-Reactive Protein ≤1 (CRP≤1); and the Patient’s Global Assessment ≤1 (PGA≤1) (Felson et al., 2011). The secondary categorical outcome was according to disease activity scores (DAS) based EULAR non-response at 4 months of follow-up. This classifies patients into good, moderate- and non-responders based both on level and changes in DAS (Supplementary Table S1) (Fransen and van Riel, 2005). Good responders were denoted as favorable outcome and moderate/non-responders were considered as unfavorable outcome. The other secondary categorical outcome was the failure to sustained MTX monotherapy defined based on persistence of MTX monotherapy. Patients with MTX monotherapy at both 12 months and 24 months reflect favorable outcome otherwise unfavorable. In addition, continuous outcomes include DAS, clinical disease activity index (CDAI) (Aletaha et al., 2005), simplified disease activity index (SDAI) (Smolen et al., 2003), joint pain and fatigue at 4 or 6 months, as well as improvement of these measures from baseline.
2.3 Selection of single-nucleotide polymorphisms
Twenty-five SNPs in nine candidate genes (marked in blue in Figure 1) were selected based on literature attending to their putative contribution in MTX metabolic pathways and mechanism of transport and action. They included, genes involved in MTX cellular influx (SLC19A1 and SLC22A11; -efflux (ABCC1), –and polyglutamate formation (GGH), in addition to the folate- (DHFR and MTHFR) and methionine pathways (MTR and MTRR), as well as de novo purine synthesis (ATIC) and finally, inflammatory cytokines (AIF-1) (Qiu et al., 2017) (Supplementary Table S2).
2.4 DNA isolation and genotyping
Genomic DNA was extracted from full blood by using GeneJET Whole Blood Genomic DNA Purification Mini Kit as per manufacturer’s instruction (ThermoFisher, Waltham, United States). Total genomic DNA was quantified and then analyzed for purity using the NanoDrop 1000 Spectrophotometer v3.7 (Thermo Fisher Scientific, Waltham, United States). Sequenom® Assay Designer 3.1 software was used for primers design, and samples were genotyped using the Sequenom iPLEX® MassARRAY platform according to the manufacturer’s instructions (OE Biotechnology Co., Ltd., Shanghai, China).
2.5 Statistical analysis
As quality control, SNPs were excluded from analysis when genotyping success rates were less than 95% or when minor allele frequency (MAF) was less than 5%. Moreover, deviation of SNPs from Hardy-Weinberg equilibrium (HWE) was assessed by a Pearson’s chi-squared test with a significance threshold of p < 0.05 and those SNPs were excluded. Demographic and clinical categorical variables were summarized by frequency and percentage while continuous variables were summarized by median and interquartile range (IQR, 25–75th percentile). Differences in the categorical variables between subgroups were analyzed by Pearson’s chi-squared test and Mann–Whitney U test for continuous variables. To analyze the associations between SNPs and categorical clinical outcomes, both univariate and multivariate analyses were performed using three genetic models. These included the genotype model denoted as major homozygote vs. heterozygote vs. minor homozygote, the recessive model denoted as major homozygote + heterozygote vs. minor homozygote and the dominant model which was major homozygote vs. heterozygote + minor homozygote). To implement univariate analysis, the odds ratio (OR) (95% confidence interval) of each SNP was calculated from cross-tabulation and Pearson’s chi-squared test, or Fisher exact test when expected case number ≤5, was performed. Haldane-Anscombe correction was adopted if the count number in a subgroup was zero when calculating OR. To conduct multivariate analysis, a logistic regression model was fitted with covariates age, gender, BMI and smoking in addition to the SNP variable, and the corresponding adjusted OR for each genotype was calculated from the estimated regression coefficient of the SNP variable in the model. To identify the association between SNPs and continuous clinical outcomes, Spearman’s rank correlation was used to calculate correlation coefficient ρ and p-values for each SNP encoded by counts of minor allele (major homozygote: 0 as reference, heterozygote: 1, minor homozygote: 2). Univariate and multivariate analysis by ordinary least squares linear regression (OLS) was employed to further explore the association between SNP and continuous endpoints, with age, gender, BMI and smoking as covariates in multivariate models. Results were considered statistically significant when two-sided p-values were less than 0.05. The adjustment of significant results for multiple testing by Bonferroni correction (significant level of 0.0022 corrected based on a group of hypotheses tested on 23 SNPs) was reported in all analyses on association between SNPs and outcomes. OLS is the only analysis where significance remains after Bonferroni correction. For other analyses we also chose to report the raw p-values (without multiple testing adjustment) because of the exploratory nature of the analyses and to better highlight nuances in the results that would be less visible after correction.
3 Results
3.1 Patient inclusion
Out of 238 patients from the ARCTIC trial (Haavardsholm et al., 2016), 14 patients were excluded as blood-samples were missing, leaving 224 patients for the present study. Next, 16 plus 4 plus 4 patients were excluded due to missing clinical data, inadequate information on treatments or transition from MTX monotherapy at 4 months, and genotyping failure respectively. This resulted in inclusion of 200 patients in the present analysis (Figure 2).
FIGURE 2. Flow chart for patient’s inclusion and exclusion criteria. N, number of patients; MTX, methotrexate; SNPs, Single nucleotide polymorphisms.
3.2 Baseline characteristics
The demographic and clinical characteristics of the recruited are summarized in Table 1 and show that 78 (39%) were males and 122 (61%) were females with median age of 54 years. Seventy (35%) of the patients were never smokers, 86 (43%) were former smokers while 44 (22%) were active smokers at the time of RA diagnosis. The median body mass index (BMI) was 25.1. Moreover, the median of disease duration indicating months since first symptoms was 5.5 months while 144 (72%) and 168 (84%) of the patients were RF and anti-cyclic citrullinated protein (anti-CCP) positive, respectively. The median baseline value of CDAI and SDAI before treatment were 21 and 22, respectively. The median PGA and physician global assessment score (PhGA) at the baseline were 48 and 36 respectively, whereas erythrocyte sedimentation rate (ESR) and CRP for patients were 20 mm/h and 7 mg/L. Moreover, median value of joint pain, fatigue, patients reported outcomes measurement information physical function T-score (PROMIS-PF) and RA impact disease score (RAID) before treatment were 45.5, 37.5, 39.6, and 4.4 respectively. The median baseline value of DAS before treatment initiation was 3.32. Median start dose for treatment was 15 mg MTX/week at baseline.
3.3 Outcomes
At 6 months, 71 (35.5%), and 129 (64.5%) and patients were classified as remission and non-remission to MTX treatment, respectively, according to ACR/EULAR Boolean remission criteria. There were no statistically significant differences in categorical variables (gender/smoking/RF-status/anti-CCP status) and disease duration variable between the remission and non-remission groups. By contrast, the non-remission group had significantly higher BMI, disease activity composite scores, tender and painful joints assessed by RAI, PGA and PhGA values at the baseline. Patients not achieving remission at 6 months also showed significantly elevated levels of ESR, CRP, joint pain, fatigue and RAID at the baseline, as compared to remission group (see Table 1 for details). Finally, the median baseline value of PROMIS-PF for the non-remission group was significantly lower (indicating worse physical condition) as compared to remission group. For another two outcomes, the comparison of variables in two subgroups is shown in Supplementary Figure S1.
3.4 Single nucleotide polymorphisms selection
As depicted in Table 2 and Figure 1, a total of 25 SNPs were selected, based on their involvement in MTX uptake and metabolism pathways. In Table 2, the gene names and chromosome location of each SNP together with genotypic distribution of minor allele frequency (MAF) and Hardy-Weinberg equilibrium (HWE) are listed. Two SNPs, rs12681874 located in GGH and rs3821353 in ATIC were excluded due to failure in genotyping success rate (<95%) and non-agreement with HWE (p < 0.05), respectively, leaving 23 SNPs for the final analysis.
3.5 Single nucleotide polymorphisms associated with categorical outcomes post methotrexate treatment
Three models, a genotype-, recessive- and dominant genetic model, were used in the present study to explore the association between SNP and classification outcomes. Pearson’s Chi-square test was used for univariate analysis of each SNP, and multivariable logistic regression model was applied with adjustment for gender, BMI, age, and smoking to assess potential confounding effects. We found that minor homozygote of two SNPs, TT in rs1476413 (C>T) in the MTHFR and AA in rs7884864 (G>A) in the ABCC1 genes were associated with ACR/EULAR non-remission at 6 months in both the recessive and the genotype minor model (Figure 3A red dots, and Table 3A). No SNPs were significantly associated with EULAR non-response at 4 months (Figure 3B; Table 3B). By contrast, the minor allele T and heterozygote CT of SNP rs1801133 in gene MTHFR was significantly associated with EULAR response at 4 months in the dominant and genotype-heterozygote models, respectively (Figure 3B blue dot and Table 3B). Furthermore, the minor homozygote GG of rs246240, the heterozygote CT of rs2274808 and the minor homozygote TT of rs1476413 were all significantly associated with failure of sustained MTX monotherapy in the recessive, genotype-heterozygote and genotype-minor models, respectively (Figure 3C red dots and Table 3C). In addition to this, the minor allele G and heterozygote TG of rs2259571 was significantly associated with sustained MTX monotherapy in dominant and genotype-heterozygote model respectively (Figure 3C blue dots and Table 3C). SNPs found not associated with clinical outcomes are shown in Supplementary Table S3–5. However, it is worth mentioning that all these associations became insignificant when significance level was adjusted by Bonferroni correction.
FIGURE 3. SNPs associated with MTX treatment inefficiency measured by categorical outcomes. Volcano figures showing SNPs significantly associated with MTX treatment inefficiency. Multivariate logistic regression (adjusted by gender, age, smoking, BMI) was implemented to calculate odds ratio and adjusted p-value (without Bonferroni correction) of each SNP. For those where case number in subgroups is zero, chi-square test was used and odds ratio was estimated by contingency table. MTX treatment inefficiency was measured by three criteria and denoted as ACR/EULAR non-remission at 6 months (A), EULAR non-response at 4 months (B) and failure to MTX monotherapy from 12 to 24 months (C). Recessive, dominant, genotype-heter and genotype-minor models were applied respectively. They are denoted as mm vs. MM + Mm (reference), mm + Mm vs. MM (reference), Mm vs. MM (reference) and mm vs. MM (reference) accordingly. The Bonferroni corrected p-value was estimated based on hypotheses tested on 23 SNPs to adjust significant results for multiple testing. The raw p-values (without Bonferroni correction adjustment) were also reported because of the exploratory nature of the analyses and to better highlight nuances in the results that would be less visible after correction. MM, major homozygote; Mm, heterozygote; mm, minor homozygote.
3.6 Identification of single nucleotide polymorphisms associated with continuous clinical outcomes post methotrexate monotherapy
Disease activity is monitored by DAS, SDAI, CDAI, and patients report outcomes (PRO) capturing features such as joint pain, fatigue. Each SNP was denoted by numbers based on count of minor allele and Spearman’s rank correlation coefficient was estimated for each SNP with these continuous clinical outcomes. Generally, the coefficient for each SNP is low ranging from −0.2 to 0.2 and we further chose the associated one with coefficient >0.18 or < −0.18 for displaying. The heatmap in Figure 4 depicts a negative association between the increasing count of the minor allele G in rs1805087 in MTR gene and joint pain at 6 months, which is in line with the trend shown in scatter plot in Figure 4A. By contrast, positive associations were found between the increasing count of the minor allele A of SNP rs1801394 in the MTRR gene and four outcomes including improved CDAI and SDAI both at 4 and 6 months. This was in consistent with trends as shown in Figures 4B–E. Additionally, it should be note that the three SNPs, rs2274808 and rs9977268 in the SLC19A1, and rs1801133 in the MTHFR gene, were all significantly associated with improved fatigue at 4 months (Figure 4 heatmap and F-H). Whereas increasing number of T allele in rs2274808 and rs9977268 were positively correlated with improved fatigue at 4 months, higher number of T allele in rs1801133 were negatively correlated (Figure 4 heatmap and F-H). As the same case in the analyses for categorical outcomes, the correlations between SNPs and outcomes became insignificant with Bonferroni correction.
FIGURE 4. Identification of SNPs correlated to different continuous clinical outcomes. Heatmap showing the correlation between SNPs involved in MTX metabolism and continuous clinical outcomes by Spearman’s rank correlation test. The correlation with significant p-value (calculated by Spearman rank **raw p-value ≤ 0.01) and relatively high coefficient (absolute value of correlation coefficient >0.18) were marked by asterisk. These significant associations were further shown by scatter plots (A–H), including rs1805087 and joint pain_6 (A), rs1801394 and ΔCDAI_4 (B), rs1801394 and ΔCDAI_6 (C), rs1801394 and ΔSDAI_4 (D), rs1801394 and ΔSDAI_6 (E), rs2274808 and ΔFatigue_4 (F), rs9977268 and ΔFatigue_4 (G), rs1801133 and ΔFatigue_4 (H). The Bonferroni corrected p-value was estimated based on hypotheses tested on 23 SNPs to adjust significant results for multiple testing. The raw p-values (without Bonferroni correction adjustment) were also reported because of the exploratory nature of the analyses and to better highlight nuances in the results that would be less visible after correction. Delta values means improvements in clinical outcomes from baseline value (baseline value–current value). DAS, Disease activity score; CDAI, Clinical Disease Activity Index; SDAI, Simplified Disease Activity Index.
Next, we further investigated the association between SNP and clinical continuous outcomes by using linear regression with ordinary least squares (OLS) in both univariate model and multivariate model (Figures 5A,B). Based on Bonferroni corrected significance level, both models identified significant positive association between rs1801394_heter (GA) and improved CDAI at 4 and 6 months, which was in line with the result in Spearman’s correlation analysis. Moreover, the negative association between rs1801133_minor (TT) and improved fatigue at 4 and 6 months, as well as the negative association between rs408626_heter (AG) and CDAI/SDAI at 6 months were identified by both univariate and multivariate analysis. Of note, after adjustment by covariates, rs2259571_heter (TG) become significantly associated with disease activities, including DAS, CDAI and SDAI at 4 months. (Figure 5B).
FIGURE 5. Identifying associated SNPs with continuous clinical outcomes by linear regression. Heatmap showing the association between SNPs involved in MTX metabolism and continuous clinical outcomes by linear regression with ordinary least squares regression (OLS). Coefficient and p-value shown in heatmap (A,B) were respectively estimated by univariate and multivariate regression analysis with age, BMI, gender and smoking as covariates. Significance level was corrected by Bonferroni correction for multiple testing and significant associations were marked by asterisk (*p-value ≤ 0.05, ***p-value ≤ 0.001). Delta values means improvements in clinical outcomes from baseline value (baseline value–current value). DAS, Disease activity score; CDAI, Clinical Disease Activity Index; SDAI, Simplified Disease Activity Index.
4 Discussion
Despite the long-term use of MTX for patients with RA, our ability to predict clinical therapeutic response to MTX is limited. Recently, studies of association between MTX metabolism related- SNPs and therapeutic response demonstrated their involvement in interpatient variability in clinical response. Here a comprehensive analysis of SNPs in genes encoding proteins for MTX uptake and metabolism and their association with multiple clinical outcomes post MTX treatment was studied in patients with RA. The disclosed SNPs include rs1801394, rs408626 and rs2259571 associated with disease activity relevant continuous outcomes, as well as rs1801133 associated with improved fatigue. Additionally, potential associations between two SNPs, rs1476413, rs3784864 and ACR/EULAR non-remission; between one SNP rs1801133 and EULAR response; between four SNPs rs246240, rs2259571, rs2274808, rs1476413 and failure to MTX monotherapy after 12–24 months were identified, despite insignificant when Bonferroni correction applied.
The increased number of minor allele A in rs1801394 G>A in the MTRR was demonstrated to be associated with more improved CDAI at 4 and 6 months by Spearman’s rank correlation test. This was further supported by both univariate and multivariate linear regression analysis showing rs1801394_GA significantly associated with more improved CDAI with Bonferroni corrected p-value as compared with GG. The similar trend can be also observed in rs1801394_AA from Figures 4B,C, despite no statistical significance found due to limited number of participants. This concluded a favorable role of the allele A in this SNP associated with positive continuous outcomes after MTX treatment. Our results were in contrast to one study on 915 European Caucasian origin patients with RA where A allele was associated with less improvement than the G allele according to ΔDAS28 (p = 0.0016) and EULAR response (p = 0.004) (López-Rodríguez et al., 2018). Additionally, other studies on Chinese and Indian population identified no association between the SNP and MTX response (Ghodke-Puranik et al., 2015; Lv et al., 2018). Due to difference in the analyzed population, a clear conclusion would be precluded. More studies on Nordic population remain needed for validation.
Another disease activity associated SNP is rs2259571 TG in AIF-1 gene which was found to be significantly associated with low level of DAS, CDAI and SDAI at 4 months. Moreover, The T allele and TG genotype were both shown to be negatively associated with failure to MTX monotherapy, despite at an uncorrected significance level. AIF-1 is a cytoplasmic, inflammation-responsive protein and Pawlik and coworkers demonstrated that the patients with the AIF-1 rs2259571 GG genotype have a poorer response (DAS28 > 2.4 at 6 months) to MTX treatment by studying 221 Poland patients with RA (p = 0.03, OR = 2.44) (Pawlik et al., 2013). This was in contrast to the present results which demonstrated that G allele was associated with more sustained MTX monotherapy and low level of disease activity. Our result was also supported by machine learning analysis in the other paper by us identifying rs2259571 TG as predictor for more sustained MTX monotherapy. The inconsistency in conclusion from other study may result from the different definition of response outcome and distinct population for analysis. Our result stressed important role of rs2259571 in predicting MTX related outcome and more research is warranted for further investigation of this SNP. In addition to rs1801394 and rs2259571, rs408626 AG genotype in DHFR gene was identified to be associated with low value of CDAI and SDAI at 6 months, this is in contrast to other study where AA genotype was associated with the less reduction in relative DAS28 score (p = 0.05) (Milic et al., 2012).
The TT genotype of rs1801133 in the MTHFR gene which was negatively associated with improvement in fatigue at 4 months post MTX monotherapy. As this was shown both through Spearman’s rank correlation test and the univariate and multivariate linear regression with Bonferroni corrected p-value, suggesting for a potent role in predicting fatigue by screening for this allele variant. Moreover, this SNP was further associated with EULAR response outcome in the dominant and genotype-heterozygote model at an uncorrected significance level. This may implicate that population carrying T allele (CT and TT) or CT may have higher proportion of responders as compared with those carrying the CC genotype. Recent studies on this SNP showed inconsistent result. Whereas a study from Poland (n = 174) demonstrated a similar result as the present study that T allele was shown to be associated with improved remission outcomes (p = 0.02, OR = 0.46) (Kurzawski et al., 2007), a different study on 233 Portuguese patients found that TT were associated with increased risk for non-response to MTX defined by DAS28< =3.2 (p = 0.013, OR = 4.63) (Lima et al., 2014). Moreover, studies on Asian population including 217 Indian RA population, 110 Chinese and 67 Pakistani showed no association between the SNP and treatment outcomes (Xiao et al., 2010; Ghodke-Puranik et al., 2015; Iqbal et al., 2015). Due to conflicted results from different studies on different ethnical population, a clear conclusion was precluded and more studies on Norwegian population are needed.
The MTHFR gene was further associated with both non-remission post MTX treatment and failure to sustained MTX monotherapy through TT genotype of SNP rs1476413 C>T, at an uncorrected significance level. The MTHFR gene encodes an enzyme responsible for thymidine and methionine synthesis, which are one of targets of MTX (Figure 1) (Owen et al., 2013b). The association of rs1476413 TT with failure in responding to MTX treatment is consistent with a study by Salazar and coworkers who identified that T allele was associated with a poor response to MTX treatment in 124 Spanish patients with RA (p = 0.0086, OR = 3.56) (Salazar et al., 2014), despite difference in the outcome definition.
The SNP rs3784864 AA in ABCC1 gene in the ABC encoding gene was associated with ACR/EULAR non-remission at 6 months at an uncorrected significance. The ABC gene family are MTX transporters responsible for MTX efflux. The gene ABCC1 is located on chromosome 16p13 and the SNP rs3784864 is characterized by an intronic G to an A substitution (Haga and Burke, 2004). Our study found that population with AA was associated with non-remission. This observation is in contrast to a recent report demonstrating that MTX non-response was associated with G carriers of rs3784864 in 233 Portugal patients with RA (p = 0.015, OR = 4.24) (Lima et al., 2015). Different from Boolean ACR/EULAR Remission as the outcome used in this study, they defined non-response as DAS28 > 3.2 in two consecutive evaluations with 3 months as each interval evaluation. These different ethical backgrounds of population and outcome definitions in these studies may be the reason for this inconsistency in results.
Moreover, the SNP rs246240 GG in ABCC1 and rs2274808 CT in SLC19A1 was positively associated with failure to sustain MTX monotherapy at an uncorrected significance level, suggesting that populations with genotype rs246240 GG and rs2274808 CT might have a higher risk of MTX withdrawal before 12 months. Such a role for rs246240 has been validated by Lima and coworkers in a multivariate genotype analysis which showed that G carriers are associated with poor response to MTX in 233 Portugal population (p = 0.008, OR = 5.47) (Lima et al., 2015). By contrast, the role for rs2274808 is inconsistent with other study: A study in 248 United Kingdom population found the T allele associated with failure to MTX treatment (p = 0.009, OR = 1.76) (Owen et al., 2013a), whereas in the present study rs2274808 CT and not TT was found to be associated with early discontinuation of MTX monotherapy compared to RA patients carrying the CC genotype.
The main strengths our study include 1) high quality of data: all patients were naïve to DMARDs and initiated with MTX as first-line therapy, with the studied population being highly homogenous regarding ethnic origin and thus representative; 2) a comprehensive assessment of 23 SNPs associated with MTX metabolism and inflammatory response as well as a variety of relevant outcome measures 3) a rigorous statistical evaluation taking into consideration common potential cofounders including age, BMI, smoking and gender, with Bonferroni correction for multiple testing to avoid false positive findings. Nevertheless, several potential limitations should be noted. First, the sample size of our present study was relatively small (n = 200), thus statistical power was relatively limited to detect variants with weak effect sizes. Secondly, there is a possibility of false negative rate caused by Bonferroni correction, which may filter some truly significantly associated SNPs. To this end, we also reported uncorrected raw p-value for analyses on association for categorical outcomes and on Spearman rank correlation for continuous outcomes.
5 Conclusion
In summary, we have identified and validated several SNPs associated with a variety of clinical outcomes post MTX monotherapy in a Norwegian population with early DMARD naïve RA. Among these SNPs, rs1801394, rs408626, and rs2259571 are significantly associated with disease activity relevant continuous outcomes, as well as rs1801133 associated with improved fatigue. Moreover, rs1476413 and rs246240 had been found associated with categorical outcomes of ACR/EULAR non-remission and failure to MTX monotherapy respectively, which was identical to previous studies and thus strengthened reliability of these findings. Despite this, more investigation should be conducted to validate the results in other cohorts before these may be applied in clinical decision-making.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.
Ethics statement
The studies involving human participants were reviewed and approved by Regional Committees for Medical and Health Research Ethics (REK) (reference number: 2010/744/REK sør-øst C). The patients/participants provided their written informed consent to participate in this study.
Author contributions
SK designed the study, conducted DNA extraction and contributed to writing the manuscript. GL designed the study, analyzed the data and contributed to writing the manuscript. FG designed the study, conducted the DNA extraction and revised the manuscript. JS assisted with data transfer and revised the manuscript. GG, TK, and NS were involved in ARCTIC project and revised the manuscript. MZ supervised the data analysis and revised the manuscript. SL was involved in ARCTIC project, designed the study and revised the manuscript. BS and EH led the project, designed and supervised the study, and contributed to writing the manuscript. All authors contributed to the article and approved for the submitted version.
Funding
This study was funded by the UiO-MED and UiO-IMB, Grants to GL and BSS, The Norwegian Research Council, grant #ES612356, and Throne Holst Foundation, grant # BSS 2019/2022. Principal funding recipient is BSS. The ARCTIC trial was supported by the Norwegian Research Council, Norwegian South-Eastern Health Region, Norwegian Women’s Public Health Association, Norwegian Rheumatism Association, investigator initiated research grants from AbbVie, UCB Pharma, Pfizer Inc, MSD Norway, and Roche Norway.
Acknowledgments
The authors thank all the patients for participating in the ARCTIC study, as well as investigators and medical staff for their contributions to patient recruitment, clinical measurement and data collection in the ARCTIC study.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphar.2022.1075603/full#supplementary-material
References
Aletaha, D., Nell, V. P. K., Stamm, T., Uffmann, M., Pflugbeil, S., Machold, K., et al. (2005). Acute phase reactants add little to composite disease activity indices for rheumatoid arthritis: Validation of a clinical activity score. Arthritis Res. Ther. 7 (4), R796–R806. doi:10.1186/ar1740
Aletaha, D., and Smolen, J. S. (2018). Diagnosis and management of rheumatoid arthritis: A review. JAMA 320 (13), 1360–1372. doi:10.1001/jama.2018.13103
Chabner, B. A., Allegra, C. J., Curt, G. A., Clendeninn, N. J., Baram, J., Koizumi, S., et al. (1985). Polyglutamation of methotrexate. Is methotrexate a prodrug? J. Clin. . 76 (3), 907–912. doi:10.1172/JCI112088
D’Cruz, L. G., McEleney, K. G., Tan, K. B. C., Shukla, P., Gardiner, P. V., Connolly, P., et al. (2020). Clinical and laboratory associations with methotrexate metabolism gene polymorphisms in rheumatoid arthritis. J. Pers. Med. 10 (4), E149. doi:10.3390/jpm10040149
Eichelbaum, M., Ingelman-Sundberg, M., and Evans, W. E. (2006)., 57. United States, 119–137. doi:10.1146/annurev.med.56.082103.104724Pharmacogenomics and individualized drug therapyAnnu. Rev. Med.
Emery, P., Burmester, G. R., Bykerk, V. P., Combe, B. G., Furst, D. E., Barre, E., et al. (2015). Evaluating drug-free remission with abatacept in early rheumatoid arthritis: Results from the phase 3b, multicentre, randomised, active-controlled AVERT study of 24 months, with a 12-month, double-blind treatment period. Ann. Rheum. Dis. 74 (1), 19–26. doi:10.1136/annrheumdis-2014-206106
Felson, D. T., Smolen, J. S., Wells, G., Zhang, B., van Tuyl, L. H. D., Funovits, J., et al. (2011). American College of Rheumatology/European League against Rheumatism provisional definition of remission in rheumatoid arthritis for clinical trials. Arthritis Rheum. 63 (3), 573–586. doi:10.1002/art.30129
Fransen, J., and van Riel, P. L. C. M. (2005). The disease activity score and the EULAR response criteria. Clin. Exp. Rheumatol. 23 (5), S93–S99.
Ghodke-Puranik, Y., Puranik, A. S., Shintre, P., Joshi, K., Patwardhan, B., Lamba, J., et al. (2015). Folate metabolic pathway single nucleotide polymorphisms: A predictive pharmacogenetic marker of methotrexate response in Indian (asian) patients with rheumatoid arthritis. Pharmacogenomics 16 (18), 2019–2034. doi:10.2217/pgs.15.145
Gispen, J. G., Alarcon, G. S., Johnson, J. J., Acton, R. T., Barger, B. O., and Koopman, W. J. (1987). Toxicity of methotrexate in rheumatoid arthritis. J. Rheumatol. 14 (1), 74–79.
Goekoop-Ruiterman, Y. P., de Vries-Bouwstra, J. K., Allaart, C. F., van Zeben, D., Kerstens, P. J., Hazes, J. M., et al. (2007). Comparison of treatment strategies in early rheumatoid arthritis: A randomized trial. Ann. Intern Med. 146 (6), 406–415. doi:10.7326/0003-4819-146-6-200703200-00005
Gosselt, H. R., Verhoeven, M. M. A., Bulatovic-Calasan, M., Welsing, P. M., de Rotte, M. C. F. J., Hazes, J. M. W., et al. (2021). Complex machine-learning algorithms and multivariable logistic regression on par in the prediction of insufficient clinical response to methotrexate in rheumatoid arthritis. J. Pers. Med. 11 (1), 44. doi:10.3390/jpm11010044
Haavardsholm, E. A., Aga, A. B., Olsen, I. C., Lillegraven, S., Hammer, H. B., Uhlig, T., et al. (2016). Ultrasound in management of rheumatoid arthritis: ARCTIC randomised controlled strategy trial. BMJ Clin. Res. ed. 354, i4205. doi:10.1136/bmj.i4205
Haga, S. B., and Burke, W. (2004). Using pharmacogenetics to improve drug safety and efficacy. JAMA 291 (23), 2869–2871. doi:10.1001/jama.291.23.2869
Hider, S. L., Bruce, I. N., and Thomson, W. (2007). The pharmacogenetics of methotrexate. Rheumatol. Oxf. 46 (10), 1520–1524. doi:10.1093/rheumatology/kem147
Hørslev-Petersen, K., Hetland, M. L., Junker, P., Podenphant, J., Ellingsen, T., Ahlquist, P., et al. (2014). Adalimumab added to a treat-to-target strategy with methotrexate and intra-articular triamcinolone in early rheumatoid arthritis increased remission rates, function and quality of life. The OPERA study: An investigator-initiated, randomised, double-blind, parallel-group, placebo-controlled trial. Ann. Rheum. Dis. 73 (4), 654–661. doi:10.1136/annrheumdis-2012-202735
Inoue, K., and Yuasa, H. (2014). Molecular basis for pharmacokinetics and pharmacodynamics of methotrexate in rheumatoid arthritis therapy. Drug Metab. Pharmacokinet. 29 (1), 12–19. doi:10.2133/dmpk.dmpk-13-rv-119
Iqbal, M. P., Ali, A. A., Mehboobali, N., and Iqbal, K. (2015). Short Communication: Lack of association between MTHFR gene polymorphisms and response to methotrexate treatment in Pakistani patients with rheumatoid arthritis. Pak J. Pharm. Sci. 28 (5), 1789–1792.
Klareskog, L., van der Heijde, D., de Jager, J. P., Gough, A., Kalden, J., Malaise, M., et al. (2004). Therapeutic effect of the combination of etanercept and methotrexate compared with each treatment alone in patients with rheumatoid arthritis: Double-blind randomised controlled trial. Lancet 363 (9410), 675–681. doi:10.1016/S0140-6736(04)15640-7
Kooloos, W. M., Huizinga, T. W. J., Guchelaar, H. J., and Wessels, J. A. M. (2010). Pharmacogenetics in treatment of rheumatoid arthritis. Curr. Pharm. Des. 16 (2), 164–175. doi:10.2174/138161210790112764
Kurzawski, M., Pawlik, A., Safranow, K., Herczynska, M., and Drozdzik, M. (2007). 677C>T and 1298A>C MTHFR polymorphisms affect methotrexate treatment outcome in rheumatoid arthritis. Pharmacogenomics 8 (11), 1551–1559. doi:10.2217/14622416.8.11.1551
Lima, A., Bernardes, M., Azevedo, R., Medeiros, R., and Seabra, V. (2015). Pharmacogenomics of methotrexate membrane transport pathway: Can clinical response to methotrexate in rheumatoid arthritis Be predicted? Int. J. Mol. Sci. 16 (6), 13760–13780. doi:10.3390/ijms160613760
Lima, A., Monteiro, J., Bernardes, M., Sousa, H., Azevedo, R., Seabra, V., et al. (2014). Prediction of methotrexate clinical response in Portuguese rheumatoid arthritis patients: Implication of MTHFR rs1801133 and ATIC rs4673993 polymorphisms. BioMed Res. Int. 2014, 368681. doi:10.1155/2014/368681
López-Rodríguez, R., Ferreiro-Iglesias, A., Lima, A., Bernardes, M., Pawlik, A., Paradowska-Gorycka, A., et al. (2018). Replication study of polymorphisms associated with response to methotrexate in patients with rheumatoid arthritis. Sci. Rep. 8 (1), 7342–7348. doi:10.1038/s41598-018-25634-y
Lv, S., Fan, H., Li, J., Yang, H., Huang, J., Shu, X., et al. (2018). Genetic polymorphisms of TYMS, MTHFR, ATIC, MTR, and MTRR are related to the outcome of methotrexate therapy for rheumatoid arthritis in a Chinese population. Front. Pharmacol. 9, 1390. doi:10.3389/fphar.2018.01390
Maciejewski, M., Sands, C., Nair, N., Ling, S., Verstappen, S., Hyrich, K., et al. (2021). Prediction of response of methotrexate in patients with rheumatoid arthritis using serum lipidomics. Sci. Rep. 11 (1), 7266. doi:10.1038/s41598-021-86729-7
Malik, F., and Ranganathan, P. (2013). Methotrexate pharmacogenetics in rheumatoid arthritis: A status report. Pharmacogenomics 14 (3), 305–314. doi:10.2217/pgs.12.214
Martínez, G., Feist, E., Martiatu, M., Garay, H., and Torres, B. (2020). Autoantibodies against a novel citrullinated fibrinogen peptide related to smoking status, disease activity and therapeutic response to methotrexate in cuban patients with early rheumatoid arthritis. Rheumatol. Int. 40 (11), 1873–1881. doi:10.1007/s00296-020-04580-x
Milic, V., Jekic, B., Lukovic, L., Bunjevacki, V., Milasin, J., Novakovic, I., et al. (2012). Association of dihydrofolate reductase (DHFR) -317AA genotype with poor response to methotrexate in patients with rheumatoid arthritis. Clin. Exp. Rheumatol. 30 (2), 178–183.
Moscow, J. A., GongM., , He, R., Sgagias, M. K., Dixon, K. H., Anzick, S. L., et al. (1995). Isolation of a gene encoding a human reduced folate carrier (RFC1) and analysis of its expression in transport-deficient, methotrexate-resistant human breast cancer cells. Cancer Res. 55 (17), 3790–3794.
Nam, J. L., VillEnEuvE, E., Hensor, E. M. A., Conaghan, P. G., Keen, H. I., Buch, M. H., et al. (2014). Remission induction comparing infliximab and high-dose intravenous steroid, followed by treat-to-target: A double-blind, randomised, controlled trial in new-onset, treatment-naive, rheumatoid arthritis (the IDEA study). Ann. rheumatic Dis. 73 (1), 75–85. doi:10.1136/annrheumdis-2013-203440
Owen, S. A., Hider, S. L., Martin, P., Bruce, I. N., Barton, A., and Thomson, W. (2013a). Genetic polymorphisms in key methotrexate pathway genes are associated with response to treatment in rheumatoid arthritis patients. Pharmacogenomics J. 13 (3), 227–234. doi:10.1038/tpj.2012.7
Owen, S. A., Lunt, M., Bowes, J., Hider, S. L., Bruce, I. N., Thomson, W., et al. (2013b). MTHFR gene polymorphisms and outcome of methotrexate treatment in patients with rheumatoid arthritis: Analysis of key polymorphisms and meta-analysis of C677T and A1298C polymorphisms. pharmacogenomics J. 13 (2), 137–147. doi:10.1038/tpj.2011.42
Pawlik, A., Herczynska, M. M., Dziedziejko, V., Malinowski, D., Kurzawski, M., Drozdzik, M., et al. (2013). Effect of allograft inflammatory factor-1 gene polymorphisms on rheumatoid arthritis treatment with methotrexate. Postepy Hig. Med. Dosw (Online) 67, 637–642. doi:10.5604/17322693.1058897
Qiu, Q., Huang, J., Shu, X., Fan, H., Zhou, Y., and Xiao, C. (2017). Polymorphisms and pharmacogenomics for the clinical efficacy of methotrexate in patients with rheumatoid arthritis: A systematic review and meta-analysis. Sci. Rep. 7, 44015. doi:10.1038/srep44015
Ranganathan, P., and McLeod, H. L. (2006). Methotrexate pharmacogenetics: The first step toward individualized therapy in rheumatoid arthritis. Arthritis rheumatism 54 (5), 1366–1377. doi:10.1002/art.21762
Saevarsdottir, S., Wallin, H., Seddighzadeh, M., Ernestam, S., Geborek, P., Petersson, I. F., et al. (2011). Predictors of response to methotrexate in early DMARD naive rheumatoid arthritis: Results from the initial open-label phase of the SWEFOT trial. Ann. Rheum. Dis. 70 (3), 469–475. doi:10.1136/ard.2010.139212
Salazar, J., Moya, P., Altes, A., Diaz-Torne, C., Casademont, J., Cerda-Gabaroi, D., et al. (2014). Polymorphisms in genes involved in the mechanism of action of methotrexate: Are they associated with outcome in rheumatoid arthritis patients? Pharmacogenomics 15 (8), 1079–1090. doi:10.2217/pgs.14.67
Salliot, C., and van der Heijde, D. (2009). Long-term safety of methotrexate monotherapy in patients with rheumatoid arthritis: A systematic literature research. Ann. Rheum. Dis. 68 (7), 1100–1104. doi:10.1136/ard.2008.093690
Smolen, J. S., and Aletaha, D. (2006). What should be our treatment goal in rheumatoid arthritis today? Clin. Exp. Rheumatol. 24 (6), S-7-13.
Smolen, J. S., Breedveld, F. C., Schiff, M. H., Kalden, J. R., Emery, P., Eberl, G., et al. (2003). A simplified disease activity index for rheumatoid arthritis for use in clinical practice. Rheumatology 42 (2), 244–257. doi:10.1093/rheumatology/keg072
Smolen, J. S., Landewe, R. B. M., Bijlsma, J. W. J., Burmester, G. R., Dougados, M., Kerschbaumer, A., et al. (2020). EULAR recommendations for the management of rheumatoid arthritis with synthetic and biological disease-modifying antirheumatic drugs: 2019 update. Ann. Rheum. Dis. 79 (6), 685–699. doi:10.1136/annrheumdis-2019-216655
Wu, C.-Y., Yang, H. Y., Luo, S. F., and Lai, J. H. (2021). From rheumatoid factor to anti-citrullinated protein antibodies and anti-carbamylated protein antibodies for diagnosis and prognosis prediction in patients with rheumatoid arthritis. Int. J. Mol. Sci. 22 (2), E686. doi:10.3390/ijms22020686
Xiao, H., Xu, J., Zhou, X., Stankovich, J., PanF., , Zhang, Z., et al. (2010). Associations between the genetic polymorphisms of MTHFR and outcomes of methotrexate treatment in rheumatoid arthritis. Clin. Exp. Rheumatol. 28 (5), 728–733.
Keywords: methotexate, rheumatoid arthritis, single nucelotide polymorphisms, clincal outcomes, association
Citation: Kolan SS, Li G, Grimolizzi F, Sexton J, Goll G, Kvien TK, Sundlisæter NP, Zucknick M, Lillegraven S, Haavardsholm EA and Skålhegg BS (2022) Identification of SNPs associated with methotrexate treatment outcomes in patients with early rheumatoid arthritis. Front. Pharmacol. 13:1075603. doi: 10.3389/fphar.2022.1075603
Received: 20 October 2022; Accepted: 02 November 2022;
Published: 17 November 2022.
Edited by:
Vijay Suppiah, University of South Australia, AustraliaReviewed by:
Ingrid Fricke-Galindo, Instituto Nacional de Enfermedades Respiratorias-México (INER), MexicoDazhi Fan, Foshan Women and Children Hospital, China
Copyright © 2022 Kolan, Li, Grimolizzi, Sexton, Goll, Kvien, Sundlisæter, Zucknick, Lillegraven, Haavardsholm and Skålhegg. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Bjørn Steen Skålhegg, b.s.skalhegg@medisin.uio.no
†These authors have contributed equally to this work and share first authorship