# Genetic components of human pain sensitivity: a protocol for a genome-wide association study of experimental pain in healthy volunteers

Annina Schmid, Kaustubh Adhikari, Luis Miguel Ramirez-Aristeguieta, Juan-Camilo Chacón-Duque, Giovanni Poletti, Carla Gallo, Francisco Rothhammer, Gabriel Bedoya, Andres Ruiz-Linares, David Bennett

► To cite this version:

Annina Schmid, Kaustubh Adhikari, Luis Miguel Ramirez-Aristeguieta, Juan-Camilo Chacón-Duque, Giovanni Poletti, et al. Genetic components of human pain sensitivity: a protocol for a genome-wide association study of experimental pain in healthy volunteers. BMJ Open, 2019, 9 (4), pp.e025530. 10.1136/bmjopen-2018-025530. [https://doi.org/10.1136/bmjopen-2018-025530](https://doi.org/10.1136/bmjopen-2018-025530). [hal-03268942](https://HAL.inhal-03268942)

HAL Id: hal-03268942

https://hal.science/hal-03268942

Submitted on 16 Nov 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# BMJ Open Genetic components of human pain sensitivity: a protocol for a genome-wide association study of experimental pain in healthy volunteers

Annina B Schmid,¹ Kaustubh Adhikari,²,³,⁴ Luis Miguel Ramirez-Aristeguieta,⁵ Juan-Camilo Chacón-Duque,²,⁶ Giovanni Poletti,⁷ Carla Gallo,⁷ Francisco Rothhammer,⁸ Gabriel Bedoya,⁹ Andres Ruiz-Linares,¹⁰,¹¹ David L Bennett¹

**To cite:** Schmid AB, Adhikari K, Ramirez-Aristeguieta LM, et al. Genetic components of human pain sensitivity: a protocol for a genome-wide association study of experimental pain in healthy volunteers. *BMJ Open* 2019;9:e025530. doi:10.1136/bmjopen-2018-025530

► Prepublication history and additional material for this paper are available online. To view these files, please visit the journal online (http://dx.doi.org/10.1136/bmjopen-2018-025530).

ABS and KA contributed equally.

Received 24 July 2018
Revised 18 January 2019
Accepted 5 March 2019

## ABSTRACT

**Introduction** Pain constitutes a major component of the global burden of diseases. Recent studies suggest a strong genetic contribution to pain susceptibility and severity. Whereas most of the available evidence relies on candidate gene association or linkage studies, research on the genetic basis of pain sensitivity using genome-wide association studies (GWAS) is still in its infancy. This protocol describes a proposed GWAS on genetic contributions to baseline pain sensitivity and nociceptive sensitisation in a sample of unrelated healthy individuals of mixed Latin American ancestry.

**Methods and analysis** A GWAS on genetic contributions to pain sensitivity in the naïve state and following nociceptive sensitisation will be conducted in unrelated healthy individuals of mixed ancestry. Mechanical and thermal pain sensitivity will be evaluated with a battery of quantitative sensory tests evaluating pain thresholds. In addition, variation in mechanical and thermal sensitisation following topical application of mustard oil to the skin will be evaluated.

**Ethics and dissemination** This study received ethical approval from the University College London research ethics committee (3352/001) and from the bioethics committee of the Odontology Faculty at the University of Antioquia (CONCEPTO 01–2013). Findings will be disseminated to commissioners, clinicians and service users via papers and presentations at international conferences.

## Strengths and limitations of this study

► We propose a genome-wide association study of both baseline pain sensitivity as well as nociceptive sensitisation.

► The study will be conducted in an admixed population in Colombia, with variable proportions of European, Native American and sub-Saharan African ancestry.

► Phenotypic data will be collected by a single trained examiner, thus producing high-quality measures.

► Because we are focusing on highly reproducible pain phenotypes obtained at a single centre, our initial sample size is limited to 1500–2000 participants. This will allow us to identify only genetic variants that have intermediate and large (but not small) effect sizes.

factors in a range of chronic pain conditions.² Importantly, a growing number of studies in patient populations suggest that genetics is an important contributory factor to pain susceptibility and severity.²⁻⁴ Of note, twin studies using either clinical pain outcomes or experimental pain models suggest that sensitivity to pain has a heritability of up to 55%.⁵⁻¹⁰ Interestingly, the heritability varies greatly depending on the clinical pain outcome or the sensory modality tested in experimental pain models.⁷ There is evidence that different sensory modalities may have distinct genetic components contributing to their variance in humans,¹⁰ and this is consistent with animal models which underline the distinct neurobiology mediating different sensory modalities.¹¹

A powerful technology for identifying genetic determinants of human complex phenotypes is genome-wide association studies (GWAS). However, GWAS analyses

## INTRODUCTION

Pain constitutes a major component of the global burden of diseases, with lower back and neck pain representing the single leading cause for years lived with disability followed closely by migraine and other musculoskeletal disorders.¹ Pain is a multidimensional experience involving a highly complex interaction of physical, biochemical, physiological, cognitive, emotional, behavioural and sociocultural factors. Many studies have identified genetic

<figure><img src="image_2.png" /></figure>

Check for updates

© Author(s) (or their employer(s)) 2019. Re-use permitted under CC BY. Published by BMJ.

For numbered affiliations see end of article.

**Correspondence to**
Professor David L Bennett;
david.bennett@ndcn.ox.ac.uk
of pain remain limited, mainly due to the high numbers of individuals required to enable adequate power, and the complexities of accurately phenotyping traits that ultimately represent a personal subjective perception.² In contrast to disease cohorts where pain variation may be influenced by the severity of the disease process or its treatment, experimental pain studies measuring baseline pain sensitivity have the advantage of studying one stimulus in a standard condition (eg, controlling intensity, location and stimulus duration). To date, several candidate gene studies have been performed to determine genetic influences on human pain sensitivity using experimental pain models,⁵⁻¹² ¹³ however, to our knowledge, a full GWAS has not yet been reported. One genome-wide study evaluated association for single nucleotide polymorphisms (SNPs) revealed by exome sequencing in a subset of twins cohort that were identified as having particularly high or low heat pain sensitivity. Using pathway analysis, there was significant enrichment for variants in genes of the angiotensin pathway.¹⁴ Whereas a direct link between experimental pain sensitivity and clinical pain severity is often not present,¹⁵ there is some evidence that findings from experimental pain models can be predictive of clinically relevant pathological pain such as postoperative pain.¹⁶ Irrespective of the association between experimental and pathological pain, understanding the genetic influences on experimental pain sensitivity will provide important biological insights into the mechanisms underlying pain sensitivity.

An important drawback of most genetic studies in the pain field so far is that they have been performed mostly in populations of European ancestry, thus, they have explored only a small fraction of human phenotypic and genetic diversity.¹⁷¹⁸ This is important in the study of pain as recent studies point to variable pain thresholds for European Americans, African Americans and Latinos.¹⁹⁻²¹ However, it is not clear if this variation in pain thresholds relates to differences in neurobiological mechanisms or other factors such as social or cultural²² parameters. In addition, recent studies suggest that a phenotype associated with increased sensitisation of the nociceptive system due to temporal summation may render people vulnerable to developing clinically relevant pain.²³ However, most studies focus on pain phenotyping in the naïve state rather than sensitised state; only a handful of studies have investigated a genetic component underlying nociceptive sensitisation.⁷ ²⁴ We, therefore, intend to use an algogen (allyl isothiocyanate [AITC], an agonist of the ligand gated ion channel TRPA1) in order to sensitise the nociceptive system to replicate the changes that could occur in pathological pain states.

This protocol proposes a GWAS on genetic contributions to baseline pain sensitivity and nociceptive sensitisation in a sample of unrelated healthy individuals of mixed European/Native American/African ancestry. We will evaluate baseline cutaneous pain thresholds as well as the variation in sensitisation following mustard oil (AITC) application, a controlled model of tissue

injury. We hypothesise that we will identify SNPs associated with experimental pain stimuli in the naïve and sensitised state. Elucidating the genetic basis of pain variation has the potential to reveal targets for future analgesic development. This can be translated into improved pain management potentially tailored to an individual's pain risk or resilience, including sensitivity differences between different human populations.

## METHODS AND ANALYSIS

This GWAS follows the Strengthening the Reporting of Genetic Association Studies guidelines.²⁵ A flow chart of the study procedure is detailed in figure 1.

### Participants

Healthy participants aged 18–40 will be recruited in Medellín, Colombia via public noticeboards at local Universities, distribution of flyers and through the local print media. In addition, we are inviting previous participants from the Consortium for the Analysis of the Diversity and Evolution of Latin America (CANDELA)²⁶ GWAS to participate in this project.

Recruiting healthy young participants has advantages in GWAS studies. They are less likely to have undetected illnesses or other problems that may influence their biological pathway of pain sensitivity. Young people will also have less overall accumulated exposure or risks from environmental (external) factors which may affect their pain sensitivity. Such factors increase the overall variability of participants’ pain perception response and reduce the power of detecting genetic causes. Since most traits are affected by a combination of genetic and environmental factors, many studies including CANDELA tend to use young participants for genetic variant discovery.²⁷

Participants will be excluded if they have chronic pain or any chronic medical condition (eg, diabetes, neurodegenerative, musculoskeletal or psychiatric conditions). Participants currently taking analgesics, anti-inflammatoryies, opioids, antihistamines, antidepressants or antiepileptic drugs will be excluded. Women who are pregnant or in their menstrual phase (self-report) will be excluded from the study. Participants will be advised to not smoke or consume coffee within 1 hour of testing, and to avoid psychoactive substances or alcohol within 8 hours prior to testing. Further exclusion criteria include current or past self-inflicted injuries, as well as dermatomal, traumatic or infectious conditions affecting the arm, and a history of severe allergic reactions to any kind of medication, materials, food or insect bites. Participants with moderate to severe anxiety (≥25 on the Hamilton Anxiety Rating Scale²⁸) or severe depression (>15 on the 16-item Quick Inventory of Depressive Symptomatology Self-Reported (QIDS-SR16)²⁹ will be excluded from the study. Recruitment started in January 2013 and is predicted to take approximately 5–7 years.

### Procedure

Participants will attend a single appointment at the quantitative sensory testing (QST) laboratory at the Universidad
<figure><img src="image_3.png" /><figcaption>**Figure 1** Study procedure. The graph details the study procedure from recruitment to data analysis. Yellow boxes represent new procedures; green boxes indicate data generation and collection; blue boxes indicate a procedural step. AITC, allyl isothiocyanate; BMI, body mass index; CANDELA, Consortium for the Analysis of the Diversity and Evolution of Latin America; CDT, cold detection thresholds; GWAS, genome-wide association study; HPT, heat pain threshold; MPT, mechanical pain threshold; PCs, principal components; PPT, pressure pain threshold; QST, quantitative sensory testing; TSL, thermal sensory limen; VDT, vibration detection threshold; VAS, Visual Analogue Scale; WDT, warm detection thresholds; WUR, wind-up ratio.</figcaption></figure>
de Antioquia, Medellín. Following informed consent, age and self-reported gender will be recorded and participants will answer questions regarding their self-reported ancestry (see online supplementary appendix 1). Height and weight will be measured and body mass index (BMI) calculated. Since psychological factors such as anxiety can influence pain perception during experimental pain testing,³⁰ participants will complete the Spanish version of the Hamilton Anxiety Rating Scale and the QIDS-SR16. The QIDS-SR16 has acceptable internal consistency and moderate to strong concurrent validity compared with other depression scores³¹ and its Spanish version shows adequate test-retest reliability and high internal consistency.³² The Hamilton Anxiety Rating Scale has shown to have high inter-rater and test-retest reliability³³ and good construct validity.³⁴

## Evaluation of sensory function in the naïve state

We will determine sensory function in the naïve state and following nociceptive sensitisation. Baseline sensory function will be evaluated using specific static and dynamic QST. These include cold detection threshold and warm detection thresholds, thermal sensory limen and heat pain thresholds (HPT) using a ThermoTester (Q-sense, Medoc, Israel, 30×30 mm thermode size). Recording of thermal thresholds will strictly follow published QST guidelines.³⁵ Mechanical pain thresholds (MPT) will be evaluated using a 20 piece von Frey hair set (Touch Test, North Coast, USA) which exerts differing forces (9.8, 13.7, 19.6, 39.2, 58.8, 78.5, 98.1, 147.1, 255.0, 588.4, 980.7, 1765.2, 2942.0 mN). The von Frey hairs will be applied at a rate of 2s on, 2s off in ascending order starting from 9.8mN baseline stimulus until participants first perceive the stimulus as sharp (pricking). Subsequently, the hairs will be applied in descending order until the stimulus is perceived as blunt. The geometric mean of five series of ascending and descending stimuli is defined as the MPT. Wind-up ratio will be determined with numerical pain ratings on a Visual Analogue Scale (VAS 0–100) for a single stimulus followed by the average pain rating for a train of 10 stimuli applied at 1 Hz within the same 1 cm² using a 255 mN von Frey hair. This will be repeated five times and the ratio will be established as the mean rating of the trains of stimuli divided by the mean rating of the single stimuli. Vibration detection thresholds (VDT) will be determined by recording the mean of 3 disappearance thresholds with a Rydel-Seiffer tuning fork. Pressure pain thresholds (PPT) will be recorded in triplicate with a manual algometer (Wagner Instruments, Greenwich, Connecticut, USA) and their mean used for analysis.

The side to be tested will be randomised and patients will first be familiarised with the sensory tests on the forearm on the control side, before performing the actual measurements on the test arm. All tests will be performed halfway over the volar side of the forearm except for VDT (ulnar styloid) and PPT (thenar muscles).

<figure><img src="image_4.png" /><figcaption>**Figure 2** Method to determine area of flare, punctuate hyperalgesia and allodynia. (A) An acetate template is used to mark a star with eight spokes containing eight points at 1 cm increments on the volar forearm. (B) A small cotton swab soaked in 30% mustard oil is applied in the center of the star and (C) held in place with a Tegaderm for 10 min. During this time, pain scores are recorded every 30 s. (D and E) After removal of the mustard oil, the skin flare will be marked and the area calculated. (F) The area of brush-evoked and punctuate hypersensitivity will be determined with a brush and a 98.1 mN von Frey hair respectively (pictured) by testing potential hypersensitivity at each point on the eight spokes.</figcaption></figure>

## Mustard oil evoked nociceptive sensitisation

After the baseline sensory measures, an acetate template will be used to mark a star with eight spokes each containing eight points at 1 cm increments on the volar forearm (figure 2). The skin temperature will be standardised by placing the 32°C warm thermode over the center of the star for 5 min before starting. We will then apply a sensitisation paradigm using mustard oil (AITC (Sigma), diluted at 30% in olive oil) as previously performed.³⁶ AITC, the active component of mustard oil, activates the ion channel TRPA1 and evokes skin flare and nociceptive sensitisation.³⁷ A small cotton swab soaked in mustard oil will be applied to a 0.64 cm² area on the volar forearm and held in place with a Tegaderm (3M) for
<figure><img src="image_5.png" /></figure>

10 min.³⁶ During this time, pain scores will be recorded every 30 s using an electronic VAS ranging from 0 to 100. After 10 min, the mustard oil will be removed and the area of the skin flare will be recorded to the nearest 0.5 cm at each spoke.⁷ Eight triangular shapes will be created by joining the points on adjacent spokes and the total area will be calculated by adding all triangular segments. The area of mustard oil application will be subtracted from the total area to determine the area of secondary flare (flare area).

After mapping of the flare, the area of brush-evoked hypersensitivity will be determined with a brush (Nr 5 Senselab, Somedic, Sweden) by applying 1 cm long strokes at each of the points on the eight spokes, starting from the outside and moving towards the sensitised centre. The area of punctuate hypersensitivity will be determined with a 98.1 mN filament (Bailey Instruments, UK) following the same procedure.³⁶ As for the flare, the primary area of mustard oil application will be subtracted from both hypersensitive areas such that the recorded areas represent secondary hyperalgesia/hypersensitivity.

Following mustard oil sensitisation, the MPT and HPT will be repeated using the same methods as described above. All postsensitisation tests will be performed within 5 min of mustard oil removal.

## Reliability of naïve and sensitised sensory function protocol

To determine intratester reliability, we repeated the sensory function protocol performed by the same investigator in n=12 healthy volunteers on two different occasions within 2–6 weeks. Intraclass correlation coefficients (3.1) revealed good to excellent agreement for all sensory testing variables (table 1).³⁸

<table><caption>Table 1 Intratester reliability of sensory function protocol</caption><thead><tr><th></th><th>Intraclass<br>correlation<br>coefficients</th><th>95% CI</th><th>P value</th></tr></thead><tbody><tr><td>CDT</td><td>0.728</td><td>0.277 to 0.914</td><td>0.003</td></tr><tr><td>WDT</td><td>0.764</td><td>0.351 to 0.927</td><td>&lt;0.0001</td></tr><tr><td>TSL</td><td>0.638</td><td>0.161 to 0.878</td><td>0.005</td></tr><tr><td>HPT</td><td>0.752</td><td>0.339 to 0.922</td><td>0.002</td></tr><tr><td>MPT</td><td>0.928</td><td>0.767 to 0.979</td><td>&lt;0.0001</td></tr><tr><td>WUR</td><td>0.634</td><td>0.113 to 0.880</td><td>0.012</td></tr><tr><td>VDT</td><td>0.956</td><td>0.860 to 0.987</td><td>&lt;0.0001</td></tr><tr><td>PPT</td><td>0.734</td><td>0.305 to 0.915</td><td>0.002</td></tr><tr><td>VAS</td><td>0.893</td><td>0.667 to 0.970</td><td>&lt;0.0001</td></tr><tr><td>Flare area</td><td>0.610</td><td>0.095 to 0.869</td><td>0.015</td></tr><tr><td>Brush-evoked allodynia</td><td>0.756</td><td>0.365 to 0.922</td><td>0.001</td></tr><tr><td>Punctuate hyperalgesia</td><td>0.615</td><td>0.094 to 0.871</td><td>0.013</td></tr><tr><td>Postsensitisation MPT</td><td>0.941</td><td>0.808 to 0.983</td><td>&lt;0.0001</td></tr><tr><td>Postsensitisation HPT</td><td>0.758</td><td>0.339 to 0.924</td><td>0.002</td></tr></tbody></table>

CDT, cold detection threshold; HPT, heat pain threshold; MPT, mechanical pain threshold; PPT, pressure pain threshold; TSL, thermal sensory limen; VDT, vibration detection threshold; VAS, Visual Analogue Scale; WDT, warm detection threshold; WUR, wind-up ratio.

## Genotyping

Each participant will donate blood or saliva (Oragene OG-500, Genotek, Canada) for DNA extraction. DNA samples will be genotyped on the Illumina HumanOmniExpress chip containing ~700 000 markers. In volunteers who already participated in the CANDELA GWAS,³⁹ genotype data from blood samples genotyped on the same chip are already available and will be reused.

Whole-genome genotype data from the Illumina array will undergo quality control⁴⁰ to exclude any markers or samples that fail stringent thresholds. Quality metrics provided by the genotype calling algorithm in the Illumina GenomeStudio software,⁴¹ such as the GenTrain score, cluster separation score, and excess heterozygosity rates will be used to filter poorly genotyped SNPs. Subsequent SNP-level and sample-level quality control thresholds such as missingness will be applied. Sex mismatch between records and genetic data of X and Y chromosomes will be checked. Only samples and SNPs that pass all criteria will be retained for analysis. Details of the currently used quality control protocol for CANDELA genotyped samples are provided in online supplementary appendix 2.

## Statistical analysis

### Sample size calculation

The power for GWAS of experimental pain phenotypes for varying sample and effect sizes was estimated following the formulae described in Visscher *et al.*⁴² Estimated power is shown for a range of effect sizes for experimental pain phenotypes taken from existing experimental pain studies. The statistical software R V.3.4.1⁴³ was used to perform the calculations and produce the figures. The codes are published on https://github.com/kaustubhad/gwas-power.

In whole-genome SNP-based GWAS studies, the association analysis is usually conducted with a multivariate linear regression model, where the trait values are regressed onto an SNP genotype (with additive coding) and other covariates which commonly include age, gender, BMI and genetic principal components (PCs). The p value threshold³⁹ ⁴² for genome-wide significant associations is commonly 5×10⁻⁸, while the threshold for a suggestive significant association is commonly 10⁻⁵. Formulae to calculate power in GWAS with genome-wide and suggestive significance thresholds are presented in online supplementary appendix 3, and power calculated for the current GWAS setting is shown in figure 3.

Figure 3A shows estimated power (in percentage) as a heatmap under the standard GWAS settings of using whole-genome genotyping data and a p value significance threshold of 5×10⁻⁸. Sample size (n) varies from 100 to 5000, while the proportion of trait variance explained by the marker (q², in percentage) varies from 0.01% to 6%. As sample size increases, power increases quickly for a range of trait variance values to reach 100%.

Figure 3B shows estimated power (in percentage) under the same settings but a suggestive p value significance
<figure><img src="image_6.png" /><figcaption>**Figure 3** Estimated power (in percentage) under the standard genome-wide association studies (GWAS) settings of using whole-genome genotyping data. (A) Estimated power (in percentage) as a heatmap, setting the significance threshold at $5×10^{-8}$, the commonly used threshold for genome-wide significance in GWAS studies. (B) Estimated power with the significance threshold set at $10^{-5}$, the commonly used threshold for suggestive significance. In panels A-B, the x-axis denotes a range of sample sizes (n) in a GWAS, the y-axis represents the proportion of trait variance ($q^2$) explained by a marker. Power of detecting the marker at a specific (n, $q^2$) combination is represented by a colour gradient. Contour lines for power at 10% intervals are also shown. Panels C-D shows power curves for the expected sample sizes for this study. (C) Expected power at genome-wide and suggestive significance thresholds for a sample size of n=1500. (D) Estimated power for a sample size of n=2000. In Panels C-D, the x-axis denotes the proportion of trait variance ($q^2$) explained by a marker, and y-axis represents estimated power (in percentage). The two curves correspond to the two commonly used GWAS thresholds. In each panel, the point for 80% power is indicated with a green triangle, so that the necessary parameter configurations can be read from the graph. In panels A-B, the contour corresponding to 80% power is also marked in green.</figcaption></figure>

threshold of $10^{-5}$. As expected, power is higher at similar sample and effect sizes for this less stringent threshold.

Simplified power estimates are shown as power curves in figure 3C,D for the expected sample sizes for this study. Figure 3C shows expected power at genome-wide and suggestive significance thresholds for a sample size

of n=1500 at varying effect sizes, while figure 3D shows estimated power for a sample size of n=2000.

The range of trait variance has been taken from Doehring *et al*,44 which provides estimates for the proportion of trait variance explained by an SNP for several experimental pain phenotypes and multiple markers. The
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 4** (A) Distributions of trait variance explained by a single marker from Doehring *et al*⁴⁴ for traits included in our study and those not included. (B) Allele frequency distributions of loci associated with experimental pain in previously published cohorts, for Europeans and Colombians.</figcaption></figure>

values ranged from 0.02% to 6%. Some of the traits were the same as the traits investigated here, while some other traits were different. Nevertheless, the distributions of trait variance for the two groups of traits are very similar, as seen in **figure 4A**.

Power of a GWAS depends on the allele frequency of the SNPs through their effect on the test statistic. While the majority of GWAS studies are conducted in European-origin individuals, including the experimental pain study used to determine sample size here,⁴⁴ our population of interest is an admixed Latin American population. Therefore, we wanted to assess the distribution of allele frequencies in Europeans versus Colombians for SNPs studied for or associated with experimental pain in various studies.⁴⁴ ⁴⁵ Minor allele frequencies were obtained for all such reported SNPs from the 1000 Genomes project database¹⁷ for Western Europeans (from Britain (GBR), Utah residents from Northern and Western Europe (CEU), Spain (IBS), and Tuscany in Italy (TSI)) and Colombians (from Medellin in Colombia, (CLM), where this study will be performed). Allele frequency distributions for both Europeans and Colombians are shown in **figure 4B**. The two distributions are quite similar, with the Colombian distribution slightly more spread out. This is somewhat expected as the Colombians have on average 60% European (Spanish) ancestry.²⁶ Having a well spread out distribution of allele frequencies is important in a GWAS as low-frequency alleles have lower power for a given sample and effect size.⁴⁶ Here, the comparison to European allele frequency distribution suggests that the current Colombian cohort will have nearly equivalent power to any European-based cohort. In contrast to European-only cohorts though, our cohort will have the advantage that alleles present in other continental populations such as sub-Saharan Africans or Native Americans that are not present in Europeans would also be detectable in a GWAS, and could be followed up in replication cohorts of specific ethnicities.

The CANDELA project includes genotpyes of ~2000 patients from Medellin. We anticipate to contact and phenotype 50%–75% of these participants as well as

contacting an additional 500 participants to bring the initial sample to 1500–2000.

## Data analysis plan

The cleaned genetic data will first be merged with reference samples worldwide, such as the 1000 Genomes Project,¹⁷ Simons Genome Diversity Project,⁴⁷ Estonian Biocentre Human Genome Diversity Panel,⁴⁸ and additional European and Native American samples that are particularly relevant for Latin American populations.⁴⁹ The merged dataset will be checked for genetic outliers, through genetic PCs and continental ancestry proportions (using supervised Admixture⁵⁰), and for unexpected genetic similarities. These steps can often detect any sample misplacement or contamination, which might be reflected in sex mismatch, unexpected genetic similarities or inflated heterozygosity rates.⁴⁰ ⁵¹ Genetic ancestry estimates will be compared with self-reported ancestry information (see online supplementary appendix 1), particularly for genetic outliers or samples showing unexpected results. Participants self-reporting for ethnicities rare in Colombians, such as East Asian or South Asian, would also be excluded as outliers. The authors have extensive experience in conducting association analysis in admixed populations, including several GWAS publications on a wide range of phenotypes which contain detailed protocols on how to conduct such analyses.²⁷ ³⁹ ⁵² ⁵³ Further details of the currently used quality control protocol for CANDELA samples⁵³ are provided in online supplementary appendix 2.

The genetic data allow estimation of the narrow-sense heritability of any quantitative trait, which is the fraction of trait variance that is explained by the genetic data. Estimates of heritability, obtained using the software GCTA,⁵⁴ will provide an idea of which traits have more of a biological basis versus which are more environmentally determined, and thus which traits would be more amenable to genetic analysis for discovery of associated genetic variants. Note however that relatively precise estimate (low SDs) of heritability by this method requires several thousand samples,⁵² so the currently proposed
sample size might be underpowered to estimate heritability accurately.

To facilitate better identification of associated loci, the genotype data will be imputed to approximately 10 million loci using the 1000 Genomes phase 3 imputation reference panel<sup>52</sup> by first haplotype phasing using SHAPEIT2<sup>55</sup> and then imputation using Impute2.<sup>56</sup> Quality control of the imputed genotypes will be performed using recommended thresholds on imputation quality score, concordance metrics and proportion of high-probability calls. Details of the currently used imputation protocol for CANDELA samples<sup>52</sup> are provided in online supplementary appendix 2.

GWAS studies will be conducted in Plink2<sup>57</sup> to perform single-locus association studies for each trait individually, across the whole genome in an additive multivariate linear regression model.<sup>42</sup> Covariates will be used in the regression to adjust for any other sources of trait variability, such as basic variables like age, sex, and BMI, and genetic PCs will be used to control for population substructure.<sup>52 58</sup>

The number of genetic PCs to be included in the regression depends on the sample composition, such as variation in ancestry and presence/absence of genetic outliers. It would be determined by inspecting the proportion of variance explained by each PC (displayed on a scree plot) and by checking PC scatter plots.

In addition to being used as exclusion criteria, anxiety and depression scores could be used as covariates in GWAS. The exact set of covariates to be used will be determined based on initial diagnostic analyses such as correlation analysis.

These single-locus association results, obtained as p values, will be visualised via the Manhattan plot. Commonly used p value thresholds for selecting associated loci are 5×10<sup>-8</sup> for genome-wide significance and 10<sup>-5</sup> for suggestive significance.<sup>42</sup>

An extension of this additive multivariate linear regression model, still within the single-trait single-locus setting, called the mixed linear model analysis which better controls for any cryptic relatedness or population substructure, will also be performed in GCTA.<sup>54</sup>

There are several extensions of the single-trait single-locus association studies that increase power for detecting associated loci: combining several related traits that may share a biological basis, using multivariate Wald tests as implemented in MultiPhen<sup>59</sup>; or gene-based tests that combine signals across all loci in a gene to increase signal strength and reduce the burden of multiple testing, such as set-based models implemented in Plink2<sup>57</sup> or fastBAT implemented in GCTA.<sup>60</sup> The admixed nature of the sample might be used in detecting associations by the method of admixture mapping,<sup>61</sup> though the potential of success of this method in detecting associated variants depends on the extent of stratification of the variant's allele frequency across continents. These analyses might help detect additional loci that are underpowered in classical GWAS due to smaller effect sizes.

## Handling of missing data

It might not be possible to record some traits in some individuals, even though the completeness of the first 100 samples suggests that missingness will be low. The single-trait methods used in traditional GWAS analyses automatically exclude individuals from the analysis of a trait who have missing values for that trait. The same applies to individuals having missing genotypes for any particular SNP. However, genotyping success rate using the Illumina HumanOmniExpress chip in the CANDELA cohort is very high (>99.8%), so the number of excluded individuals in any analysis would be very low overall.

Some multivariate analyses such as PCs when applied on the set of phenotypes require having recorded values of all phenotypes for an individual. Instead of using the subset of individuals who have the complete set of phenotypes recorded, which would incur some loss in sample size, the missing phenotype data for each individual will be imputed following standard statistical procedures as implemented in the R package 'mice'.<sup>62</sup> When the proportion of missing data is small, imputation is preferable in such multivariate analyses than sample exclusion, and is routinely applied to genetic data such as while calculating genetic PCs.<sup>58</sup>

## DISCUSSION

This GWAS including a well-defined cohort of healthy participants will provide important insights into the genetic aspects underlying experimental pain sensitivity in the naïve and sensitised state. This may allow further exploration of potential biological mechanisms underlying pain sensitivity. Future studies will be required to extrapolate these findings to patient populations with chronic pain.

## Patient and public involvement

No patient involvement is performed during this study.

## Ethics and dissemination

Findings will be disseminated to commissioners, clinicians and service users via papers and presentations at international conferences such as the biennial World Congress of International Association for the Study of Pain. We will also post our findings to the publicly available painnetworks.org database.

## Author affiliations

<footnote n="1"><sup>1</sup>Nuffield Department of Clinical Neurosciences, Oxford University, Oxford, UK</footnote>

<footnote n="2"><sup>2</sup>Department of Genetics, Evolution and Environment, University College London, London, UK</footnote>

<footnote n="3"><sup>3</sup>Department of Cell and Developmental Biology, University College London, London, UK</footnote>

<footnote n="4"><sup>4</sup>School of Mathematics and Statistics, Faculty of Science, Technology, Engineering and Mathematics, The Open University, Milton Keynes, UK</footnote>

<footnote n="5"><sup>5</sup>GST Lab. Faculty of Odontology, Universidad de Antioquia, Medellin, Colombia</footnote>

<footnote n="6"><sup>6</sup>Department of Earth Sciences, Natural History Museum, London, UK</footnote>

<footnote n="7"><sup>7</sup>Unidad de Neurobiología Molecular y Genética, Laboratorios de Investigación y Desarrollo, Facultad de Ciencias y Filosofía, Universidad Peruana Cayetano Heredia, Lima, Peru</footnote>

<footnote n="8"><sup>8</sup>Instituto de Alta Investigación, Universidad de Tarapacá, Arica, Chile</footnote>
⁹GENMOL (Genética Molecular), Universidad de Antioquia, Medellín, Colombia

¹⁰Ministry of Education Key Laboratory of Contemporary Anthropology and Collaborative Innovation Center of Genetics and Development, School of Life Sciences and Human Phenome Institute, Fudan University, Shanghai, Shanghai, China

¹¹CNRS, EFS, ADES, Aix-Marseille Université, Marseille, France

**Acknowledgements** We thank Peter Kamerman for his input to parts of the study design, and Macarena Fuentes-Guajardo for helping with the phenotyping photographs in figure 2. We would like to thank the CANDELA team for their contributions in organisation, recruitment, data collection, sample processing and genotyping: Alvaro Alvarado, Ana Carolina Orozco, Caio C. Silva de Cerqueira, Claudia Jaramillo, David Balding, Diana Rogel Diaz, Francisco de Ávila Becerril, Francisco M. Salzano, Francisco Quispealaya, Hugo Villamil-Ramírez, Ilich Jafet Moreno, Javier Mendoza-Revilla, Javier Rosique, Jorge Gómez-Valdés, Joyce De la Piedra, Lavinia Schuler-Faccini, Mónica Ballesteros Romero, Major Eugênio Correa de Souza Junior, Malena Hurtado, María Pizarro, María Teresa Del Solar, Marcelo Zagonel de Oliveira, Mari-Wyn Burley, María-Cátira Bortolini, Martha Granados Riveros, Miguel Ángel Contreras Sieck, Paola Everardo, Paola Everardo Martínez, Paola León-Mimila, Ricardo Cebrecos, Ricardo Gunski, Rodrigo Barquera Lozano, Rolando Gonzalez-José, Rosilene Paim, Ruth Rojas, Samuel Canizales-Quinteros, Sergeant João Felisberto Menezes Cavalheiro, Tábita Hunemeier, Valeria Villegas, Vanessa Granja, Vanessa Sarabia, Victor Acuña-Alonzo, Virginia Ramallo, Wendy Hart, William Arias and William Flores.

**Contributors** DLB and AR-L conceptualised the study. Funding was acquired by DLB, ABS, GB and AR-L. ABS, J-CC-D, KA and LMR-A have contributed to the study methodology. AR-L, KA, J-CC-D, GB, FR, CG and GP contributed to the design and execution of the CANDELA aspect of this study. The first draft of the protocol was prepared by ABS, DLB, KA and AR-L and all authors provided critical evaluation and approved the final manuscript.

**Funding** Work leading to this publication was funded by grants from: the Newton Fund Institutional Links Grant from the British Research Council (grant number 216398412), the Excellence Initiative of Aix-Marseille University—A*MIDEX (a French 'Investissements d'Avenir' programme) and from a Wellcome Trust Senior Scientist Fellowship to DLB (ref. no. 095698z/11/z and 202747/Z/16/Z). KA was supported by a Wellcome Investigator Award WT107055AIA (to C.D. Stern). ABS was supported by the National Institute for Health Research (NIHR) Biomedical Research Centre (BRC) Oxford.

**Disclaimer** The views expressed are those of the author(s) and not necessarily those of the NHS, the NIHR or the Department of Health.

**Competing interests** None declared.

**Patient consent for publication** Not required.

**Ethics approval** This study received ethical approval from the University College London research ethics committee (3352/001) and from the bioethics committee of the Odontology Faculty at the University of Antioquia (CONCEPTO 01-2013).

**Provenance and peer review** Not commissioned; externally peer reviewed.

**Open access** This is an open access article distributed in accordance with the Creative Commons Attribution 4.0 Unported (CC BY 4.0) license, which permits others to copy, redistribute, remix, transform and build upon this work for any purpose, provided the original work is properly cited, a link to the licence is given, and indication of whether changes were made. See: https://creativecommons.org/licenses/by/4.0/.

## REFERENCES

1. GBD 2015 Disease and Injury Incidence and Prevalence Collaborators. Global, regional, and national incidence, prevalence, and years lived with disability for 310 diseases and injuries, 1990-2015: a systematic analysis for the Global Burden of Disease Study 2015. *Lancet* 2016;388:1545-602.

2. Zorina-Lichtenwalter K, Meloto CB, Khoury S, et al. Genetic predictors of human chronic pain conditions. *Neuroscience* 2016;338:36-62.

3. Bennett DL, Woods CG. Painful and painless channelopathies. *Lancet Neurol* 2014;13:587-99.

4. Denk F, McMahon SB. Neurobiological basis for pain vulnerability: why me?. *Pain* 2017;158(Suppl 1):S108-14.

5. Belfer I, Diatchenko L. *Pain Genetics: basic to translational science*. Wiley: Oxford, 2014.

6. Nielsen CS, Stubhaug A, Price DD, et al. Individual differences in pain sensitivity: genetic and environmental contributions. *Pain* 2008;136:21-9.

7. Norbury TA, MacGregor AJ, Urwin J, et al. Heritability of responses to painful stimuli in women: a classical twin study. *Brain* 2007;130(Pt 11):3041-9.

8. Hocking LJ, Morris AD, Dominiczak AF, et al. Heritability of chronic pain in 2195 extended families. *Eur J Pain* 2012;16:1053-63.

9. Angst MS, Phillips NG, Drover DR, et al. Pain sensitivity and opioid analgesia: a pharmacogenomic twin study. *Pain* 2012;153:1397-409.

10. Nielsen CS, Knudsen GP, Steingrimsdóttir ÓA. Twin studies of pain. *Clin Genet* 2012;82:331-40.

11. Mogil JS. Pain genetics: past, present and future. *Trends Genet* 2012;28:258-66.

12. Kim H, Mittal DP, Iadarola MJ, et al. Genetic predictors for acute experimental cold and heat pain sensitivity in humans. *J Med Genet* 2006;43:e40.

13. Kim H, Neubert JK, San Miguel A, et al. Genetic influence on variability in human acute experimental pain sensitivity associated with gender, ethnicity and psychological temperament. *Pain* 2004;109:488-96.

14. Williams FM, Scollen S, Cao D, et al. Genes contributing to pain sensitivity in the normal population: an exome sequencing study. *PLoS Genet* 2012;8:e1003095.

15. Hübscher M, Moloney N, Leaver A, et al. Relationship between quantitative sensory testing and pain or disability in people with spinal pain-a systematic review and meta-analysis. *Pain* 2013;154:1497-504.

16. Werner MU, Mjöbo HN, Nielsen PR, et al. Prediction of postoperative pain: a systematic review of predictive experimental pain studies. *Anesthesiology* 2010;112:1494-502.

17. Auton A, Brooks LD, Durbin RM, et al. A global reference for human genetic variation. *Nature* 2015;526:68-74.

18. Popejoy AB, Fullerton SM. Genomics is failing on diversity. *Nature* 2016;538:161-4.

19. Rahim-Williams FB, Riley JL, Herrera D, et al. Ethnic identity predicts experimental pain sensitivity in African Americans and Hispanics. *Pain* 2007;129:177-84.

20. Rahim-Williams B, Riley JL, Williams AK, et al. A quantitative review of ethnic group differences in experimental pain response: do biology, psychology, and culture matter? *Pain Med* 2012;13:522-40.

21. Bagwath Persad LA, Kamerman PR, Wadley AL. Predictors of Cold and Pressure Pain Tolerance in Healthy South African Adults. *Pain Med* 2017;18:pnw291-37.

22. Zborowski M. Cultural components in responses to pain.1. *J Soc Issues* 1952;8:16-30.

23. Weissman-Fogel I, Granovsky Y, Crispel Y, et al. Enhanced presurgical pain temporal summation response predicts post-thoracotomy pain intensity during the acute postoperative phase. *J Pain* 2009;10:628-36.

24. Tegeder I, Adolph J, Schmidt H, et al. Reduced hyperalgesia in homozygous carriers of a GTP cyclohydrolase 1 haplotype. *Eur J Pain* 2008;12:1069-77.

25. Little J, Higgins JP, Ioannidis JP, et al. STrengthening the REporting of Genetic Association studies (STREGA)--an extension of the STROBE statement. *Eur J Clin Invest* 2009;39:247-66.

26. Ruiz-Linares A, Adhikari K, Acuña-Alonzo V, et al. Admixture in Latin America: geographic structure, phenotypic diversity and self-perception of ancestry based on 7,342 individuals. *PLoS Genet* 2014;10:e1004572.

27. Adhikari K, Fontanil T, Cal S, et al. A genome-wide association scan in admixed Latin Americans identifies loci influencing facial and scalp hair features. *Nat Commun* 2016;7:10815.

28. Hamilton M. The assessment of anxiety states by rating. *Br J Med Psychol* 1959;32:50-5.

29. Rush AJ, Trivedi MH, Ibrahim HM, et al. The 16-Item Quick Inventory of Depressive Symptomatology (QIDS), clinician rating (QIDS-C), and self-report (QIDS-SR): a psychometric evaluation in patients with chronic major depression. *Biol Psychiatry* 2003;54:573-83.

30. Ignatiadis M, Polyzos A, Stathopoulos GP, et al. A multicenter phase II study of docetaxel in combination with gefitinib in gemcitabine-pretreated patients with advanced/metastatic pancreatic cancer. *Oncology* 2006;71(3-4):159-63.

31. Reilly TJ, MacGillivray SA, Reid IC, et al. Psychometric properties of the 16-item Quick Inventory of Depressive Symptomatology: a systematic review and meta-analysis. *J Psychiatr Res* 2015;60:132-40.

32. Gili M, Lopez-Navarro E, Homar C, et al. Psychometric properties of Spanish version of QIDS-SR16 in depressive patients. *Actas Esp Psiquiatr* 2014;42:292-9.
33. Shear MK, Vander Bilt J, Rucci P, et al. Reliability and validity of a structured interview guide for the Hamilton Anxiety Rating Scale (SIGH-A). *Depress Anxiety* 2001;13:166-78.

34. Clark DB, Donovan JE. Reliability and validity of the Hamilton Anxiety Rating Scale in an adolescent sample. *J Am Acad Child Adolesc Psychiatry* 1994;33:354-60.

35. Rolke R, Magerl W, Campbell KA, et al. Quantitative sensory testing: a comprehensive protocol for clinical trials. *Eur J Pain* 2006;10:77-88.

36. Kremeyer B, Lopera F, Cox JJ, et al. A gain-of-function mutation in TRPA1 causes familial episodic pain syndrome. *Neuron* 2010;66:671-80.

37. Koltzenburg M, Lundberg LE, Torebjörk HE. Dynamic and static components of mechanical hyperalgesia in human hairy skin. *Pain* 1992;51:207-19.

38. Shrout PE, Fleiss JL. Intraclass correlations: uses in assessing rater reliability. *Psychol Bull* 1979;86:420-8.

39. Adhikari K, Reales G, Smith AJ, et al. A genome-wide association study identifies multiple loci for variation in human ear morphology. *Nat Commun* 2015;6:7500.

40. Anderson CA, Pettersson FH, Clarke GM, et al. Data quality control in genetic case-control association studies. *Nat Protoc* 2010;5:1564-73.

41. Guo Y, He J, Zhao S, et al. Illumina human exome genotyping array clustering and quality control. *Nat Protoc* 2014;9:2643-62.

42. Visscher PM, Wray NR, Zhang Q, et al. 10 Years of GWAS Discovery: Biology, Function, and Translation. *Am J Hum Genet* 2017;101:5-22.

43. TTeam RC. R: A language and environment for statistical computing. Vienna, Austria. 2018 https://www.R-project.org/2018

44. Doehring A, Küsener N, Flühr K, et al. Effect sizes in experimental pain produced by gender, genetic variants and sensitization procedures. *PLoS One* 2011;6:e17724.

45. Duan G, Guo S, Zhang Y, et al. The Effect of SCN9A Variation on Basal Pain Sensitivity in the General Population: An Experimental Study in Young Women. *J Pain* 2015;16:971-80.

46. Nm L, Lange C. *The Fundamentals of Modern Statistical Genetics*: Springer, 2011.

47. Mallick S, Li H, Lipson M, et al. The Simons Genome Diversity Project: 300 genomes from 142 diverse populations. *Nature* 2016;538:201-6.

48. Pagani L, Lawson DJ, Jagoda E, et al. Genomic analyses inform on migration events during the peopling of Eurasia. *Nature* 2016;538:238-42.

49. Reich D, Patterson N, Campbell D, et al. Reconstructing Native American population history. *Nature* 2012;488:370-4.

50. Alexander DH, Novembre J, Lange K. Fast model-based estimation of ancestry in unrelated individuals. *Genome Res* 2009;19:1655-64.

51. Chacón-Duque JC, Adhikari K, Fuentes-Guajardo M, et al. Latin Americans show wide-spread Converso ancestry and imprint of local Native ancestry on physical appearance. *Nat Commun* 2018;9:5388.

52. Adhikari K, Fuentes-Guajardo M, Quinto-Sánchez M, et al. A genome-wide association scan implicates DCHS2, RUNX2, GLI3, PAX1 and EDAR in human facial variation. *Nat Commun* 2016;7:11616.

53. Adhikari K, Mendoza-Revilla J, Sohail A, et al. A GWAS in Latin Americans highlights the convergent evolution of lighter skin pigmentation in Eurasia. *Nat Commun* 2019;10:358.

54. Yang J, Lee SH, Goddard ME, et al. GCTA: a tool for genome-wide complex trait analysis. *Am J Hum Genet* 2011;88:76-82.

55. O'Connell J, Gurdasani D, Delaneau O, et al. A general approach for haplotype phasing across the full spectrum of relatedness. *PLoS Genet* 2014;10:e1004234.

56. Howie B, Fuchsberger C, Stephens M, et al. Fast and accurate genotype imputation in genome-wide association studies through pre-phasing. *Nat Genet* 2012;44:955-9.

57. Chang CC, Chow CC, Tellier LC, et al. Second-generation PLINK: rising to the challenge of larger and richer datasets. *Gigascience* 2015;4:7.

58. Price AL, Patterson NJ, Plenge RM, et al. Principal components analysis corrects for stratification in genome-wide association studies. *Nat Genet* 2006;38:904-9.

59. O'Reilly PF, Hoggart CJ, Pomyen Y, et al. MultiPhen: joint model of multiple phenotypes can increase discovery in GWAS. *PLoS One* 2012;7:e34861.

60. Bakshi A, Zhu Z, Vinkhuyzen AA, et al. Fast set-based association analysis using summary data from GWAS identifies novel gene loci for human complex traits. *Sci Rep* 2016;6:32894.

61. Seldin MF. Admixture mapping as a tool in gene discovery. *Curr Opin Genet Dev* 2007;17:177-81.

62. Buuren S, Groothuis-Oudshoorn K. mice: Multivariate Imputation by Chained Equations in R. *J Stat Softw* 2011;45.