Abstract
Background
Bilateral sagittal split ramus osteotomy (BSSRO) is a common orthognatic surgical procedure. Sensory disturbances in the inferior alveolar nerve, including hypoesthesia and dysesthesia, are frequently observed after BSSRO, even without distinct nerve injury. The mechanisms that underlie individual differences in the vulnerability to sensory disturbances have not yet been elucidated.
Methods
The present study investigated the relationships between genetic polymorphisms and the vulnerability to sensory disturbances after BSSRO in a genome-wide association study (GWAS). A total of 304 and 303 patients who underwent BSSRO were included in the analyses of hypoesthesia and dysesthesia, respectively. Hypoesthesia was evaluated using the tactile test 1 week after surgery. Dysesthesia was evaluated by interview 4 weeks after surgery. Whole-genome genotyping was conducted using Illumina BeadChips including approximately 300,000 polymorphism markers.
Results
Hypoesthesia and dysesthesia occurred in 51 (16.8%) and 149 (49.2%) subjects, respectively. Significant associations were not observed between the clinical data (i.e., age, sex, body weight, body height, loss of blood volume, migration length of bone fragments, nerve exposure, duration of anesthesia, and duration of surgery) and the frequencies of hypoesthesia and dysesthesia. Significant associations were found between hypoesthesia and the rs502281 polymorphism (recessive model: combined X2 = 24.72, nominal P = 6.633 × 10−7), between hypoesthesia and the rs2063640 polymorphism (recessive model: combined X2 = 23.07, nominal P = 1.563 × 10−6), and between dysesthesia and the nonsynonymous rs2677879 polymorphism (trend model: combined X2 = 16.56, nomina P = 4.722 × 10−5; dominant model: combined X2 = 16.31, nominal P = 5.369 × 10−5). The rs502281 and rs2063640 polymorphisms were located in the flanking region of the ARID1B and ZPLD1 genes on chromosomes 6 and 3, whose official names are “AT rich interactive domain 1B (SW11 -like)” and “zona pellucida-like domain containing 1”, respectively. The rs2677879 polymorphism is located in the METTL4 gene on chromosome 18, whose official name is “methyltransferase like 4”.
Conclusions
The GWAS of sensory disturbances after BSSRO revealed associations between genetic polymorphisms located in the flanking region of the ARID1B and ZPLD1 genes and hypoesthesia and between a nonsynonymous genetic polymorphism in the METTL4 gene and dysesthesia.
Keywords
Background
Neuropathic pain in the orofacial region is a clinical manifestation of trigeminal nerve injury following oral surgery. Neuropathic pain subsequent to nerve damage at a central or peripheral site remains a major problem for both patients and clinicians because the pain is usually extremely intense and often refractory to various conventional pain therapies. Moreover, remarkable individual differences in the vulnerability to neuropathic pain exist. Many studies have been performed to reveal the mechanisms that underlie neuropathic pain, but only a few genetic studies have focused on neuropathic pain [1,2], possibly because individual differences across patients with neuropathic pain are usually affected by various factors other than genetic factors.
Sensory disturbances, including hypoesthesia and dysesthesia, often appear as a prodromal symptom of neuropathic pain. Sensory disturbances or neuropathic pain in the inferior alveolar nerve are inevitably caused by a primary lesion or dysfunction of the nerve. The symptoms, however, are subject to individual differences in daily clinical practice and may be related to genetic factors. Bilateral sagittal split ramus osteotomy (BSSRO) is commonly conducted to correct jaw deformities, such as mandibular prognathism. Sensory disturbances in the inferior alveolar nerve, including hypoesthesia and dysesthesia, are frequently observed in the lower lip and mental area after BSSRO, even without distinct nerve injury. Symptom frequency 1 or 2 weeks after BSSRO is reported in 25–56% of patients [3–5]. Considering that almost all patients who undergo BSSRO are young and healthy and the degree of surgical invasiveness, surgical site, and surgical procedures are highly consistent across cases, environmental factors appear to have relatively little impact on individual differences in the vulnerability to sensory disturbances or neuropathic pain after BSSRO.
Innovative techniques have been used to investigate the genetic factors related to various human traits. A wide array of information on the entire human genome has accumulated, and the results of genome-wide association studies (GWASs) have been reported [6,7]. A marked increase in the rate of discovery of genes associated with various diseases has also occurred [8].
The present GWAS investigated the relationships between genetic polymorphisms and the vulnerability to sensory disturbances after BSSRO.
Results
Clinical data overview and SNP data management for GWAS
Hypoesthesia and dysesthesia occurred in 51 (16.8%) and 149 (49.2%) of the 304 and 303 patients, respectively (Table 1). Logistic regression analysis revealed no significant associations between the clinical data and frequency of hypoesthesia or dysesthesia after BSSRO (data not shown).
Expression frequency of hypoesthesia and dysesthesia after BSSRO
After filtering the markers by genotype call frequency, “Cluster sep”, and minor allele frequencies in the first quality control assessment of the genotyping data, 243,501 markers were selected. These merged genotype data from five different BeadChips consisted of single nucleotide polymorphism (SNP) markers on the autosome or sex chromosome, and no mitochondrial marker was included. Furthermore, 272 markers were excluded based on the Hardy-Weinberg equilibrium test (P ≤ 2 × 10−7). As a result, a total of 243,229 SNP markers (including 4,822 nonsynonymous SNPs) were selected for the subsequent association study (Additional file 1: Figure S1 and Additional file 2: Figure S2).
GWAS identified several loci associated with sensory disturbances in the inferior alveolar nerve after BSSRO
The GWAS was performed to detect any signals associated with hypoesthesia or dysesthesia after BSSRO as three-stage analyses for two independent patterns: (1) a normal GWAS procedure that targeted all of the SNPs that were available (Additional file 1: Figure S1) and (2) a GWAS procedure that targeted only nonsynonymous SNPs that tended to affect the function of the protein encoded by the relevant gene (Additional file 2: Figure S2).
In the first analysis that targeted all of the SNPs, six, five, and 22 SNPs were selected as the top candidates associated with hypoesthesia for the trend, dominant, and recessive models for each minor allele, respectively, after the final stage (Table 2). Seven, four, and nine SNPs were selected as the top candidates associated with dysesthesia for the trend, dominant, and recessive models for each minor allele, respectively, after the final stage (Table 3). Among these, two SNPs, rs502281 and rs2063640, showed significant associations with hypoesthesia after the final stage in the recessive model (rs502281: X2 = 16.44, Q = 0.0196; rs2063640: X2 = 14.38, Q = 0.0291; Table 2). None of the SNPs showed significant associations with dysesthesia after the final stage in any of the models (Table 3).
Top candidate SNPs selected after final stage analysis in 3-stage GWAS targeting all SNPs (hypoesthesia)
CHR, chromosome number; Position, chromosomal position (bp); Q, Q value for FDR correction of multiple comparison; Related gene, the nearest gene from the SNP site;*, Significant after FDR correction(Q < 0.05).
Top candidate SNPs selected after final stage analysis in 3-stage GWAS targeting all SNPs (dysesthesia)
CHR, chromosome number; Position, chromosomal position (bp); Q, Q value for FDR correction of multiple comparison; Related gene, the nearest gene from the SNP site;†modified from the Illumina annotation file.
In the second analysis that targeted nonsynonymous SNPs, four, three, and 14 SNPs were selected as the top candidates associated with hypoesthesia for the trend, dominant, and recessive models for each minor allele, respectively, after the second stage (Table 4). Three, five, and two SNPs were selected as the top candidates associated with dysesthesia, respectively, after the second stage (Table 5). Among these, none of the SNPs showed significant associations with hypoesthesia after the final stage in any of the models (Table 4). One SNP, rs2677879, showed significant associations with dysesthesia after the final stage in the trend and dominant models (trend model: X2 = 6.585, Q = 0.0309; dominant model: X2 = 6.669, Q = 0.0491; Table 5). Statistical power analyses revealed that the expected power (1 minus type II error probability) was only 19.5% and 15.1% for the Cohen's conventional “small” effect size of 0.10 [9] and 90.8% and 84.6% for the medium effect size of 0.30, with a total of 120 valid samples in each stage. The degrees of freedom were set at 1 and 2, respectively, for the nominal type I error probability of 0.05. The estimated effect sizes were 0.26 and 0.28 to achieve 80% power for this type I error probability using our samples. The degrees of freedom were set at 1 and 2, respectively.
Top candidate SNPs selected after second stage analysis in 3-stage GWAS targeting nonsynonymous SNPs (hypoesthesia)
CHR, chromosome number; Position, chromosomal position (bp); Q, Q value for FDR correction of multiple comparison; Related gene, the nearest gene from the SNP site; NA, data not available.
Top candidate SNPs selected after second stage analysis in 3-stage GWAS targeting nonsynonymous SNPs (dysesthesia)
CHR, chromosome number; Position, chromosomal position (bp); Q, Q value for FDR correction of multiple comparison; Related gene, the nearest gene from the SNP site; *, Significant after FDR correction (Q < 0.05).
Candidate loci revealed by the GWAS were located around/within the gene regions of ARID1B, ZPLD1, and METTL4
Figures 1 and 2 present the genome-wide associations between polymorphism markers and the susceptibility to hypoesthesia evaluated by the Semmes-Weinstein pressure aesthesiometer test after BSSRO for all of the samples in each model for each chromosome. Significant associations were found between hypoesthesia and the rs502281 SNP (recessive model: combined X2 = 24.72, nominal P = 6.633 × 10−7; Table 2; Additional file 3: Table S1) and rs2063640 SNP (recessive model: combined X2 = 23.07, nominal P = 1.563 × 10−6; Table 2; Additional file 3: Table S2) and between dysesthesia and the rs2677879 SNP (trend model: combined X2 = 16.56, nominal P = 4.722 × 10−5; dominant model: combined X2 = 16.31, nominal P = 5.369 × 10−5; Table 5; Additional file 3: Table S3) in two independent patterns of analyses with all of the samples. According to the annotation information supplied by the manufacturer of the whole-genome genotyping arrays (Illumina, San Diego, CA), the rs502281 and rs2063640 SNPs are located within the gene flanking region of ARID1B and ZPLD1 on chromosomes 6 and 3 (Table 2; Figure 1), whose official names are “AT rich interactive domain IB (SWI1-like)” and “zona pellucida-like domain containing 1”, respectively, based on the National Center for Biotechnology Information database [10]. The rs2677879 SNP is located within the gene region of METTL4 on chromosome 18 (Table 5; Figure 2), whose official name is “methyltransferase like 4”, based on the same database.
Discussion
The present study explored genome-wide associations between common genetic variations and sensory disturbances after BSSRO. There are occasional reports in the literature about the relationship between individual genetic polymorphisms and neuropathic pain [11,12]. One study investigated the association between catechol-O-methyltransferase gene polymorphisms and pain sensitivity and musculoskeletal pain attributed to temporomandibular disorders [13]. Another study focused on the association between HLA gene polymorphisms and postherpetic neuralgia, also known as intractable chronic pain disorder [14]. Although a GWAS was previously conducted in patients with neuropathic pain induced by administration of paclitaxel for breast cancer [15], no other such studies have been performed to determine the development of postoperative peripheral neuropathy. BSSRO is among the most frequent surgical procedures in the area of oral surgery, and its procedures are well standardized. Because patient candidates for BSSRO are relatively healthy and young, they are a good population for studies of postoperative peripheral neuropathy. We conducted a GWAS to investigate the onset of sensory disturbances after BSSRO.
The results of the present study showed that hypoesthesia and dysesthesia occurred in 16.8% (51 of 304) and 49.2% (149 of 303) of the patients, respectively. Our incidence rate for hypoesthesia tended to be lower than previously reported incidences that ranged from 25% to 56% [3–5]. One reason for this may be the fact that BSSRO is performed by a limited number of skilled surgeons at our hospital, although several other reasons may explain the lower incidence of hypoesthesia. Hypoesthesia and dysesthesia are classified into vulnerability of the peripheral nerve to external stress and property of emergence of neuropathic pain following nerve injury, respectively. Thereby, the candidate genes, which were found in the present study, should be associated with these two aspects.
The GWAS identified ARID1B, ZPLD1, and METTL4 as candidates that may be associated with the onset of sensory disturbances. The ARID1B gene, which is located in 6q25.3, encodes a protein that is a member of the ARID family of DNA-binding proteins and a subunit of human SWI/SNF-related complexes. The SWI/SNF complexes are known to use energy generated by an integral adenosine triphosphatase subunit to remodel chromatin. These complexes are involved in maintaining normal cellular functions and restricting the access of regulatory factors to nucleosomal DNA [16]. The ARID1B gene has been suggested to be associated with the occurrence of Coffin-Siris syndrome [17], a multiple congenital anomaly/mental retardation syndrome characterized by mild to moderate mental retardation, moderate to severe hypotonia, epilepsy, and congenital malformation, including a coarse facial appearance and incompletely formed fifth fingers and toes. Haploinsufficiency of the ARID1B gene is speculated to be a common potential cause of intellectual disability and speech impairment. The nervous system may be involved in the intractability and chronicity of neuropathic pain [18,19], but it is unclear whether ARID1B is associated with pain mechanism. According to the HapMap database [20], however, the rs2817461 and rs502281 SNPs identified in the present study are located upstream (approximately 200 kbp) from the ARID1B gene. Further studies are needed to examine the effects of these SNPs on ARID1B gene expression and function.
The functions of ZPLD1 remain unclear, but one report investigated the involvement of ZPLD1 in cerebral cavernous malformations [21]. The ZPLD1 gene may be involved in the development of cerebral cavernous malformations at the mRNA expression level. Additionally, a high incidence of epilepsy is found in patients with cerebral cavernous malformations [22], suggesting the involvement of ZPLD1 in the nervous system. ZPLD1 is also reportedly associated with childhood obesity [23]. However, it is unclear whether ZPLD1 is associated with pain mechanism. According to the HapMap database, the re2063640 SNP identified as a candidate in the present study is located in a relatively downstream region (approximately 4 kbp) that is close to the ZPLD1 gene. This SNP may exert an effect on the gene expression level of ZPLD1, but this needs to be clarified in future studies.
The METTL4 gene is located on the chromosome region 18p11.32. Detailed information on the functions of its gene product, however, is unavailable. No studies of which we are aware have reported associations between METTL4 and specific diseases. Based on the molecular structure of METTL4, it may affect methylation, which plays a major role in various epigenetic regulatory mechanisms. DNA methylation, recognized as the most common type of epigenetic modifications, is involved in gene silencing and plays an important role in gene regulation, development, and tumorigenesis. It has also been shown to be associated with the pathophysiology of various nervous and mental disorders. A mutation in MeCP2, a methyl-CpG binding protein, reportedly causes Rett syndrome, characterized by mental retardation and autism [24]. With regard to acquired mental disorders, abnormal DNA methylation is found in the brains of patients with schizophrenia and depression. Using microarray technology, Mill et al. comprehensively analyzed DNA methylation in the frontal lobe in patients with schizophrenia and bipolar (manic-depressive) disorder and found changes in the DNA methylation of genes involved in brain development and stress responses [25]. According to the dbSNP database [26], the rs2677879 SNP, a candidate identified in the present GWAS of nonsynonymous polymorphisms, leads to amino acid substitution from Gin to Lys, likely causing functional changes in the protein. Although the precise functions of METTL4 are poorly understood, a representative METTL, METTL11A, reportedly exhibited catalytic activity as a histone methyltransferase [27]. Although future studies are needed, the action of METTL4 might be involved in methyltransferase activity and thus cause the methylation of genomic DNA close to related genes, which could result in the modulation of neural transmission related to sensory disturbances.
The genes identified in the present study are different from those previously reported to be associated with neuropathic pain. Future studies that involve larger numbers of patients may identify previously reported gene polymorphisms and determine the functional relationships between the present study and peripheral neuropathy. We did not consider the patients' personalities (i.e., psychological factors) in the present study, which should be addressed in future studies.
Conclusion
The present GWAS determined the onset of sensory disturbances after BBSRO and identified three gene polymorphisms in or near the region of the ARIBD1, ZPLD1, and METTL4 genes. Elucidating the relationship between neuropathic pain and genetic factors will elucidate the risk factors for neuropathic pain in individual patients, thereby allowing the selection of tailored treatments.
Methods
Patients
Enrolled in the study were 304 healthy patients (American Society of Anesthesiologists Physical Status I; age, 15–50 years; 114 males and 190 females) who were scheduled to undergo BSSRO for mandibular prognathism at Tokyo Dental College Suidoubashi Hospital (Table 6). The study protocol was approved by the Institutional Review Board, Tokyo Dental College, Chiba, Japan, and the Institutional Review Board, Tokyo Institute of Psychiatry (currently Tokyo Metropolitan Institute of Medical Science), Tokyo, Japan. Written informed consent was obtained from all of the patients or parents when the patients were younger than 20 years old and any accompanying image. Patients who presented with distinct nerve injury during BSSRO were excluded from the study.
Clinical data
Anesthesia and surgery
Four experienced, skilled surgeons were selected. These surgeons were board-certified in the oral surgery specialty. General anesthesia was induced with target-controlled infusion (TCI) of propofol using a TCI pump (TE-371, Terumo, Tokyo, Japan). Vecuronium (0.1 mg/kg) was administered to facilitate nasotracheal intubation. After the induction of anesthesia, 10 ml of venous blood was sampled for the preparation of DNA specimens. General anesthesia was maintained with propofol at a target blood concentration of 4–6 μg/ml. Vecuronium was administered at a rate of 0.08 mg/kg/h The lungs were ventilated with oxygen-enriched air. Local anesthesia was performed on the right side of the surgical field with 8 ml of 2% lidocaine that contained 12.5 μg/ml epinephrine, and right mandibular ramus osteotomy was performed. Local anesthesia was then performed on the left side, and left mandibular ramus osteotomy was performed. The bilateral mandibular bone segments were fixed in appropriate positions (Figure 3). Whenever systolic blood pressure or heart rate exceeded +20% of the preinduction value during surgery, intravenous (i.v.) fentanyl (1 μg/kg) was administered. At the end of surgery, a rectal diclofenac sodium suppository (50 mg) and dexamethasone (8 mg, i.v.) were administered to prevent orofacial edema and postoperative pain. Oral mecobalamin (1.5 mg/day) was administered for 4 weeks after the operation.
Illustration of bilateral sagittal split ramus osteotomy, which sagittally splits the mandibular ramus into inside and outside bone fragments.
Evaluation of sensory disturbances
Sensory disturbances were determined postoperatively by the presence of hypoesthesia or dysesthesia in the mental nerve area. Hypoesthesia was evaluated by tactile-threshold tests 1 week after the operation. The 1 week time-point was chosen for assessment to avoid testing during the time when postoperative pain was severe. The tactile-threshold test was performed using a Semmes-Weinstein pressure aesthesiometer (Research Design, Houston, TX, USA; (Figure 4) [28]. The Semmes-Weinstein pressure aesthesiometer consisted of 20 filaments with different diameters. The end of each filament was mounted into an individual Lucite rod. The amount of force applied through the long axis of each filament to achieve a noticeable bend was determined. The magnitude of these forces ranged from 0.0045 g to 447 g. This test was performed by two experienced dentists.
Touch stimulation was performed using the method of Bell [29]. The Semmes-Weinstein pressure aesthesiometer was perpendicularly lowered to a test region for 1-1.5 s and then lifted for 1-1.5 s. Stimulation was applied three times with 1.65–4.08 manufacturer's filament marking and calculated force (Fmg) and once with 4.17–6.65 Fmg at each point. All of these filaments, with the exception of the largest (6.65 Fmg), bent when they reached the specified pressure. Stimulation began with the 1.65 Fmg filament (i.e., the thinnest filament), and the stimulation force was increased until the patient perceived the stimulation. Tactile sensitivity was recognized to be positive when the patient perceived any stimulation, even if the stimulation was not perceived as a normal tactile sensation.
Based on the running courses of the labial inferior ramification and mental ramification, measurements were performed at two points [3]: (1) the vermilion border at one-third the distance between the oral angles and (2) the midpoint of the perpendicular line from point (1) to the lower margin of the mentum.
The worst among the values obtained at the four total test-points on both sides was regarded as the representative value. This value was evaluated by the interpretation scale reported by Bell [29]. In this scale, sensory function is classified into five grades. In the present study, the patients who were classified into grades that were worse than the second grade (2.83 Fmg) were regarded as hypoesthesic.
A patient who spontaneously recognized any abnormal sensations was regarded as dysesthesic. The evaluation of dysesthesia was based on the definition of the International Association for the Study of Pain. Subjective symptoms were assessed by interview 4 weeks after the operation. The patients were asked to select words from the McGill Pain Questionnaires [30] to describe their pain (i.e., temporal, brightness, thermal, dullness, traction pressure, constrictive pressure, etc.). The time-point of 4 weeks was chosen for assessment to avoid testing during the time of Wallerian [31] degeneration and retrograde degeneration after nerve damage.
Whole-genome genotyping
Genomic DNA was extracted from whole-blood samples using standard procedures. The extracted DNA was dissolved in TE buffer (10 mM Tris-HCl and 1 mM ethylenediaminetetraacetic acid, pH 8.0). The DNA concentration was adjusted to 100 ng/μl using a NanoDrop ND-1000 Spectrophotometer (NanoDrop Technologies, Wilmington, DE, USA).
Whole-genome genotyping was performed using Infinium assay II utilizing an iScan system (Illumina) according to the manufacturer's instructions, with a total of 361 samples including those of the patients enrolled in the study. Genotyping was conducted basically the same way as a previous report [32]. Five kinds of BeadChips were used for genotyping 40, 67, 6, 120, and 128 samples, respectively: HumanHap300 (total markers: 317,503), HumanHap300-Duo (total markers: 318,237), Human610-Quad v1 (total markers: 620,901), Human1M v1.0 (total markers: 1,072,820), and Human 1M-Duo v3 (total markers: 1,199,187). Some BeadChips include a number of probes that are specific to copy number variation markers, but most were for SNP markers on the human autosome or sex chromosome. Approximately 300,000 SNP markers were commonly included in all of the BeadChips.
Quality control
The data for the genotyped samples were analyzed using BeadStudio or GenomeStudio with the Genotyping module v3.3.7 (Illumina) to evaluate the quality of the results. The genotype data from all five of the BeadChips were merged to analyze all of the samples simultaneously (i.e., only the markers common to all of the BeadChips were included in the analysis, and the others were automatically excluded). In the data-cleaning process, the samples with a genotype call rate of less than 0.95 were excluded from further analyses. Markers with a genotype call frequency of less than 0.95, “Cluster sep” (i.e., an index for genotype cluster separation) of less than 0.1, and minor allele frequencies of less than 0.05 were excluded from the subsequent association study.
Statistical analysis
Prior to the GWAS, associations between the clinical data and hypoesthesia or dysesthesia expression frequency after BSSRO were analyzed. Clinical data included gender, age, body weight, body height, loss of blood volume, migration length of bone fragments, duration of anesthesia, and duration of surgery (Table 6). A logistic regression analysis was performed using SPSS (12.0J for Windows, SPSS Japan, Tokyo, Japan).
The Fisher's exact test was performed for all of the genotype frequency data to investigate the deviation of the distributions from those in the theoretical Hardy-Weinberg equilibrium, which sometimes reflects genotyping errors or population stratification of the samples. Markers with P values (df = 1) greater than approximately 2 × 10−7 (0.05/300,000) were considered for the GWAS.
A multistage GWAS was conducted for the patients who underwent painful cosmetic surgery to investigate the association between genetic variations and sensory disturbances after BSSRO. Among 361 subjects, one subject did not meet the quality control criteria in our preliminary analysis, and 57 and 58 subjects lacked clinical data for hypoesthesia and dysesthesia, respectively. Therefore, genotype data for a total of 360 subjects were used for our three-stage GWAS (120 subjects for each of the first-, second-, and final-stage analyses). Clinical data for a total of 304 and 303 subjects were used for our three-stage GWAS of hypoesthesia (104, 98, and 102 subjects for the first-, second-, and final-stage analyses, respectively) and dysesthesia (105, 96, and 102 subjects for the first-, second-, and final-stage analyses, respectively), respectively. The subjects were recruited within several years and randomly categorized into three independent groups to minimize bias in the clinical data, indicating that the samples and clinical data were not used in chronological order for our first-, second-, and final-stage analyses. In our preliminary analysis that used merged markers between different BeadChips with BeadStudio or GenomeStudio, 295,036 SNPs (including 6,016 nonsynonymous SNPs) were selected for the analyses.
For the GWAS, the Cochran-Armitage trend test was performed to explore markers that might confer susceptibility to hypoesthesia evaluated by the Semmes-Wemstem pressure aesthesiometer test or dysesthesia after BSSRO. The patients were divided into two groups based on the presence or absence of symptoms, and a linear trend analysis of the increased rate of subjects with an increased number of variant risk alleles was performed for all markers. Moreover, dominant and recessive genetic models for each minor allele were used for the analyses because of the previously insufficient knowledge about the genetic factors associated with sensory disturbances after BSSRO. The association study included both female and male subjects for autosomal markers, although male genotypes were excluded from the analysis of X chromosome markers. All of the statistical analyses were performed using gPLINK v. 2.050, PLINK v. 1.07 PLINK [33], and Haploview v. 4.1 [34]. Single-nucleotide polymorphism annotations were created based on an annotation file within Human 1M-Duo v3 supplied by the manufacturer of the BeadChips. For calculation of Q-values, SFDR (Stratified False Discovery Rate) software [35] or PLINK v. 1.07 was used. Power analyses were performed using G*Power v. 3.0.5 [36].
The GWAS procedure is summarized in the Additional file 1: Figure S1 and Additional file 2: Figure S2. In the first-stage analysis of 104 and 105 subjects for hypoesthesia and dysesthesia, respectively, the SNPs that had statistical P values of less than 0.05 were selected as the candidate SNPs for the second-stage analysis among the SNP that passed the quality control criteria within the 295,036 SNPs (6,016 nonsynonymous SNPs). For these SNPs, the second-stage analysis was conducted. Again, the SNPs that had P values of less than 0.05 were considered potential candidates and selected for further final-stage analysis. Linkage disequilibrium (LD)-based SNP pruning was also conducted in this stage utilizing PLINK v. 1.07 software, and SNPs that were in approximate linkage equilibrium with an SNP were excluded based on the following process: (i) consider a window of 50 SNPs, (ii) calculate LD between each pair of SNPs in the window, (iii) remove one of a pair of SNPs if the LD is greater than 0.8, and (iv) shift the window five SNPs forward and repeat the procedure. In the final stage, the association study was conducted to determine whether the possible associations between the SNPs selected in the second stage and phenotypic traits would be strictly replicated. In this stage, the Q values of the false discovery rate were calculated to correct for multiple testing, in addition to P values based on previous reports [37,38]. The SNPs with Q < 0.05 in the analysis were considered genome-wide significant.
Two independent patterns of the GWAS were conducted to effectively explore candidate SNPs that showed statistically strong association with the phenotypic traits and those that could functionally impact neighboring genes. In the first pattern, a normal GWAS procedure targeted all of the SNPs that were available (Additional file 1: Figure S1). In the second pattern, the GWAS procedure targeted only nonsynonymous SNPs that tended to affect the function of the protein encoded by the relevant gene (Additional file 2: Figure S2).
A log quantile-quantile (QQ) P-value plot as a result of the GWAS for the combined samples was subsequently drawn to check the pattern of the generated P-value distribution, in which the observed P values against the values expected from the null hypothesis of uniform distribution, calculated as -log10 (P value), were plotted for each model. Many of the plots were mostly concordant with the expected line (y = x), especially over the range of 0 < -log10 (P value) < 4, indicating no apparent population stratification of the samples used in the study, although the plots for the recessive model, especially for hypoesthesia, apparently deviated over the range of -log10 (P value) > 3 (Additional file 4: Figure S3 and Additional file 5: Figure S4).
Additional files
Additional file 1: Figure S1
Schematic illustration of the multistage GWAS that targeted all of the SNPs that were available. Potent candidate SNPs associated with (
Additional file 2: Figure S2
Schematic illustration of the multistage GWAS that targeted only nonsynonymous SNPs. Potent candidate SNPs associated with (
Additional file 3 Table S1
Frequencies of hypoesthesia in patients with the ARID1B (rs502281) genotype.
Additional file 4: Figure S3
Log quantile-quantile (QQ) P value plot for all of the samples as a result of the GWAS for hypoesthesia evaluated by the Semmes-Weinstein pressure aesthesiometer test after BSSRO in (
Additional file 5: Figure S4
Log quantile-quantile (QQ) P value plot for all of the samples as a result of the GWAS for dysesthesia after BSSRO in (
Competing interests
The authors declare that they have no competing interests.
Authors' contributions
DK conceived the study, analyzed the data, generated the figures, and contributed to writing the manuscript. DN conceived the study, performed the molecular genetic studies, analyzed the data, generated the figures, and contributed to writing the manuscript. YT and TK performed most of the operations on the patients in this study and analyzed the data. SK performed the molecular genetic studies and conceived the study. Kl and KF participated in conceiving the design, analyzed the data, and edited the manuscript. All the authors read and approved the final manuscript.
Footnotes
Abbreviations
Acknowledgements
We acknowledge Mr. Michael Arends for assistance with editing the manuscript. We are grateful to the volunteers for their participation in the study and anesthesiologists and surgeons in the Department of Oral Health and Clinical Science, Division of Dental Anesthesiology, Orofacial Pain Center, Suidoubashi Hospital, Tokyo Dental College, for collecting the clinical data. This work was supported by grants from the Ministry of Health, Labour and Welfare of Japan (H21-3jigan-ippan-011, H22-lyaku-015), Ministry of Education, Culture, Sports, Science and Technology of Japan (20390162, 22790518, 23390377, and 25116532), Smoking Research Foundation, Naito Foundation, Astellas Foundation for Research on Metabolic Disorders, and Mitsubishi Foundation.