A longitudinal genome-wide association study of anti-tumor necrosis factor response among Japanese patients with rheumatoid arthritis

This is the first GWAS investigating genetic biomarkers of response to anti-TNF therapy
in Japanese patients with RA, utilizing a longitudinal approach to examine associations
between genome-wide SNPs and repeated measures of anti-TNF response at 3 and 6 months.
We found borderline significant association (p 1×10
^?6
) at three non-correlated regions within our study population, with the associated
SNPs mapping to or close to the following genes: MAP3K7, BACH2 (6q15), WDR27 (6q27) and GFRA1 (10q25.3). Each of these regions harbored numerous SNPs demonstrating evidence of
association with change in DAS28. Furthermore, the 6q15 locus was significantly associated
with response at 6 months (p?=?2.5×10
^?8
), and the association at the GFRA1 locus represents a replication of a previously reported association among Caucasian
patients. We therefore considered these regions worthy of being reported so that they
may be investigated further in larger datasets.

The MAP3K7 gene encodes transforming growth factor beta-activated kinase 1 (TAK1) which is a
key regulator in multiple inflammatory signaling pathways 44], 45], including the p38 MAPK and nuclear factor kappa B signaling pathways. TAK1 deficiency
leads to reduced pro-inflammatory cytokine production in cultured RA synoviocytes
46]. It is thus an excellent candidate that may influence the effect of anti-TNF agents,
as proposed 47], and is already a candidate therapeutic target to block pro-inflammatory pathways
in RA 48], 49]. Transcription factor BACH2, on the other hand, appears to be a key negative regulator
of effector T cell differentiation, promoting immune homeostasis 50]. In the mouse, it appears to be a super-enhancer repressing a network of genes critical
for T cell function 51]. Of interest, variants at the BACH2 locus have been associated with multiple autoimmune diseases, including RA 52]–57]. The GFRA1 protein is a member of the Glial cell line-derived neurotrophic factor
(GDNF) receptor family and mediates activation of the RET tyrosine kinase receptor.
GDNF is produced by astrocytes in response to pro-inflammatory cytokines including
TNF? 58] and appears to suppress interleukin-17 (IL-17)-mediated inflammation via the NF-kappa
B pathway 59]. The function of the WDR27 protein has not been established.

The association with the GFRA1 gene was previously identified by Plant et al.16]. SNP rs7070180 mapping to an intron of the GFRA1 gene, was associated with anti-TNF response in a cohort of 566 Caucasian RA patients
in the UK (p?=?2.24×10
^?4
) and in a meta-analysis including a cohort of 379 additional patients (p?=?6.42×10
^?5
). However, this locus was not reported among the major findings of that study as
SNP rs7070180 failed to genotype in one of the cohorts. While this SNP was not genotyped
or imputed in our data, several other SNPs within the 3’ untranslated region (UTR)
of the GFRA1 gene were associated with anti-TNF response among our patients. There are no reports
of MAP3K7 being associated with anti-TNF response in RA, although it has been proposed as a
good candidate for pathway pharmacogenetics relating to TNF inhibitors 47]. Among the genes of the p38 MAPK network that have been investigated 7], 9], 13], suggestive evidence of an association with MAP2K6 was reported 9], though not replicated 13]. Other associations for p38 MAPK candidate genes were reported in a sample of 1,102 patients using a generously non-stringent
significance threshold of p 0.1 7].

A major strength of the present study is the use of repeated measures of anti-TNF
response at 3 months and 6 months after treatment was started. Previous studies included
clinical response from a single time point in standard linear or logistic regression
models 5]–9], 11], 12], 16]–19]. However, assessment of response at a single point in time may not adequately reflect
a patient’s response to therapy, as response may fluctuate over time. Hence, using
response data from at least two time points is more reliable and clinically relevant.
The longitudinal approach enables the use of repeated measures of response from different
time-points for each patient, thus increasing the power to detect an association as
we have demonstrated, while taking into account within-patient correlation. Patients
with missing data at one time point were still included in the analyses, as the GEE
uses all available data. However, while the association with the MAP3K7 locus achieved genome-wide significance for SNP rs284511 when using anti-TNF response
at only 6 months, the lack of association with this SNP at 3 months led to the GEE
model only detecting a borderline significant association.

We did not identify any other overlap between our results and previous findings, possibly
due to differences in ethnicity, response variable, i.e., two time points vs a single
time point, or duration of anti-TNF treatment (3–12 months in previous studies). Further,
the lack of consistent findings between previous studies may also have been the result
of differences in a number of factors including study design, clinical outcomes examined
(DAS28ESR vs DAS28CRP), specific anti-TNF medications used, concomitant disease modifying
anti-rheumatic drugs (DMARDs), sex ratios, and small sample sizes in some cases 9], 12], 19]. Another important difference between studies is the heterogeneity in phenotypes
introduced by differences in time from baseline to assessment of response ranging
from 3 to 12 months in different studies. As seen in our data, treatment duration
is significantly associated with response, and should be given due consideration when
using multiple datasets for combined analyses to minimize phenotypic heterogeneity
or when comparing results between studies. Further, the low significance thresholds
used to identify previously reported associations may have led to false positive associations
being included. We also cannot exclude the possibility that our findings may include
false positives until they can be replicated in independent datasets or that our sample
size was not adequately powered to detect some of the previously reported findings.

The present study has a number of limitations. First, the sample size of 444 patients
is modest compared to previous GWAS, which combined data from different populations
of European ancestry to achieve large sample sizes 16], 18]. Nonetheless, it represents the largest reported sample size for studies of anti-TNF
response among Japanese or any East Asian RA population 18], 24]. Given the lack of pharmacogenomics studies examining associations with anti-TNF
response among Japanese RA patients, these results represent an important contribution
to the field. Although power may be limited due to the modest sample size, this may
in part be compensated for by the ethnic homogeneity of the patient population and
the use of response data from two time points as described earlier. Second, the response
phenotype in our study included response to three different anti-TNF agents, which
may have introduced some bias in responder status because a patient who was a non-responder
to one drug, might have responded well to another drug. For example, of the 18 patients
who had switched to a new anti-TNF drug at baseline due to inefficacy of the previous
drug, 11 had a good response to the second line of anti-TNF agent. However, all patients
starting a new anti-TNF agent at study baseline – i.e., those who switched and those
who were anti-TNF naïve – were followed for 6 months to assess response as is routinely
done in clinical practice. Hence, none of them were switched to a new agent for the
duration of the study. To mitigate misclassification bias in the response phenotype
that might have arisen from inclusion of three different anti-TNF agents in the analysis,
we adjusted for the type of anti-TNF agent used. Third, as there is no gold standard
measure to evaluate treatment response in RA, we used the DAS28CRP3 as a surrogate
to assess disease activity although it may not be a perfect measure of response. The
possibility that such a complex phenotype may be associated only with modest genetic
effects has been raised 4], 16]. A comparison of the component variables, i.e., tender and swollen joint counts and
CRP levels, showed similar trajectories from baseline to 6 months to the composite
DAS28 score. This suggests that no major differences in associations would be expected
by focusing on components of the DAS28 as the outcome in our dataset. In order to
more closely capture variations in patient response over time, we chose to use repeated
measures of the change in DAS28 from baseline. We did not categorize the outcome into
European League Against Rheumatism (EULAR) responses as this would have led to a reduction
in power. Last, the association between a genetic predictor and clinical response
could be confounded by factors that influence response. As far as possible, we adjusted
for likely confounders associated with change in DAS28 in our dataset. We did not,
however, adjust for variations in drug dosage. We had previously reported that response
to anti-TNF therapy is influenced by sex in the long term 60], but that these sex differences were not observed during the first 6 months of treatment.
Thus, the lack of an association with sex in the present dataset, which was followed
for only 6 months, is in agreement with our previous findings.