Genome-wide association study meta-analysis of chronic widespread pain: evidence for involvement of the 5p15.2 region

Background and objectives Chronic widespread pain (CWP) is a common disorder affecting ∼10% of the general population and has an estimated heritability of 48–52%. In the first large-scale genome-wide association study (GWAS) meta-analysis, we aimed to identify common genetic variants associated with CWP. Methods We conducted a GWAS meta-analysis in 1308 female CWP cases and 5791 controls of European descent, and replicated the effects of the genetic variants with suggestive evidence for association in 1480 CWP cases and 7989 controls. Subsequently, we studied gene expression levels of the nearest genes in two chronic inflammatory pain mouse models, and examined 92 genetic variants previously described associated with pain. Results The minor C-allele of rs13361160 on chromosome 5p15.2, located upstream of chaperonin-containing-TCP1-complex-5 gene (CCT5) and downstream of FAM173B, was found to be associated with a 30% higher risk of CWP (minor allele frequency=43%; OR=1.30, 95% CI 1.19 to 1.42, p=1.2×10−8). Combined with the replication, we observed a slightly attenuated OR of 1.17 (95% CI 1.10 to 1.24, p=4.7×10−7) with moderate heterogeneity (I2=28.4%). However, in a sensitivity analysis that only allowed studies with joint-specific pain, the combined association was genome-wide significant (OR=1.23, 95% CI 1.14 to 1.32, p=3.4×10−8, I2=0%). Expression levels of Cct5 and Fam173b in mice with inflammatory pain were higher in the lumbar spinal cord, not in the lumbar dorsal root ganglions, compared to mice without pain. None of the 92 genetic variants previously described were significantly associated with pain (p>7.7×10−4). Conclusions We identified a common genetic variant on chromosome 5p15.2 associated with joint-specific CWP in humans. This work suggests that CCT5 and FAM173B are promising targets in the regulation of pain.


INTRODUCTION
Chronic widespread pain (CWP) is a common disorder, affecting about 10% of the general population. 1 The prevalence of CWP increases with age for both men and women, but is more common in women at any age. 1 CWP represents a major underestimated health problem and is associated with substantial impairment and a reduced quality of life. It has been related to a number of physical and affective symptoms such as fatigue, psychological distress and somatic symptoms. 1 2 Chronic musculoskeletal pain is one of the most common conditions seen in rheumatology clinics and accounts for 6.2% of the total healthcare costs in The Netherlands every year. 3 Further research is needed to be able to understand the causal mechanisms and optimal treatment for CWP patients.
CWP causally relates to an initial local pain stimulus, such as an acute injury or athletic injuries, or another pain state such as low back pain or local pain due to osteoarthritis (OA) or rheumatic arthritis (RA). [4][5][6] However, most injured subjects do not develop CWP, and only a proportion of patients with OA or RA develop CWP. We therefore hypothesise that several discrete stimuli may initiate CWP via a common final pathway that involves the generation of a central pain state through the sensitisation of second order spinal neurons.
CWP is a complex trait since both environmental and genetic factors play a role in the aetiology. Heritability estimates of twin studies suggest that 48-52% of the variance in CWP occurrence is due to genetic factors, implying a strong genetic component. 7 A number of studies have examined genetic variants for CWP. These candidate gene studies examined polymorphisms in genes involved in both the peripheral and the central nervous system. 8 In particular, genes involved in neurotransmission ( pathway of dopamine and serotonin [9][10][11][12][13][14][15][16][17][18][19], and genes important for the hypothalamic-pituitaryadrenal axis have been considered. 20 A number of genetic variants in these candidate genes were found to be associated with CWP, individual pain sites or experimental pain. However, no consistent significant associations have been demonstrated. The most studied gene in relation to pain is catechol-O-methyltransferase (COMT), an enzyme that degrades neurotransmitters including dopamine. The variant allele of rs4680 (or V158M) results in reduced enzymatic activity due to its effect on thermostability, 21 and has been associated with reduced opioid activity in response to painful stimuli resulting in increased pain sensitivity. 22 But also for COMT, no consistent results have been observed in genetic association studies. 13 [23][24][25][26][27][28][29] Overall, the results have been conflicting, which is likely due to the modest sample sizes used and paucity of replication. In general, candidate studies are biased by previous knowledge of the aetiology of the disease under study. Since knowledge about the pathophysiology of CWP is poor, the chances of success using this approach are low. Therefore our objective was to identify genetic variants involved in CWP by means of a large-scale hypothesis-free genome-wide association study (GWAS) meta-analysis including 2788 cases and 13 780 controls. To our knowledge, this is the first study presenting a large-scale GWAS meta-analysis of chronic pain. The prevalence of CWP is approximately two times higher in women than in men and there is strong evidence that women tolerate less thermal and pressure pain than men. 30 Therefore only women were included in this study to reduce heterogeneity and thereby increase power.

MATERIALS AND METHODS
We performed a meta-analysis (stage 1) of GWAS data of 1308 female Caucasian CWP cases and 5791 female Caucasian controls, derived from five studies, and focused our follow-up efforts on the single-nucleotide polymorphisms (SNPs) with suggestive evidence of association ( p<1×10 −5 ) with CWP (stage 2). The study outline is summarised in figure 1.

Phenotype
CWP was defined as subjects having pain in the left side of the body, in the right side of the body, above the waist, below the waist, and in the axial skeleton (following the Fibromyalgia Criteria of the American College of Rheumatology 2 ). Controls were defined as subjects not having CWP. Subjects using analgesics (ATC code: N02 31 ) were excluded from the control group. Detailed descriptions of the study specific inclusion criteria are presented in supplementary table S1.

Study design summary
We combined the summary statistics of GWAS in a meta-analysis comprising 1308 CWP female Caucasian cases and 5791 female Caucasian controls (stage 1). We focused our follow-up efforts on the SNPs with suggestive evidence of association ( p<1×10 −5 ) with CWP in 1480 CWP cases and 7989 controls available for replication (stage 2).

Subjects
A full detailed description of all study cohorts is presented in table 1 and in the supplementary methods section. For the stage 1 analysis, we included studies from The Netherlands (the Erasmus Rucphen Family study (ERF study), 32 Rotterdam Study I, II and  III (RS-I, RS-II and RS-III) 33 ), and the UK (TwinsUK 34 35 ). All studies were approved by their institutional ethics review committees and all participants provided written informed consent. For our stage 2 analysis, we sought follow-up samples with pre-existing GWAS in silico data (stage 2a) as well as de novo genotyping (stage 2b). The studies are from the UK (the British 1958 Birth Cohort (1958BC), 23 36-38 the Chingford Study (CHINGFORD), 39 46 47 ). All Genotypes of the stage 2a studies (1958BC, AGES, DSDBAC, FOA, GARP and SHIP) were obtained from SNP arrays and imputed data. Where unavailable, proxy SNPs were selected based on high linkage disequilibrium (LD). The stage 2b studies (CHINGFORD, EPIFUND and HCS) performed de novo genotyping, using both Sequenom iPLEX and TaqMan-based assays (supplementary methods). Genotyping platforms, calling algorithms, quality control before imputation, imputation methods and analysis software used were all study-specific (see supplementary tables S4 and S5). The explicit number of follow-up SNPs genotyped in the different studies and whether the original or a proxy SNP was used is summarised in supplementary table S6.

GWAS analysis in the stage 1 studies
CWP was analysed as a binary trait (cases vs controls) using logistic regression under an additive model with adjustment for age and body mass index (see supplementary table S7). To adjust for population substructure, we included the four most important PCs as covariates in the regression analysis of RS-I, RS-II and RS-III. These PCs were derived from a multidimensional scaling analysis of identity-by-state distances, using PLINK software. 48 Detailed descriptions of the GWAS methods are provided in supplementary table S8).

Stage 1: GWAS meta-analysis
p Values for association were combined using the Meta-Analysis Tool for genome-wide association scans (METAL). 49 The genomic control method 50 as implemented in METAL was used to correct for any residual population stratification or relatedness not accounted for by the four most important PCs. A p value <5×10 −8 was considered genome-wide significant while a p value <1×10 −5 was considered suggestive. 51 Power calculations were

SNP selection for replication
We aimed to select SNPs for replication (stage 2) that were enriched for signals of association with CWP. All SNPs with suggestive evidence for association in the stage 1 analyses were selected and separated into independent loci by taking the most significantly associated SNP and eliminating all SNPs that have a HapMap CEU pairwise correlation coefficient r 2 >0.8 with that SNP using the PLINK software.

Meta-analysis of stage 1 and stage 2 results
We combined the stage 1 and stage 2 association results to derive a combined meta-analysis for the suggestively associated loci. METAL was used to conduct a fixed-effects meta-analysis as in stage 1. Estimated heterogeneity variance and forest plots were generated using comprehensive meta-analysis (http:// www.meta-analysis.com).

Functional analysis of associated SNPs
To determine whether the associated SNPs have any regulatory effect on gene expression levels, we checked their effect (and the effect of the linked SNPs) on the expression levels of their neighbouring genes. We used the 1000 genomes data in the SNAP software 52 53 to identify those SNPs having LD thresholds of r 2 >0.1. We searched two publicly available eQTL databases: the NCBI GTEx (Genotype-Tissue Expression) eQTL browser (http://www.ncbi.nlm.nih.gov/gtex/GTEX2/gtex.cgi) and the expression Quantitative Trait Loci database (http:// eqtl.uchicago.edu/cgi-bin/gbrowse/eqtl/). We used SIFT 54 to predict whether the coding non-synonymous variant causing an amino acid substitution affects protein function.

RNA expression analyses in mice
For functional follow-up, two independent mouse models of inflammatory pain were studied. The first model was based on carrageenan injections; female C57Bl/6 mice received an intraplantar injection of 20 μl λ-carrageenan (2% (w/v), Sigma Aldrich, Zwijndrecht, the Netherlands) in saline in both hind paws. 55 The second model was based on Complete Freund's Adjuvant (CFA) injections; male C57Bl/6 mice (Harlan Laboratories) received an intraplantar injection of 20 μl CFA (Sigma-Aldrich) in saline in both hind paws. 56 Controls were injected with saline only. At day 3 (after CFA injection) or day 6 (after carrageenan injection), thermal sensitivity (heat withdrawal latency time) was measured using the Hargreaves (IITC Life Science, Woodland Hills, California, USA) test as described. 57 Intensity of the light beam was chosen to induce heat withdrawal latency time of approximately 8 s at baseline. After measurement the mice were sacrificed and the lumbar (L2-L5) spinal cord and the dorsal root ganglions (DRG) (L2-L5) were isolated. These areas of spinal cord and DRG were selected because pain transmission from the hind paws is mediated via primary sensory neurons that have their cell bodies in the lumbar DRG, and transmit the signal to the lumbar spinal cord through sensory fibres in the dorsal roots. Total RNA was isolated and mRNA levels of Cct5 and Fam173b were measured in the spinal cord and the DRG. For more details, see the supplementary methods section.
All experiments were performed in accordance with international guidelines and approved by the experimental animal committee of the University Medical Center Utrecht (carrageenan experiment) or the UK Home Office Animals (Scientific Procedures) Act 1986 (CFA experiment). Mice used for the carrageenan experiment were bred and maintained in the animal facility of the University of Utrecht (The Netherlands).

Systemic review of genetic variants previously described
We systematically searched for associations earlier reported with pain in the HugeNavigator PhenoPedia database. 58 We used the search term 'pain' and checked all publications for genes and SNPs associated with pain at least twice. Genes and SNPs associated with drug therapy, facial pain, migraine and postoperative pain were excluded. For all reported SNPs, we examined their association with CWP in our stage 1 meta-analysis. The significance threshold was set at p<8×10 −4 using Bonferroni correction for 65 independent genetic loci. Again, power calculations were performed using CaTS software (http://www.sph.umich.edu/csg/ abecasis/CaTS/). With an α level of 8×10 −4 , power calculations showed that we had approximately 80% power to detect an OR of 1.22 for SNPs with a minor allele frequency of 20% or higher.

Meta-analysis of GWAS replication
For the 10 independent SNPs with suggestive evidence, we pursued in silico replication data in six studies (stage 2a: 1203 CWP cases and 5032 controls) and performed de novo genotyping in subjects from three additional studies (stage 2b: 277 CWP cases and 2957 controls) (a detailed description of the studies is presented in table 1 and supplementary methods). The summary results of the stage 1 and 2 meta-analysis are presented in table 2. After combining the results of stage 1 and stage 2, the top SNP was rs13361160 (OR=1.17, 95% CI 1.10 to 1.24, p=4.7×10 −7 , I 2 =28.4%). Figure 3 shows a forest plot of the association of rs13361160 with CWP across the stage 1 and stage 2 studies. The overall effect in the replication studies (stage 2 studies) was in a consistent direction but not significant (OR=1.06, 95% CI 0.98 to 1.16, p=0.16). In the combined analysis, moderate heterogeneity was observed (I 2 =28.4%). Supplementary table S1 shows the different pain assessment methods used in the different studies to define CWP. Since four out of five stage 1 studies included joint-specific pain only (ERF, RS-I, RS-II and RS-III), we performed a sensitivity analysis in which stage 2 cohorts using non-joint pain were excluded (1958BC, DSDBAC, EPIFUND, HCS and SHIP). This resulted in a combined OR of 1.23 (95% CI 1.14 to 1.32, p=3.4×10 −8 , I 2 = 0%). An overview of the results of the combined meta-analysis and the separate stage 1 and stage 2 analyses is presented in table 3.

Functional analysis of rs13361160 and rs2386592
The SNPs rs13361160 and rs2386592 (r 2 =0.97) are annotated to the 5p15.2-region and located 81 kb upstream of CCT5 and 57 kb downstream of FAM173B (figure 4). We tested whether rs13361160 and rs2386592 and their linked SNPs (r 2 >0.1) affected gene expression levels of CCT5 or FAM173B. In total, we identified 130 SNPs in LD with our top SNPs, of which two SNPs were located in the coding region: one synonymous SNP rs1042392 in the CCT5 gene (r 2 =0.16, D 0 =0.85) and one non-synonymous SNP rs2438652 in the FAM173B gene (r 2 = 0.17, D 0 =1.0) (see supplementary table S9). The minor allele of rs2438652 causes a threonine-to-methionine substitution (T75M) which is thought to be functionally neutral. SNPs rs13361160 and rs2386592 were not recorded as influencing the expression levels of CCT5 and FAM173B, however the linked intronic SNP rs2445871 (r 2 =0.14 for both) had a direct eQTL effect on FAM173B expression levels in liver tissue. 59

RNA expression analysis in mice
We studied gene expression levels of the two nearest genes, Cct5 and Fam173b, in the lumbar spinal cord and the DRG in two independent mouse models of chronic inflammatory pain. In both the carrageenan treated group and the CFA treated group, mice had shorter heat withdrawal latency times than mice injected with saline only, confirming enhanced pain sensitivity ( p<0.001) (see supplementary figure S1).
The results from the multivariate analysis using the two genes (Cct5 and Fam173b examined as dependent variables), the different treatments (saline, carrageenan and CFA) and the different tissues (DRG and spinal cord) confirmed that there is a significant treatment effect for Cct5 (F(2,25) 5). These findings indicate that in spinal cord but not in DRG, both Fam173b and Cct5 expression levels were up-regulated in response to two different inducers of inflammatory pain. DRG Fam173b and Cct5 expression levels in CFA/carrageenan-treated mice were indistinguishable from saline-treated mice.

Candidate SNPs previously associated with chronic pain
We examined whether genetic variants previously described for association with pain were associated with CWP in our large stage 1 meta-analysis. We identified a total of 44 genes, of which 136 SNPs had been reported at least twice with any pain phenotype (excluding facial pain, migraine, postoperative pain and response to drug therapy), and we examined the association of these 136 SNPs with CWP in the GWAS stage 1 meta-analysis. Out of 136 candidate SNPs, we were able to check 92 common SNPs (MAF>5%) in 65 independent genetic loci (see supplementary table S10). Five SNPs had a too low MAF (<=5%) and 39 SNPs were not genotyped or imputed in our meta-analysis. None of the earlier reported SNPs passed the significance threshold ( p<8×10 −4 ). Interestingly, the strongest associated SNPs are located in three genes that have been reported to be associated with pain phenotypes most frequently: COMT, GCH1 (GTP cyclo-hydrolase 1) and OPRM1 (mu opioid receptor). The effects of the SNPs in GCH1 are in the same direction as reported earlier [60][61][62] : individuals having the minor allele for rs10483639, rs4411417 or rs752688 have 15% less pain than those exhibiting the common alleles. The effect of the SNP rs599548 in OPRM1 is also in the same direction as reported earlier 63 : those having the minor allele for rs599548 have 19% more pain than those exhibiting the major allele. The two COMT SNPs are in weak LD with the wellknown amino acid changing variant rs4860, but previously have not been reported to be significantly associated with pain. 23 64 We have found a protective effect for the minor allele of rs2020917 (those having a minor allele have 15% less pain) and an adverse effect for the minor allele of rs5993883: those having the minor allele of rs5993883 have 14% more pain.

DISCUSSION
In this study, we identified a genetic variant near CCT5 and FAM173B to be associated with CWP. Chronic pain coincided with higher RNA expression of Cct5 and Fam173b in the lumbar spinal cord of mouse models of inflammatory pain. This finding indicates that both genes in the 5p15.2 region are regulated in the context of inflammatory pain.
Interestingly, Bouhouche et al 65 reported a human pedigree in which a CCT5 mutation caused hereditary sensory neuropathy (Online Mendelian Inheritance in Man (OMIM) ID=610150), a syndrome characterised by a sensory deficit in the distal portion of the lower extremities, chronic perforating ulcerations of the feet and progressive destruction of underlying bones. Symptoms can include pain and numbness, tingling in the hands, legs or feet, and extreme sensitivity to touch. CCT5 is a subunit of the chaperonin containing t-complex polypeptide 1 (TCP-1) which assists in protein folding and assembly in the brain. 66 CCT5 interacts with the serine/threonine-protein phosphatase 4 catalytic subunit PP4C. [67][68][69] Zhang et al 70 confirmed that protein phosphates like PPP4C may have a regulatory effect on the central sensitisation of nociceptive transmission in the spinal cord. Interestingly, sensitisation is thought to contribute to chronic inflammatory pain. 71 Since the function of the FAM173B gene is not yet known, it is difficult to postulate the mechanism by which this gene could influence CWP. Further research into the genes in this locus is needed to ascertain whether either or both CCT5 and FAM173B are driving the observed association.
By combining the effects across the different stage 2 studies, moderate heterogeneity was observed in the meta-analysis. This heterogeneity might be caused by different pain assessment methods used by the stage 2 cohorts. In particular, four cohorts asked the participants about joint pain specifically, while the other five also included non-joint pain. When the non-joint pain phenotype were excluded, the heterogeneity across the cohorts reduced to 0% and the overall p value for rs13361160 now reached genome-wide significance by combining the stage 1 and stage 2 effects. This might suggest that indeed phenotype heterogeneity was introduced by including non-joint pain. In general, it is anticipated that pain is a very complex trait, with different aetiological pathways introducing phenotypic heterogeneity. A limitation of our study is that we were not able to examine possible phenotype subgroups, such as individuals with RA, a chronic systemic inflammatory disorder that principally affects the synovial joints. Stratifying these groups of individuals might serve to increase power to find genetic loci. We here decided to analyse all CWP cases together, based on the hypothesis that several discrete stimuli need to initiate CWP via a common final pathway that involves the generation of a central pain state through the sensitisation of second order spinal neurons. In addition, the prevalence of RA is very low (about 0.5-1%), 72 and the earlier defined GWAS hits for RA (ie, the HLA locus) 73 were not in our top list. So, we assume the results were not dominated by this small number of individuals with RA.
It would be helpful to dissect the phenotype of pain into quantitative sub-phenotypes, for example by measuring pain sensitivity and pain thresholds for temperature or pressure, 74 or by examining functional MRIs. 75 The use of quantitative and possibly more objective pain measurements in response to painful stimuli (rather than reported pain) will be of pivotal importance for future pain research. Because we have focused on the clinical pain definition using questionnaires and pain homunculus, we accept that we may have missed true pain susceptibility alleles. However, this study represents the largest genome-wide meta-analysis looking into the genetics of human CWP to date. The experiments in two independent mouse models of chronic inflammatory pain showed that the expression of Cct5 and Fam173b was higher in the lumbar spinal cord of mice with chronic inflammatory pain but not in DRG. In the spinal cord, the expression profiles of both genes were up-regulated in response to two different inducers of inflammatory pain. These findings indicate that both genes in the 5p15.2 region are co-regulated in the spinal cord during inflammation-induced pain in both independent pain models, thereby possibly contributing to the neurobiology of pain. In the lumbar DRG, containing the cell bodies of the primary sensory neurons that detect pain signals from the hind paws, Cct5 and Fam173b gene expression levels did not change by inflammation. Because of these complementary results from the two independent tissues (spinal cord and DRG), we hypothesise that the 5p15.2 region is likely to play a role in spinal central pain processing and not in regulating primary sensory neuron responses.  In both analyses the effect estimates of the models refer to the minor allele (=effect allele). BMI, body mass index; MAF, minor allele frequency.
In the study of candidate genes previously reported to be associated with a pain phenotype, we showed that none of the 92 studied variants were significantly associated with CWP in our GWAS meta-analysis. This can be explained by the fact that many of the previous reported loci were studied in relative modest sample sizes and in a large variety of pain phenotypes. 76 Power calculations show that we had approximately 80% power to detect an OR as low as 1.22 for SNPs with an allele frequency of 20% or higher. So, even in this large meta-analysis, power was still modest to detect small ORs and we therefore cannot exclude smaller effect sizes of the tested variants, resulting in lack of reproducibility. 77 This lack of reproducibility of SNPs in candidate genes in large GWAS meta-analyses has been shown before for other phenotypes such as bone mineral density (BMD). 78 It is interesting to note that among the candidate SNPs, the strongest associated ones were located in the three most studied pain genes, COMT, GCH1 and OPRM1. The directions of the effects of these SNPs were the same as reported earlier, which would support true associations.
In conclusion, our study reports a GWAS meta-analysis on CWP. We identified the genetic variant rs13361160 at the 5p15.2 locus, located 81 kb upstream of the CCT5 gene and 57 kb downstream of the FAM173B gene, to be associated with CWP. We showed an increase in expression levels of Cct5 and Fam173b in the spinal cord of inflammatory pain models of  mice, and since these genes both seem to influence the central mechanism of sensitisation, they may represent a novel pathway involved in pain sensation.