Genetic Variations in Pattern Recognition Receptor Loci Are Associated with Anti-TNF Response in Patients with Rheumatoid Arthritis

Objectives To determine whether genetic variation within genes related to the Toll-like receptor, inflammasome and interferon-γ pathways contributes to the differences in treatment response to tumour necrosis factor inhibitors (anti-TNF) in patients with rheumatoid arthritis (RA). Methods In a retrospective case-case study, we assessed 23 functional single nucleotide polymorphisms (SNPs) in 15 genes. We included 538 anti-TNF naïve Danish RA patients from the nationwide DANBIO database. Multivariable logistic regression analyses were performed to detect associations (p-value<0.05) between genotypes and European League Against Rheumatism (EULAR) treatment responses. False Discovery Rate corrections for multiple testing (q-value) and stratified analyses were performed to investigate association with individual therapies and IgM-rheumatoid factor (RF) status. Results Six of twenty successfully genotyped polymorphisms were nominally associated with EULAR treatment response. Three of these were in weak to moderate linkage disequilibrium with polymorphisms previously reported associated with anti-TNF treatment response. TLR5(rs5744174) variant allele carriers (odds ratio(OR) = 1.7(1.1–2.5),p = 0.010,q = 0.46) and TLR1(rs4833095) homozygous variant carriers (OR = 2.8(1.1–7.4),p = 0.037,q = 0.46) had higher odds for a positive treatment response. NLRP3(rs10754558) variant allele carriers (odds ratio(OR) = 0.6(0.4–1.0),p = 0.045,q = 0.46) were more likely to have a negative treatment response. The association in TLR5(rs5744174) remained significant after correction for multiple comparisons among patients negative for RF (OR = 6.2(2.4–16.3),p = 0.0002,q = 0.024). No other association withstood correction for multiple testing. Post hoc analyses showed that change in Patient Global score on a visual analogue scale (VAS) and change in pain VAS were the main factors responsible for the association. Conclusions We reproduced previously reported associations between genetic variation in the TLR10/1/6 gene cluster, TLR5, and NLRP3 loci and response to anti-TNF treatment in RA. Changes in VAS pain and patient global scores were the main contributors to the association found for TLR5. Furthermore, we identified other candidate genes that require replication in independent cohorts.

The association in TLR5(rs5744174) remained significant after correction for multiple comparisons among patients negative for RF (OR = 6.2 (2.4-16.3),p = 0.0002,q = 0.024). No other association withstood correction for multiple testing. Post hoc analyses showed that change in Patient Global score on a visual analogue scale (VAS) and change in pain VAS were the main factors responsible for the association.

Introduction
Predictive biomarkers capable of stratifying patients into responders and non-responders to treatment with tumour necrosis factor inhibitors (anti-TNF) will enable selection of the optimal treatment for the individual patients and thereby improve patient care. Such biomarkers, however, are scarce and none are applicable in a clinical setting. Some single nucleotide polymorphisms (SNPs) may show associations with anti-TNF response and might be useful for prediction but despite several studies addressing this issue, only few associations have been replicated (e.g. the PTPRC [1] and PDE3A-SLCO1C1 loci [2]). Overall, associations between genetic polymorphisms and anti-TNF response are characterized by small effect sizes.
However, together with e.g. expression profiling, epigenetic, para-clinical, or clinical markers the growing number of polymorphisms associated with anti-TNF treatment response may be potentially useful for prediction [3]. Furthermore, associations between genetic variation and treatment outcome may provide insight into aberrant molecular pathways and thus form the basis for developing new treatment strategies.
We recently published a candidate gene study [4] focussing on genes involved in the NF-κB mediated signalling pathway. Polymorphisms associated with anti-TNF response were found in the genes encoding NACHT, LRR and PYD domains-containing protein 3 (NLRP3/NALP3) (rs4612666) and interferon-γ (IFN-γ) (rs2430561). The association with NLRP3 is supported by a recent study of British patients [5].
The NLRP3 inflammasome belongs to a group of intracellular innate immune sensors, which generally sense stimuli associated with infection and stressed tissue [6]. Caspase recruitment domain-containing protein 8 (CARD8) is involved in NLRP3-inflammasome formation [6] and polymorphisms in this gene have been associated with response to anti-TNF drugs [5]. Upon inflammasome activation, caspase-1 proteolytically activates the pro-inflammatory cytokines interleukin (IL)-1β and IL-18. Like the NLRP3-inflammasome, Toll-like receptors (TLRs) are innate immune sensors. They sense exogenous and endogenous antigens and activate pathways (e.g. NF-κB, AP-1) that increase production of pro-inflammatory cytokines. A recent well-powered study found polymorphisms in the TLR signalling pathways (TLR2, TLR4, MyD88, CHUK) associated with response to anti-TNF treatment [7]. Also, polymorphisms in the genes encoding TLR-1 and TLR-5 (rs5744174) have been shown to be associated with increased IFN-γ secretion [8,9], and, for the latter, with reduced expression of IL-6 and IL-1β [10].
Interferon-γ (IFN-γ) is a central cytokine for both adaptive and innate immunity that mediates downstream signalling through binding to the heterodimeric IFN-γ receptor (IFNGR1 and IFNGR2), which in turn activates the JAK-STAT (Janus kinase Signal Transducers and Activator of Transcription) pathway [11]. STAT4 regulates proliferation, survival, and differentiation of lymphocytes and genetic variation in the corresponding gene has been associated with RA [12].
Thus, our aim in this study was to assess additional functional polymorphisms related to the TLR, inflammasome and interferon-γ pathways and thereby extend the analysis of associations found in our previous study [4]. We analysed 23 functional polymorphisms in CARD8,  IFNGR1, IFNGR2, IL12B, IL12RB1, IL12RB2, IL17A, IL18, JAK2, NLRP1, NLRP3, TBX21,  TIRAP, TLR1 and TLR5 in a Danish cohort of 538 RA patients treated with anti-TNF and prospectively monitored for treatment effect using validated clinical scoring methods.

Ethics statement
The study was conducted in accordance with the Declaration of Helsinki and was approved by the Regional Ethics Committee of Central Denmark Region (M-20100153 and S-20120113) and the Danish Data Protection Agency (J. 2010-41-4719). The Regional Ethics Committee of Central Denmark Region gave exemption from obtaining informed consent because samples were taken for other reasons and data were analysed anonymously.

Patients and Samples
We included 538 anti-TNF naïve RA patients. All patients initiated their first anti-TNF treatment, had clinical variables registered at baseline and follow-up in DANBIO (The National Danish Registry for Biological Treatment of Rheumatic Diseases), and had available blood samples. The baseline (pre-treatment) visit was defined as a visit 0-30 days before onset of anti-TNF treatment and follow-up as a visit within 60-180 days after treatment onset (at the contact closest to 120 days if more than one visit was registered). A detailed description of the DANBIO registry and the additional patient data were published previously [4,13].

Candidate gene analysis
PubMed was searched for functional polymorphisms in genes closely related to the up-and down-stream signalling molecules of the NLRP3 inflammasome and interferon-γ [4]. A total of 23 SNPs were ultimately chosen primarily based on evidence of biological effect and secondly based on documented association with autoimmune disease (S1 Table). Expected minor allele frequencies (MAFs) ranged from 14% to 44% except for one polymorphism (rs11810249, MAF: 4.4%).
The polymorphisms were genotyped by PCR-based KASP genotyping assay by LGC Genomics (Middlesex, United Kingdom-www.lgcgenomics.com) on DNA extracted from blood (Maxwell 16 LEV Blood DNA Kit, Promega, Madison, Wisconsin, USA) as described by Bank et al. [14]. Genotyping failed for rs8134145 (IFNGR2) and rs2072493 (TLR5) due to their close proximity to neighbouring genotyped SNPs, and for IL17A for unknown reasons. All other chosen assays had a call rate exceeding 95%. Repeated genotyping of 94 randomly selected samples in an inflammatory bowel disease cohort yielded >99% identical calls.

Statistical methods
Primary outcome (good, moderate or none) was defined by the European League Against Rheumatism (EULAR) response criteria [15] at follow-up. We performed the analyses by comparing either EULAR good/moderate vs. none (G/M vs. N) response or EULAR good vs. moderate/none (G vs. M/N) response. We also analysed treatment response defined by the American College of Rheumatology outcome measure, ACR50 response [16] and relative change in disease activity score for 28 joints (DAS28) (relDAS28 = (baseline DAS28 − follow-up DAS28)/baseline DAS28) to enable comparison with other studies.
Multivariate logistic regression analyses were performed to assess association between genotype and anti-TNF treatment response at the 5% significance level. Adjustments were made in the multivariate analyses for sex, baseline health assessment questionnaire (HAQ) score, DAS28 and concomitant treatment with disease modifying anti-rheumatic drug (DMARD). We have also performed stratified analyses for patients positive for IgM-RF (seropositive RA patients), for the specific type of anti-TNF drug (not certolizumab and golimumab due to insufficient power) and for a combined group of patients treated with monoclonal antibodies (infliximab, adalimumab, and golimumab).
We performed correction for multiple testing using False Discovery Rate classical one-stage method set at 0.05 (q-value) [17]. Values were based on the number of test performed in each of the primary, secondary and stratified analyses, respectively (see S1-S6 tables). Values can be interpreted as the expected proportion of false positive results at least as extreme as the observed result.

Study population
Baseline characteristics and treatment response of the study population are presented in Table 1. Seronegative RA patients were significantly younger and had higher pain and global scores as well as higher tender joint count compared with the seropositive patients. No difference in response rates to anti-TNF was observed between the two groups. Both had 29-30% non-responders. Patients with RA receiving concomitant synthetic DMARDs (sDMARD) comprised 84%.
Clinical data were collected at variable time-points due to the study design. Eighty-seven per cent of the patients had baseline data collected 7 days or less before treatment onset, and 74% of the follow-up data were registered 16±4 weeks after treatment onset (S1 Fig).

Genotype associations with anti-TNF response
We assessed the association with EULAR anti-TNF response regardless of RF status among all patients. We found that variant allele carriers of the TLR5 rs5744174 (OR = 1.7, p = 0.010), IL18 rs187238 (OR = 1.5, p = 0.026) and the homozygous variant genotype of TLR1 rs4833095 (OR = 2.8, p = 0.037) had a positive treatment response (EULAR G vs. M/N) ( Table 2). With EULAR G/M vs. N as treatment response, carriers of another variant in IL18 (rs1946518; r 2 = 0.59) were more likely to have positive treatment response (OR = 1.6, p = 0.022). Furthermore, individuals heterozygous for IL12B rs6887695 (OR = 0.6, p = 0.017) and variant allele carriers of NLRP3 rs10754558 (OR = 0.6, p = 0.017) seemed to have a less favourable treatment response (EULAR G/M vs. N). Analyses of associations with treatment response outcome defined by ACR50 and relDAS28 are shown in S2 Table. Association with altered ACR50 treatment response was found for TLR1 rs4833095 and with altered relDAS28 treatment responses for NLRP3 rs10754558 and TLR5 rs5744174.
We also assessed the functionally relevant haplotypes in IFNGR2 and IL12B (S1 Table) for association with EULAR anti-TNF response (G/M vs. N and G vs. M/N) under the same multivariate model. No associations between risk-haplotypes and anti-TNF response and no interaction with smoking were found.
Stratified analyses. Next, we performed stratified analyses of polymorphisms for associations with treatment response. For seropositive RA patients (n = 407), IL12B rs6887695, IL18 rs187238 and NLRP3 rs10754558 were associated with similar effect estimates and significance levels (S3 Table) as the total cohort. For IL18 rs1946518 and TLR1 rs4833095 effect estimates were similar (p-value <0.1) and for TLR5 rs5744174 effect estimates were substantially lower and non-significant compared to the overall cohort (S3 Table).
Stratified analyses by type of anti-TNF drug (infliximab (n = 168), etanercept (n = 166), adalimumab (n = 134) and the combined group of monoclonal antibodies (infliximab, adalimumab, and golimumab) (n = 372)) are shown in S4 Table. Variant allele carriers of TLR5 rs5744174 treated with infliximab (OR = 2.7, p = 0.018) and etanercept (OR = 2.3, p = 0.022) were more likely to have a good treatment response (EULAR G vs. MN) (S4 Table) and a clear allele-dosage effect was observed. Altered treatment responses were also seen among variant allele carriers of NLRP1 (rs878329) and NLRP3 rs10754558 treated with infliximab and among variant allele carriers of IL12B rs6887695 treated with the monoclonal anti-TNF drugs (S4 Table).
Post hoc analyses. Since there was some discrepancy between associations found in the primary analyses and stratified analyses of seropositive patients, we also assessed the association in the seronegative subgroup of RA patients for TLR5 rs5744174 (S3 Table). Seronegative carriers of the variant allele of TLR5 rs5744174 had a greater chance of good treatment response (EULAR G vs. M/N) (OR = 6.2 (2.4-16.3), p = 0.0002, q = 0.024) (Fig 1).
Due to lack of association among seropositive RA patients with the TLR5 polymorphism, we explored other clinical variables of importance for this association using linear regression analyses of the association between the TLR5 polymorphism and the components of DAS28 (change (Δ) between baseline and follow-up of swollen and tender joint counts (SJC and TJC), CRP and patient global score).
In this cohort of RA patients, there were marked clinical differences between the seropositive and seronegative subgroups (fraction of smokers, erosive status, age, baseline tender joint count, pain and patient global score), and this could potentially cause confounding. However, additional multivariate regression analysis adjusting for these potential confounders did not markedly change the significance levels or effect sizes for the RF status-stratified analyses of CARD8 (rs2043211), IFNGR2 (rs17882748), IL12B (rs6887695), IL18 (rs187238, rs360719), TLR5 (rs5744174) or TLR1 (rs4833095) (S6 Table).

Discussion
In this study, we extended the analyses of polymorphisms in genes related to inflammasomeand interferon-γ pathways. Overall, we successfully genotyped 20 of 23 functional polymorphisms, and six SNPs in IL12B, IL18, NLRP3, TLR1 and TLR5 (Table 2) were nominally associated with EULAR treatment response in RA patients treated with anti-TNF. Furthermore, the TLR5 polymorphism lacked association among the seropositive patients but was strongly associated with response among the seronegative patients.
We chose to analyse at the 0.05 significance level. At this level, we had !80% chance to detect effect sizes between 1.3 and 1.5 assuming an additive genetic model, and effect sizes between 1.6 and 2.9 assuming a recessive model (MAF between 0.1-0.4). Since multiple testing increases the risk of false positive associations, we have also presented associations corrected for multiple testing (False Discovery Rate, q-value). The association found for TLR5 among seronegative RA patients was the only that withstood correction for multiple testing.
Since the present study is exploratory, replication in independent cohorts is needed, in particular, for those found in the stratified analyses. However, the selection of functional polymorphisms in genes with a biologically plausible role in RA pathogenesis increases the prior probability of true associations. The presented work does not provide evidence for causal relationships between the studied polymorphisms and the altered anti-TNF response because the associated polymorphisms could also be proxies for other genetic markers.
It should be noted that, due to the design of the study, the time-points of collected response data vary between 60 and 180 days from treatment onset. When we analysed the subset of patients with response data at 16±4 weeks (comprising 74% of the patients) in the overall cohort, associations found for IL12B rs6887695, IL18 rs1946518, and the SNPs in TLR1 and TLR5 remained nominally associated with similar or higher odds ratios (S5 Table).
The contrasting odds ratios for IL12B rs6887695 heterozygous and homozygous variant genotypes substantially lower the likelihood of a true association for this polymorphism.
The present study reproduces associations between three gene loci (TLR10/1/6 gene cluster, NLRP3 and TLR5) (Fig 2) and RA anti-TNF response previously found in independent studies [4,5,7,18]. TLR1 rs4833095 is located within a cluster of genes encoding TLR-1, TLR-6 and  TLR-10 in a 54-kb region on chromosome 4. TLR1 rs4833095 is in modest linkage disequilibrium with TLR10 rs11096957 (R 2 = 0.52, D' = 0.94), which was found associated with EULAR GM vs. N response to anti-TNF treatment in a UK cohort (N = 909) [7]. TLR1 rs4833095 has been associated with a high peripheral blood mononuclear cell (PBMC) TLR-1 cell surface expression [19]. Here, we found the homozygous variant genotype associated with a greater chance of good treatment response to anti-TNF. However, this finding relies on small patient groups and needs confirmation in an independent and preferably larger cohort.
To our knowledge, this study is the third to show an association between genetic variation in NLRP3 and anti-TNF treatment response. Mathews et al. found three polymorphisms in NLRP3 with nominal association with EULAR response with rs4925659 as the most tightly linked with rs10754558 (D' = 0.93; r 2 = 0.40) [5]. In the same cohort as used in the present study, we have previously reported the rs4612666 T allele to be associated with a lesser chance of EULAR GM vs. N response [4]. The two SNPs are not in linkage disequilibrium (D' = 0.15; r 2 = 0.06), thus not allowing haplotype analysis.
Both rs4612666 and rs10754558 in NLRP3 have been functionally investigated by Hitomi et al. [20]. In THP-1 cells they found minor alleles of these two polymorphisms to decrease and increase NLRP3 expression, respectively. Based on this, our results for these two polymorphisms seem contradictory since the high expression alleles of the two polymorphisms are associated with positive and negative treatment response, respectively. However, the in vitro experiments performed on THP-1 cells may not reflect in vivo conditions in RA patients.
We found two polymorphisms in IL18 associated with anti-TNF treatment response. IL-18 is an interleukin that is activated by the NLRP3 inflammasome, and these results therefore further support a role of this pathway for anti-TNF treatment outcome in RA.
The association found for TLR5 confirms the association between a TLR5 gene locus and EULAR anti-TNF response previously found in a Dutch cohort (n = 182) [18]. In that study, both rs5744174 and another TLR5 polymorphism (rs2072493) in linkage disequilibrium (genetic distance: 71bp; D' = 1; r 2 = 0.08) were genotyped (Fig 2). Only rs2072493 was associated with anti-TNF response. In a Swedish validation cohort (n = 269), these findings could not be confirmed. Due to different MAFs of these two polymorphisms, they have low tagging capabilities of each other (r 2 = 0.08). The lack of associations for rs5744174 in the Dutch cohort and for rs2072493 in the Swedish cohort could be caused by the low statistical power, causing false negative results. Differences in statistical methods may also play a role. Coenen et al. performed unadjusted Fisher's exact or Chi-square test for analyses of 3x3 tables (EULAR response groups / genotype counts), and analyses of our data are non-significant when analysed this way (data not shown). Further, differences in anti-TNF drugs received in the two cohorts may also be of importance. As an example, the fraction of adalimumab-treated patients is larger in the Dutch cohort (34% vs. 25%). Our data do suggest a trend for a negative adalimumab response for carriers of the TLR5 variant allele compared to a positive treatment response found in the overall cohort and among infliximab and etanercept treated patients (S4 Table).
Functional studies of TLR5 rs5744174 have shown that flagellin (a TLR5 ligand) stimulation results in lower chemokine (C-C motif) ligand (CCL20) production in HEK 293T cells [21] and reduced IL-6 and IL-1β expression in primary immune cells from healthy homozygous carriers of the variant allele [10]. Another study reported that anti-TNF treatment decreases TLR5 (and TLR4) expression among ankylosing spondylitis patients [22], which supports the hypothesis of TLR5 having a role in the differential anti-TNF response. Interestingly, post hoc analyses showed that change in patients' subjective general health assessment (Δ patient global score) and probably also change in pain score could explain most of the correlation with the multi-component EULAR response. A recent paper by Coenen et al. showed that separate components of DAS28, such as tender and swollen joint counts, have a stronger heritability in pharmacogenetic studies of anti-TNF response than DAS28 itself [23]. Although they did not show this for patient global score assessment, it may be the case for TLR5 rs5744174. In fact, it is well established that TLRs are central in chronic pain pathogenesis [24] and TLR2/3/4/5 deficient mice have reduced neuropathic pain [25]. In a future validation study of the predictive value of this polymorphism, it will therefore be very relevant to include analyses adjusting for patient global score and RF status.
Overall, the effects of the associated polymorphisms were too weak to be used as independent biomarkers for treatment response but they have a potential use as part of a larger panel of predictive biomarkers if validated in independent cohorts.
Based on the results from this study, we hypothesize that so far undefined RA sub-phenotype(s) not responding to anti-TNF treatment are characterized by increased NLRP3-inflammasome activation and increased IL-1β/IL-18 production.