Genetics in TNF-TNFR pathway: A complex network causing spondyloarthritis and conditioning response to anti-TNFα therapy

Objectives We investigated whether polymorphisms (SNPs) in the promoter region of TNFA, or in the autoinflammatory TNFRSF1A and MEFV genes, concur with HLA-B27 in enhancing the risk of Spondyloarthritis (SpA) and/or in predicting the response to anti-TNFα treatment. Methods 373 controls and 137 SpA (82 with Psoriatic Arthritis-PsA and 55 with Ankylosing Spondylitis- AS; 98/137 under TNFα inhibitor therapy) from the Veneto Region (Italy) were studied. TNFA polymorphisms (-1031T>C;-857C>T;-376G>A;-308G>A;-238G>A) and HLA-B27 were assayed by RT-PCR. Direct sequencing of MEFV (exons 2,3,5 and 10) and TNFRSF1A (exons 2,3,4 and 6) genes were performed. Results HLA-B27 was associated with AS (χ2 = 120.1; p = 0.000). Only the TNFA -1031T>C was singly associated with SpA, and the haplotype C/G, resulting from -1031T>C/-308G>A combination, was significantly associated with a reduced risk of SpA (OR: 0.67, CI: 0.46–0.97; p = 0.035). Two SNPs were identified in TNFRSF1A, the R92Q (Minor allele frequency-MAF = 0.034) and c.625+10A>G (MAF = 0.479). None of them was associated with SpA (p>0.05). The TNFRSF1A c.625+10 G allele was associated with late response to anti-TNFα therapy (p = 0.031). Twenty-one SNPs were identified in MEFV gene, 10 with a known potential functional significance. Variant alleles were extremely rare in our population (MAF<0.025) except for R202Q (MAF = 0.27). None was associated with SpA diagnosis (p>0.05). Conclusion TNFRSF1A and MEFV gene SNPs are not associated with SpA in the North-East of Italy. AS risk appears to depend not only on HLA-B27, but also on the protective TNFA haplotype -1031C/-308G. The TNFRSF1A c.625+10A>G impacts on the response to anti-TNFα therapy.


Introduction
The spondyloarthritis (SpA), a group of chronic inflammatory diseases affecting the sacro ileal joint that includes the two major subtypes Ankylosing Spondylitis (AS) and Psoriatic Arthritis (PsA), results from a complex interplay among genetic background and environmental factors, that leads to the activation of autoinflammation and dys-regulation of the immune-system. Two main groups of SpA have been defined, axial and peripheral SpA, the latter being characterized by the involvement of less than five peripheral joints occurring before, at the same time or after the onset of axial SpA [1,2].
At onset of disease, usually occurring at a young age, SpA diagnosis might be delayed because of the overlap of clinical features and absence of specific blood biomarkers. Although increased serum levels of C-reactive protein (CRP), erythrocyte sedimentation rate (ESR), cytokines, vascular endothelial growth factor (VEGF), calprotectin and matrix metalloproteinase-3 (MMP3) might be observed in SpA patients, these biomarkers have limited specificity and sensitivity, because they are common markers of systemic inflammation and their increase parallels disease activity [3,4]. The delayed SpA diagnosis, estimated to be 8-11 years, underlies a retard in treatment and imparts a tremendous symptomatic burden and loss of function in these patients during the productive years of life [5].
The delay in diagnosis, the lack of specific diagnostic and prognostic biomarkers and of complete clinical response even after therapy with tumor necrosis factor α (TNF-α) inhibitors, such as etanercept, infliximab, adalimumab and golimumab, have raised the interest in the pathogenetic mechanism involved in SpA diseases, with the aim also of identifying biochemical and genetic biomarkers that might help the diagnostic work up and the evaluation of treatment effectiveness [6].
Genes in major histocompatibility complex (MHC) and non-MHC locus have been claimed to play a role in pathogenesis of SpA. The MHC accounts for about 40-50% of the genetic risk of developing SpA, 30% of the total heritability being related to the human leukocyte antigen (HLA)-B27 haplotype [7]. Single nucleotide polymorphisms (SNPs) in genes involved in TNFα signalling, as those in the promoter region of TNFA gene, in TNFSF15 (TNF ligand superfamily, member 15), TNFR-1 (TNF receptor 1) and TRADD (TNF receptor type 1-associated death domain protein), have been identified as potentially associated with SpA [8][9][10][11][12][13]. This is of great interest because TNF-α is not only involved in the propagation and perpetuation of inflammation in SpA, but also because of the clinical efficacy of treatments based on drugs targeting the TNF-α pathway [14,15]. Taking into consideration the potential role of autoinflammation in disease pathogenesis, other potential candidates are the MEFV (Mediterranean fever) and TNFRSF1A (TNF receptor superfamily member 1A) genes, involved in the pathogenesis of the autoinflammatory disorders Familial Mediterranean Fever (FMF) and Tumor necrosis factor Receptor-Associated Periodic Syndrome (TRAPS) respectively [16,17]. The latter gene encodes the TNFR-1 and genetic variations might underlie functional alterations of the TNF-α signalling, while variations in MEFV gene were shown to be very common among FMF Turkish patients with AS as a sole clinical manifestation [18][19][20][21][22][23][24].
The aim of the present study, conducted in 137 SpA patients and 373 controls coming from the North-East of Italy, was to verify whether in addition to the established HLA-B27 genetic predisposing factor, biomarkers of inflammation, SNPs in the promoter region of TNFA, or variants of the autoinflammatory TNFRSF1A and MEFV genes, might be of help in SpA diagnosis and/or in predicting the response to anti-TNFα treatment.

Ethic statement
This case-control study was approved by the Institutional Review Board of the University-Hospital of Padova (Protocol number: 3024P/13) and all subjects gave their fully informed written consent at enrolment.

Study population
A total of 510 cases and controls, coming from the Veneto Region, a North-East Italian region, were enrolled from January 2014 to December 2017. Cases were 137 SpA patients (84 males, 53 females; mean age ± SD: 51.6 ± 12.6 years) and controls were 373 healthy blood donors (230 males, 143 females; mean age ± SD: 46.5 ± 9.8 years).
Among cases, 55 patients (40%) had a definite diagnosis of AS, according to the modified New York criteria [25] and 82 (60%) had a diagnosis of PsA according to the ClASsification criteria for Psoriatic ARthritis (CASPAR criteria) [26].
All patients underwent clinical, clinimetric and functional examinations following a standardized protocol. Detailed demographic and clinical characteristics of the enrolled patients are shown in Table 1.
At enrolment, 98 patients (71.5%) were under treatment with TNF-α inhibitors (32 Infliximab, 38 Adalimumab, 21 Etanercept, 4 Golimumab, 2 Ustekinumab and 1 Secukinumab) and 16/98 (16.3%) experienced in their clinical history, at least one switch of TNF-α inhibitor before entering the study. Although there is no consensus regarding the criteria to classify the response to TNF-α inhibitors, the European League Against Rheumatism (EULAR) stated that maximal efficacy usually is not seen before 6 months even when improvement is achieved at 3 months, and indicates monitoring every 6-12 months once the treatment target has been attained [27]. Based on these recommendations, we considered 10 months-evaluation an adequate time point to make a final classification of patients as early or late responders, being responders those who reached a Bath Ankylosing Spondylitis Disease Activity Index (BAS-DAI) score lower than/equal to 4, on a scale from 0 (no activity) to 10 (high disease activity), in agreement with Braun et al [28].

Haematological and biochemical indices
The complete blood count, ESR, CRP, glucose, uric acid, creatinine, alanine aminotransferase (ALT) and prealbumin were measured by standard laboratory methods.
All genetic analyses were performed in an exploratory cohort of 91 patients and 27 controls. To improve the estimate of significant associations found at interim analysis in the exploratory cohort, the number of examined controls was extended from 27 to 214 for the HLA-B27 and to 218 for the R202Q SNP in MEFV gene, the maximum number of available biomaterial. The

TNF-α mRNA expression analysis
Peripheral blood mononuclear cells (PBMCs) obtained from 95 controls, 19 AS and 26 PsA were isolated by gradient centrifugation (Histopaque, 1077, Sigma-Aldrich, Milano, Italy). Total RNA from PBMCs was isolated according to the manufacturer's instructions (High Pure RNA Isolation Kit, Roche, Germany). Five hundred ng of total RNA were reverse transcribed into cDNA (Random primers and Superscript TM II -Reverse Transcriptase, Invitrogen, Italy). Relative quantification of TNFα was undertaken by RT-PCR with ABI Prism 7900 HT (Applied Biosystems, CA, USA) using the primers and probe set supplied by ThermoFisher Scientific (Waltham, USA, catalogue number 4331182, assay ID Hs00174128_m1). Relative quantifications were performed in a final volume of 20 μL with 33 ng of cDNA. The reference gene, HPRT1 (Hypoxanthine-guanine phosphoribosyltransferase), was selected according to the method commonly used for internal control for quantitative gene expression analyses, and its expression was determined by commercially available HPRT1 primers and probe sets (Part Number 4326321E, Applied Biosystems, CA, USA). TNF and HPRT1 were analyzed in duplicate for each sample. PCR was run at 2 minutes at 50˚C, 10 min at 95˚C, followed by 40 cycles of 15 seconds at 95˚C and 1 min at 60˚C. To determine the relative mRNA expression levels of TNF gene we used the comparative Ct method, a mathematical model that calculates changes in gene expression as a relative fold difference between an experimental and a pool derived from 20 healthy blood donors used as calibrator sample.

Statistical analysis
The chi-square (χ 2 ) test, Fisher's exact test, binary logistic regression analysis, Student's t test for unpaired data, analysis of variance (ANOVA) and Bonferroni's adjustment of p value for multiple testing were performed using Stata software, version 13.1 (StataCorp, Lakeway Drive, TX, USA). Genotype frequencies were tested for Hardy-Weinberg equilibrium proportions using χ 2 test. Haplotypes phases and frequencies were estimated by the retrospective profilelikelihood approach using the Stata haplologit package. Odds ratio (OR) was calculated by logistic regression analysis using the additive model [31].

Biomarkers of inflammation in AS and PsA
Based on the premise that inflammation of the sacroiliac joint is the hallmark of SpA, we first verified whether this confined inflammatory process could induce alterations in haematological and biochemical inflammatory indices.
The results found at enrolment in cases and controls are reported in Table 2. Both AS and PsA were associated with haematological and biochemical signs of systemic inflammation, being the number of white blood cells (WBC), polymorphonuclear cells (PMN) and platelets (PLT), and serum CRP levels increased in patients with respect to controls. PsA was also correlated with increased levels of glucose and ALT.

HLA-B27 is associated with AS
The genetic basis of SpA is complex, but the role of HLA-B27 haplotype is well established. For this reason we first evaluated in a series of 214/373 controls and 91/137 patients the HLA-B27 haplotype, results being presented in Table 3. As expected, SpA cases carrying the HLA-B27 were more frequent than controls (χ 2 = 120.12, p< 0.0001), being this haplotype mainly associated with AS than with PsA.

TNFA -1031C/-308G haplotype reduces AS risk
A number of variants in genes encoding for inflammatory and immunomodulatory cytokines and for their receptors could be involved in SpA. In this study we investigated the presence of association between SpA and five SNPs in the promoter region of TNFA. These SNPs regulate singly or combined in haplotypes the transcription of TNFα, one of the main cytokines involved in SpA pathogenesis and target of therapy with TNFα inhibitors.
The genotypes and Minor Allele Frequency (MAF) resulting from the five studied SNPs of the TNFA gene and found in patients (AS and PsA) and controls are reported in Table 4. All polymorphisms were in Hardy-Weinberg equilibrium (p >0.01) and only the -1031T>C SNP was significantly correlated with disease diagnosis, with a lower number of AS cases carrying the TC or the CC genotypes.
The studied polymorphisms of the TNFA gene, all closely located in the promoter region, are carried as different haplotypes. The haplotypes and their frequencies in controls resulting from the pairwise combination of the five SNPs were inferred by statistical analysis (S2 Table).
To verify whether TNFA haplotypes exert any independent role over SpA diagnosis, logistic regression analyses were performed, on considering SpA as the outcome variable and TNFA haplotype combinations as predictors. Table 5. reports the OR with 95% confidence intervals (CI) with respect to the reference haplotype (Ref), which, for any combination, was the most frequent in controls.
The haplotype resulting from the combination of the rare TNFA -1031 C allele with the common TNFA -308 G allele was significantly associated with SpA with an OR lower than 1.00, suggesting its protective role. When SpA were considered singly, the TNFA -1031C/-308G haplotype was confirmed to be associated with AS (OR: 0. By selecting the two SNPs -1031T>C and -308G>A of the TNFA gene, the six following genotypes resulting from haplotypes combinations, were inferred: TG/TG, TG/TA, TG/CG, TA/TA, CG/TA, CG/CG. Since the most frequent TNFA haplotype was -1031T/-308G while the protective TNFA haplotype was -1031C/-308G, genotypes were grouped as follows: common homozygous (TG/TG), carriers of the protective haplotype (CG/TG, CG/TA and CG/ CG), and others (TG/TA and TA/TA). Subjects carrying the protective TNFA -1031C/-308G   haplotype were less frequently found among AS (10%, χ 2 = 12.674 and p = 0.002), not among PsA (26%, χ 2 = 0.213 and p = 0.899), with respect to controls (49%). TNFα expression in PBMCs obtained from a series of 95 controls was found to be significantly reduced among carriers of the protective TNFA -1031C/-308G haplotype (Fig 1;
Among the ten identified missense SNPs, nine were extremely rare in our population, being their MAF lower or equal to 0.025. Only the R202Q SNP (rs224222) had a MAF of 0.288 and the rare homozygous AA genotype was found only among cases (AS frequency = 0.09; PsA frequency = 0.10), not in controls (χ 2 = 7.4837; p = 0.006) (S4 Table). To improve the estimate of this association, the number of examined controls was extended from 27 to 218 reaching a ratio between cases and controls close to 1:2. By doing this, the association between the R202Q SNP and SpA was not confirmed (χ 2 = 1.05, p = 0.59) (S5 Table).

TNFRSF1A c.625+10A>G SNP predicts response to anti-TNFα treatment
Therapy with TNFα inhibitors although effective, may fail in a subset of SpA patients. Based on the hypothesis that a different genetics background underlying the TNFα pathway may condition the clinical response, we verified whether HLA-B27, TNFA and TNFRSF1A genetics correlate with response to TNF-α inhibitors. The analysis was performed considering a subset of 65 patients who were under treatment with TNF-α inhibitors and for whom the aforementioned genetic data were available. Thirty six patients had early remission, i.e. BASDAI score lower than/equal to 4 in a time ranging from 1 to 10 months (mean ± SD: 5.94 ± 1.19 months; median: 6 months), while 29 patients were late responders, i.e BASDAI score lower than/equal to 4 in a time ranging from 10 to 36 months (mean ± SD: 16.55 ± 7.16 months; median: 12 months).
Only the TNFRSF1A c.625+10A>G SNP was associated with response to therapy, the rare GG genotype being more frequently found among late than early responders (Table 6), the association being confirmed among patients treated with infliximab (Fisher's exact = 0.026) but not in those treated with adalimumab (Fisher's exact = 1.000). The number of patients treated with the other TNFα inhibitors was too low to support statistical analysis and these drugs were therefore not considered in separate analyses.

Discussion
The SpA are a group of chronic inflammatory diseases characterized by the sharing of distinctive pathological, clinical and radiographic aspects and a strong genetic predisposition [1]. Since no specific biomarkers of diagnosis and prognosis have been validated for the clinical assessment of these diseases [5,6], we verified whether inflammatory biochemical markers commonly used in clinical practice and genetic variants of HLA locus, TNF-TNFR pathway and MEFV play a part in increasing the risk of the most prevalent SpA subtypes, AS or PsA, and/or in predicting the response to therapy with TNF-α inhibitors.
The two SpA subtypes, as expected, were well distinguished by characteristic clinical manifestations. Family history of disease, personal history of psoriasis, dactylitis and peripheral arthritis were more frequently found among PsA, and this fit well with the CASPAR criteria [26]. Back pain and buttock pain, although commonly shared by the two disease groups, were more frequent among AS than PsA, while some extra-articular clinical manifestations were extremely rare in both groups.
In both AS and PsA, WBC, PMN, PLT and serum CRP were higher than in controls, although often comprised within the reference intervals. These minimal differences in systemic signs of inflammation between cases and controls might depend on treatment (all patients were on TNF-α inhibitors or nonsteroidal anti-inflammatory drugs-NSAIDs/disease-modifying antirheumatic drugs-DMARDS) and indicate that the systemic inflammatory process, underlying these diseases, cannot be completely switch off [32].
The main genetic determinant of SpA, HLA-B27, was confirmed in the present study to be strictly associated with AS, and less strictly with PsA [33]. Located in proximity of the HLA genes, the TNFA gene appears of potential interest in SpA. Five main polymorphisms in the promoter region of this gene have been described in Caucasians and reported to be associated with different inflammatory and autoimmune diseases [34,35]. Considering the specific SpA setting, data in the literature are not definitive and sometimes contrastive, not only because of differences in ethnicity, number and selection criteria of the studied populations, but also because these SNPs regulate TNFA gene transcription mainly acting as haplotypes, which were verified in this study. Only TNFA -1031T>C SNP was associated with SpA diagnosis, while the other studied SNPs were not, this finding being in agreement with previous data [8]. The haplotypes deriving from the pairwise combinations of the five studied SNPs were inferred and only the haplotype -1031C/-308G was significantly associated with SpA, exerting in these diseases a protective role, which was mainly for AS and confirmed also among HLA-B27 positive patients. This protective role might be explained taking into account that any single allele concurring in determining the studied haplotype was demonstrated to be associated with reduced TNFα release. In fact, the rare -1031C allele, in complete linkage disequilibrium with the -863A allele, associates with low TNFα production [36], and the -308G common allele associates with a lower spontaneous or stimulated TNFα release both in vitro and in vivo with respect to the rare -308A allele [37,38]. The hypothesis that TNFA -1031C/-308G haplotype might determine a reduced TNFα transcription was confirmed in this study by the analysis of mRNA expression levels in PBMCs obtained from healthy subjects. However, in the presence of diseases, the dampened TNFα expression associated with TNFA -1031C/-308G haplotype was not seen. It is not surprising because of the limited sample size and the complicate situations of patients, such of treatment of TNFα inhibitors. Based on the above data, it is reasonable to suggest that the TNFA -1031C/-308G haplotype exerts a protective role for SpA, but mainly for AS, in HLA-B27 positive subjects by limiting the transcription of TNFα.
In the complex scenario of TNFα effects in inflammation and autoimmunity one might bear in mind the role of TNFα receptors, which genetic variability might underlies complex autoinflammatory diseases such as TRAPS [17]. Only a few number of variants of TNFRSF1A, encoding for the TNFR-1, were detected in our studied population and none of them was associated with the diagnosis of AS or PsA. Therefore, in the absence of copy number variation or other mutations, a role for TNFRSF1A genetics in SpA could be reasonably excluded.
Other genes involved in the inflammation processes, as genes encoding cytokines (e.g. IL23R, IL17, IL2R) or proteins of the inflammasome pathway (e.g. CARD9, CARD14, CARD16) have been reported as associated with SpA in several studies [9,39]. An interesting candidate involved in inflammasome pathway is MEFV gene, encoding for the protein pyrin, critically involved in the autoinflammatory disease FMF and recently suggested as potentially involved also in SpA [18][19][20][21][22][23][24]. Several MEFV variants were identified in our studied population, half without functional significance and half with potential functional consequences for the predicted amino acid substitutions in the protein. These variants were however sporadic and none was associated with SpA diagnosis. Intriguingly the main MEFV variant that has been demonstrated to be associated with AS, namely the M694V [21,23,24], was never recorded among our series of patients. These apparent discrepant results are probably consequent to differences in ethnicity (Italian vs Turkish) and in patient's selection criteria (SpA patients in our series vs FMF patients with AS in the Turkish series).
One of the main problems in the clinical setting of SpA is the interindividual variability in response to therapy, mainly in response to TNF-α inhibitors [40][41][42]. Different drugs belonging to this category are now available and they are frequently used to treat both AS and PsA especially when NSAIDs or DMARDS are not beneficial [43,44]. In this study about 88% AS and 59% PsA were on treatment with TNF-α inhibitors, and among them 11% AS and 19% PsA required in their disease history the switch from one to at least another drug type for an incomplete response to treatment. This empirical approach is the only one actually proposed since no predictive response biomarkers are validated. Their identification is really necessary also considering the high costs for the health care system of this type of therapy [45]. Among factors that might affect variability in response to therapy, genetics should be considered besides physiological (age, sex, weight, and body fat), pathophysiological (liver, kidney, and cardiovascular function, and associated diseases), and environmental factors (tobacco and alcohol consumption and concomitant treatments) [46][47][48]. In our series carriers of rare c.625 +10 G allele of TNFRSF1A gene were mainly late responders. This polymorphism might affect the response to TNF-α inhibitors because it is associated with an altered TNF-α/TNFR1 balance. In fact, patients carrying the rare G allele have been suggested to synthesize a truncated TNFR-1 [49], and the common A allele was previously demonstrated by us to be correlated with increased serum TNF-α levels in patients with celiac disease [30].
Although a limitation of this study is represented by the lack of complete genetic data for all enrolled patients and controls, any group of data were obtained by a significant number of patients and controls. Future perspectives are confirmatory studies on genetics as predictor of the response to TNFα inhibitors. In this setting, besides TNFRSF1A, many other genes of the TNF/TNFR pathway and of related pathways should be considered and comprehensively investigated possibly through new technologies, like next generation sequencing, in large series of patients to support clinical translation of results.
In conclusion the results of this study indicate the relevant role of TNF-TNFR pathway genetics in the complex network causing SpA and conditioning response to therapy. TNFA was shown to be a predisposing factor for SpA, but mainly for AS, among Italian patients, while genetics of the autoinflammatory gene MEFV appears of no impact in this setting. The haplotype resulting from TNFA-1031C/-308G, associated with lower TNFα production, exerts a protective role in AS, while the TNFRSF1A c.625+10A>G polymorphism emerged as a potential predictor of response to TNF-α inhibitors.  Table. MEFV gene missense polymorphisms in controls and in cases (AS and PsA patients). Results from the exploratory study. Table reports the genotypes and Minor Allele Frequency (MAF) resulting from the ten missense polymorphisms, identified in patients (AS and PsA) and controls (explanatory group), through the direct sequencing of exons 2, 3, 5 and 10 of MEFV gene. (DOC) S5 Table. R202Q polymorphism in the entire studied population. Table reports