Genome-Wide Association Study Identifies a Susceptibility Locus for Comitant Esotropia and Suggests a Parent-of-Origin Effect

Abstract

Purpose: To identify genetic variants conferring susceptibility to esotropia. Esotropia is the most common form of comitant strabismus, has its highest incidence in European ancestry populations, and is believed to be inherited as a complex trait.

Methods: White European American discovery cohorts with nonaccommodative (826 cases and 2991 controls) or accommodative (224 cases and 749 controls) esotropia were investigated. White European Australian and United Kingdom cohorts with nonaccommodative (689 cases and 1448 controls) or accommodative (66 cases and 264 controls) esotropia were tested for replication. We performed a genome-wide case–control association study using a mixed linear additive model. Meta-analyses of discovery and replication cohorts were then conducted.

Results: A significant association with nonaccommodative esotropia was discovered (odds ratio [OR] = 1.41, P = 2.84 × 10⁻⁰⁹) and replicated (OR = 1.23, P = 0.01) at rs2244352 [T] located within intron 1 of the WRB (tryptophan rich basic protein) gene on chromosome 21 (meta-analysis OR = 1.33, P = 9.58 × 10⁻¹¹). This single nucleotide polymorphism (SNP) is differentially methylated, and there is a statistically significant skew toward paternal inheritance in the discovery cohort. Meta-analysis of the accommodative discovery and replication cohorts identified an association with rs912759 [T] (OR = 0.59, P = 1.89 × 10⁻⁰⁸), an intergenic SNP on chromosome 1p31.1.

Conclusions: This is the first genome-wide association study (GWAS) to identify significant associations in esotropia and suggests a parent-of-origin effect. Additional cohorts will permit replication and extension of these findings. Future studies of rs2244352 and WRB should provide insight into pathophysiological mechanisms underlying comitant strabismus.

Strabismus is the pathological misalignment of the eyes and affects up to 4% of the population. All but a few percent of cases of strabismus are comitant, in which the angle of misalignment between the two eyes remains relatively constant with changes in gaze direction.^1,2 Most individuals with comitant strabismus have esotropia (ET, inward eye deviation) or exotropia (XT, outward eye deviation). The prevalence of ET is ∼2.5% among White populations of European ancestry and 0.5% among Africans and Asians. By contrast, XT has a prevalence of ∼1% among Africans, African Americans, and White populations of European ancestry, and is only slightly higher (∼1.2%) among Asians.³ 

Comitant strabismus is associated with poor binocular vision and amblyopia (uniocular visual neglect), and individuals with amblyopia have a much higher lifetime risk of bilateral visual impairment.⁴ Comitant strabismus can also disturb normal interpersonal interactions, resulting in poor self-esteem, social anxiety, phobias,^5,6 and limited employment options,⁷ leading to subtle but pervasive losses of productivity in the population.⁸ While glasses and surgery are standard treatments for comitant strabismus, many patients will require more than one procedure over a lifetime.^9,10 The pathogenesis of comitant strabismus remains largely unknown, and while opposing theories relevant to the infantile form of ET have long been debated,^11–15 they remain unproven. 

There is a significant heritable component to comitant strabismus, and familial clustering was described as early as approximately 400 BCE by Hippocrates.¹⁶ Subsequent family-based and twin studies have supported a genetic contribution,^17–20 with the relative risk for first-degree relatives of a comitant strabismus-affected proband estimated to fall between 3 and 5.^19,21 Moreover, heritability remains significant following correction for environmental risk factors.^22,23 Despite a clear genetic contribution to comitant strabismus, there are limited data supporting Mendelian segregation. Linkage analysis of two White pedigrees of European ancestry, one with unspecified comitant strabismus and one with ET, led to reports of linkage to chromosome 7p22.1 (STBMS1 locus, Online Mendelian Inheritance in Man [OMIM]185100) under a recessive²⁴ and a dominant²⁵ inheritance model, respectively. A third study combined small Japanese ET and XT pedigrees and reported genetic loci on chromosomes 4q28.3 and 7q31.2, with stronger significance under a parent-of-origin linkage model,^26,27 and suggested that MGST2 and WNT2 are potential candidates.²⁸ Beyond this, the genetic contributions to comitant strabismus remain undefined. This is in contrast to the rare forms of congenital incomitant strabismus in which ocular misalignment is associated with an angle of deviation that changes in different directions of gaze; these are often inherited as Mendelian traits and can result from gene mutations that perturb ocular motor neuron and axon development, resulting in primary malfunction in the ocular motor output pathways.^29–32 

Given the high prevalence, significant morbidity, and familial clustering of comitant strabismus, and data supporting ET as a complex trait, we conducted a genome-wide association study (GWAS) of ET in participants of White European ancestry as a first step to define its genetic architecture and unravel its pathogenesis. 

Subjects and Methods

The study was approved by the local Institutional Review Boards of Boston Children's Hospital (Boston, MA, USA), The Cleveland Clinic (Cleveland, OH, USA), the National Research Ethics Service (Leicestershire, Northamptonshire, UK), the Rutland Research Ethics Committee (Leicestershire, Northamptonshire, UK), and the Human Research Ethics Committee, Royal Victorian Eye and Ear Hospital (East Melbourne, Victoria, Australia), Princess Margaret Hospital (Perth, Western Australia), and the Sir Charles Gairdner Hospital (Perth, Western Australia, Australia). All investigations were conducted in accordance with the principles of the Declaration of Helsinki. 

Inclusion Criteria

Inclusion criteria for study of nonaccommodative ET required one of the following: (1) manifest or intermittent nonaccommodative ET of any size; (2) manifest or intermittent partially accommodative ET of any size; (3) esophoria ≥ 10 prism diopters; or (4) status post surgery for comitant ET. Inclusion criteria for study of accommodative ET required manifest ET that reduced with hyperopic correction (with or without bifocals) to a range such that fusion may be achieved (<10 prism diopters) or be reliably demonstrated. Inclusion for the purpose of testing association among individuals with amblyopia or hyperopia (independent of whether the ET was accommodative or nonaccommodative) required: (1) for amblyopia, a two or more line difference in best-corrected visual acuity between the two eyes, or a strong fixation preference in those unable to perform recognition acuity at time of first examination, or a record of diagnosis of or management for amblyopia; (2) for hyperopia, a refractive error of ≥ +3.50 diopters in either eye at any age. Definitions of ET subtypes are provided in Supplementary Methods. 

Exclusion Criteria

Exclusion criteria for any of the studies included (1) structural ocular abnormality causing acquired vision loss; (2) structural brain abnormality as determined by neuroimaging (normal neuroimaging was documented in 82 cases within the discovery cohort and 4 cases within the replication cohort); (3) conditions causing occlusion of the eye and leading to deprivation amblyopia; (4) molecularly defined genetic syndromes or other diagnoses associated with strabismus such as trisomy 21 or craniosynostosis; (5) consecutive ET; (6) other defined nonheritable etiology of strabismus. 

Phenotyping

Phenotyping was based on a combination of (1) participant examinations by an ophthalmologist, optometrist, or an orthoptist with the exception of 148 participants who were enrolled based on documentation of ET surgery, (2) participant questionnaires, and (3) review of additional medical records when available. 

Discovery Cohort Cases

A total of 1174 participants who self-reported as White of European ancestry were enrolled into the discovery cohort: 1105 from Boston Children's Hospital, 52 from Cole Eye Institute (Cleveland Clinic), 5 from Children's Hospital of Philadelphia, and 12 self-referred. After quality control (QC) procedures the total number of participants in the discovery cohorts was 1050 (Supplementary Table S1). 

Replication Cohort Cases

A total of 856 participants who self-reported as White of European ancestry were enrolled into the replication cohort: 745 from private ophthalmologists and public hospitals in New South Wales, Western Australia, Victoria, and Tasmania, Australia; 111 from Leicester, Devon and Bradford Hospitals, Leicester, United Kingdom. After QC procedures the total number of participants in the replication cohorts was 755 (Supplementary Table S1). 

Controls

Control genotypes generated on Illumina SNP microarrays were derived from the database of Genotypes and Phenotypes (dbGaP) (Supplementary Table S1; Supplementary Methods). 

Genotyping and Imputation

Discovery and replication cohorts were genotyped on Illumina Infinium Human OmniExpress_24v1-0 array (San Diego, CA, USA). dbGaP control cohorts had been genotyped on Omni 2.5 Versions 4v_1H and 4v_1D arrays with the exception of the eMERGE cohort, which had been genotyped on OmniExpress Version 12-v1_c. The OmniExpress platform used for genotyping cases was well represented on the larger panels on which the controls were genotyped, and genotypes from the discovery, replication, and control cohorts underwent rigorous QC measures. The genotypes were then merged, providing data for 337,204 single nucleotide polymorphisms (SNPs) in common among the two discovery cohorts, two replication cohorts, and ancestry-matched controls. Only cases and controls passing the QC measures (Supplementary Fig. S1; Supplementary Methods) were included in the final association tests. The SNPs overlapping between all platforms were used as input for imputation against 1000 Genomes phase I integrated variant set release (v3) European samples (provided in the public domain by the Center of Statistical Genetics, http://csg.sph.umich.edu//abecasis/MaCH/download/1000G.2012-03-14.html). A total of 14,882,799 SNPs were successfully imputed using IMPUTE2 program (version 2) (available in the public domain, https://mathgen.stats.ox.ac.uk/impute/impute_v2).³³ Of these, 6,470,615 imputed SNPs with information metric > 0.3 and minor allele frequency (MAF) > 0.01 were used. Notably, QC measures, including differential missingness testing, control–control associations, simultaneous imputation of cases and controls, and calculation of genomic inflation lambda (λ = 1.01), were conducted, in part, to limit the possibility of false-positive results that could stem from genotyping of cases and controls separately (see Supplementary Methods for details). 

Genetic Association and Meta-Analysis

Following QC, statistical tests for association were carried out in the discovery and then in the replication cohort using the SVS suite (Golden Helix, Bozeman, MT, USA). Single marker analyses for the genome-wide data were carried out using mixed linear additive model. We used EMMAX method (Efficient Mixed-Model Association eXpedited),³⁴ incorporating a precomputed kinship matrix, to account for known relatedness, as a random effect and the first three principal components as covariates. X chromosome SNPs were first analyzed separately in males and females and subsequently combined by meta-analysis. Because we used publicly available controls from dbGaP and hence some of the controls might have strabismus, we corrected for the association tests' results assuming a prevalence of 2% in control data. In addition, we repeated the statistical tests using an independent set of control individuals from the Health and Retirement control cohort. The SNPs that had significant or suggestive levels of association initially were comparable to the results using the new control set, indicating that our findings do not represent an artifact resulting from the use of historic controls. Manhattan and quantile-quantile (Q-Q) plots were generated with SVS suite, cluster plots were generated by Genome-Studio Software (Illumina, San Diego, CA, USA), and regional association plots were generated with LocusZoom (provided in the public domain by University of Michigan, http://locuszoom.sph.umich.edu/locuszoom/).³⁵ SNPs reaching significant or suggestive levels of association in the discovery cohort were checked in the replication cohort. In addition, we performed meta-analysis between the discovery cohort and the replication cohort using SVS suite. We performed both fixed-effects (inverse-variance–based approach using the effect size for the risk alleles and standard errors) and random-effect meta-analyses by calculating the random-effects variance component. The presence of heterogeneity between the two cohorts was investigated using the I² statistic. We also investigated our association data for potential replication of previously reported GWAS hits identified for ocular phoria.³⁶ Lastly, we confirmed the genotypes at the significantly associated SNP (rs2244352) by genotyping 15 US cases, 9 cases from Australia, and 6 cases from the United Kingdom (that were part of the cases included in the GWAS) using Sanger sequencing. 

Analysis and Functional Prediction of rs2244352

To evaluate the potential functional consequences of SNP rs2244352, we first used GTEx (The Genotype-Tissue Expression [GTEx] project) portal to investigate if it acts as an expression quantitative trait locus (eQTL) in different body tissues as well as in different regions of the brain (provided in the public domain by The Broad Institute of MIT and Harvard, http://www.gtexportal.org/).³⁷ We also used the Regulome DB (provided in the public domain by Stanford University, http://www.regulomedb.org/),³⁸ which compiles eQTLs, ChIP-seq, and DNase-seq data from the ENCODE (Encyclopedia of DNA Elements) project,³⁹ the Roadmap Epigenomics Mapping Consortium⁴⁰ and HaploReg v4.1 (provided in the public domain by the Broad Institute, http://compbio.mit.edu/HaploReg)⁴¹ to predict the possible function of rs2244352 as well as those of other SNPs in linkage disequilibrium. Potential chromatin interactions between the genomic region containing rs2244352 and its neighboring regions were investigated using published Hi-C chromatin interactions.⁴² 

To investigate if a parent-of-origin effect of the rs2244352 [T] allele among the nonaccommodative ET discovery cohort exists, we identified probands heterozygous for the [T] allele for whom parental DNA was available (202 informative trios), and genotyped the parents at rs2244352 using a Taqman assay (Life Technologies, Carlsbad, CA, USA). We compared the parental inheritance pattern at rs2244352 to those of 408 informative trios from unrelated syndromic cranial dysinnervation disorder and structural heart birth defect cohorts that underwent whole genome sequencing through the Gabriella Miller Kids First Pediatric Research Program. Unfortunately, parental DNA was not available from parents of the replication cohort. 

Results

Discovery and Replication Cohorts

A total of 2030 individuals diagnosed with comitant ET and self-reported to be White of European ancestry were enrolled (1052 females and 978 males, 2–91 years of age); 1174 from the United States served as the initial cohort for discovery and 856 from Australia and the United Kingdom served as the initial cohort for replication prior to QC procedures. Hypothesizing that accommodative and nonaccommodative forms of ET have different underlying etiologies, participants were divided into two subgroups: the “nonaccommodative ET” group, and the “accommodative ET” group (see Subjects and Methods section). 

Genome-Wide Association Results

The nonaccommodative ET discovery cohort included 826 cases from the United States and 2991 controls from the Genetic Variation in Refractive Error Substudy of the National Eye Institute Age-Related Eye Disease Study (AREDS) and the Health and Retirement Study (HRS). SNP rs2244352 [T] (the minor [T] allele is also the reference allele) reached genome-wide significance in the nonaccommodative discovery cohort with P = 2.84 × 10⁻⁰⁹ (Fig. 1; Table), generating an odds ratio (OR) of 1.41 assuming a 2% prevalence of ET in the control population. The quantile-quantile (Q-Q) plot for the nonaccommodative ET GWAS is shown in Supplementary Figure S2 and the cluster plot of rs2244352 is shown in Supplementary Figure S3. 

Figure 1

 

Association results of the nonaccommodative esotropia discovery cohort. (A) Manhattan plot of association results from the discovery data for the nonaccommodative esotropia cohort. The y-axis represents the −log10 (P values). The bolded upper horizontal line indicates the threshold for genome-wide significance (P < 5 × 10⁻⁸) and the lower line indicates the threshold for suggestive association (P < 1 × 10⁻⁵). The SNP rs2244352 had a P = 2.84 × 10⁻⁹. (B) Regional association plot from the discovery nonaccommodative esotropia cohort. The top-ranked SNP, rs2244352, is shown as a solid purple diamond. The color scheme indicates linkage disequilibrium between rs2244352 and other SNPs in the region using the r² value calculated from the 1000 Genomes project. The y-axis is the −log10 (P value) computed for the displayed SNPs.

Table

 

Top SNPs at Each Locus Found to Pass a Suggestive Level of Significance (P ≤ 1 × 10⁻⁵) in the Nonaccommodative Esotropia Discovery Cohort

 

We investigated the association of the 105 autosomal SNPs that reached P ≤ 1 × 10⁻⁵ (including rs2244352) in the nonaccommodative ET discovery cohort in a replication cohort composed of 689 cases from the United Kingdom and Australia and 1448 controls from a subset of the Fuchs' Endothelial Corneal Dystrophy (FECD) cohort. SNP rs2244352 showed evidence of replication (P = 0.006, OR = 1.23) (Table; Supplementary Table S2). Of note, rs2244352 was a genotyped variant and there were no missing calls in cases and controls in the discovery or replication cohorts. We then conducted a meta-analysis of the association results between the nonaccommodative ET discovery and replication cohorts; the meta-analysis P value at rs2244352 was 9.58 × 10⁻¹¹, OR = 1.33 (Table). Finally, given that the nonaccommodative cohort was heterogeneous, we tested whether any of the three largest subtypes (manifest, intermittent, and infantile ET) drove the association signal at rs2244352. P values at rs2244352 were 3.11 × 10⁻⁶ for manifest ET alone, 8.17 × 10⁻⁸ when intermittent and manifest ET were combined, and 3.23 × 10⁻⁷ when infantile and manifest ET were combined, despite the infantile ET cohort being larger than the intermittent (Supplementary Table S1). This suggests that intermittent ET was a stronger driver than infantile ET in our cohort. 

To test for associations with accommodative ET, the discovery cohort included 224 cases from the United States and 749 controls from the HRS that did not overlap with the nonaccommodative controls. SNP rs912759 [T], an intergenic variant located within a ∼2.3-Mb gene desert on chromosome 1, was the most significant, with P = 6.53 × 10⁻⁰⁷, OR = 0.59 (Supplementary Table S3). rs912759 was a genotyped variant and there were no missing calls in cases and controls in the discovery or replication cohorts. The Q-Q plot for the accommodative GWAS, the cluster plot, and regional association results for rs912759 are shown in Supplementary Figures S2, S3, and S4, respectively. We investigated the top SNPs from the accommodative ET discovery cohort (P ≤ 1 × 10⁻⁰⁵) in a replication cohort of 66 cases from Australia and the United Kingdom and 264 controls from the eMERGE cohort (Supplementary Table S1). The discovery association results for rs912759 replicated (P = 0.008, OR = 0.59) with consistent effect size and direction. When meta-analysis of the combined accommodative ET cohorts was conducted, rs912759 exceeded genome-wide significance (overall P = 1.89 × 10⁻⁰⁸, OR = 0.59) (Supplementary Table S3; Supplementary Fig. S4). 

Although we hypothesize that accommodative and nonaccommodative ET do not have the same genetic etiology, we assumed that some loci might overlap and hence we combined both cohorts in a meta-analysis, which maximized our overall sample size. The meta-analysis did not reveal any new significant associations and reduced the significance and increased the heterogeneity at the loci with significant associations: rs2244352, the SNP associated with nonaccommodative ET, had an overall fixed effect P = 2.07 × 10⁻¹⁰, random effect P = 1.7 × 10⁻⁰⁵, and I² statistic = 0.31; rs912759, the SNP associated with accommodative ET, had an overall fixed effect P = 0.0077, random effect P = 0.079 and I² statistic = 0.88. Supplementary Table S4 compares the results at rs2244352 and rs912759 between the accommodative and nonaccommodative ET cohorts. 

We also investigated whether specific associations existed for amblyopia or hyperopia within our discovery and replication cohorts, independent of whether the ET was accommodative or nonaccommodative. For amblyopia, we identified and analyzed data from 386 and 45 cases within the discovery and replication cohorts, respectively, matched to 1700 and 180 controls. Neither the primary association test in the discovery cohort nor the meta-analysis identified SNPs that reached statistical significance (Supplementary Table S5). For hyperopia, we identified and analyzed data from 361 and 172 cases within the discovery and replication cohorts, respectively, matched to 1700 and 690 controls. Again, no SNP in the primary or meta-analysis showed significant association (P < 5 × 10⁻⁸) (Supplementary Table S5). We did not replicate the associations on chromosome 15q14 or 8q12 reported for hyperopia, perhaps because our cohort was much smaller than in the original study.⁴³ We were also underpowered to test for specific associations with ocular phoria, but we did not find an association with the recently reported locus on chromosome 6p22³⁶ for any of the phenotypes we tested in the current study. 

Functional Assessment of Significantly Associated Variants

rs2244352, associated with nonaccommodative ET, is located within DNase hypersensitivity clusters in multiple tissues, including fetal brain.³⁹ DNA footprinting and Chip-Seq experiments predict that the [T] allele alters the sequence of several regulatory motifs (Fig. 2A), and the binding of POLR2A, TAF1, NRF1, E2F1, and USF1 transcription factors in a variety of cell lines.^39,44 Moreover, interrogation of Hi-C chromatin interaction data⁴⁵ reveals that WRB (within which rs2244352 resides) and other eQTL genes in the region (LCA5L, SH3BGR, HMNG1, and BRWD1) fall within the boundaries of a topologically activated domain (TAD) (Fig. 2B). 

Figure 2

 

DNA-footprinting prediction of the effect of the risk allele [T] on regulatory motifs and two-dimensional heatmap of chromatin interaction in the neighborhood of SNP rs2244352. (A) DNA footprinting predicts that the [T] allele alters the sequence of c-Ets-1 and PU-1 regulatory motifs.³⁸ Sequence logos were created by CENTIPEDE.^38,44 (B) UCSC genome browser shot of chromosome 21 (NCBI Build 37; provided in the public domain by University of California-Santa Cruz, https://genome.ucsc.edu/) created by 3D Genome Browser (provided in the public domain by Pennsylvania State University, http://promoter.bx.psu.edu/hi-c)⁷⁸ showing WRB and surrounding genes (PSMG1, BRWD1, HMGN1, LCA5L, and SH3BGR) included within a topologically associating domain (TAD) as defined by chromatin conformation capture data (Hi-C profiles). The two-dimensional heatmap represents normalized Hi-C interaction frequencies.

rs2244352 was assigned the score [1b] in the Regulome Database. This is the second highest of 14 possible scores supporting a functional role for a SNP, and is the highest score assigned to any SNPs in WRB, LCA5L, and SH3BGR, the genes with which rs2244352 is in high linkage disequilibrium (LD) (Supplementary Table S6). We interrogated the GTEx portal⁴⁶ to investigate if the rs2244352 risk allele [T] acted as an eQTL in different body tissues. The rs2244352 risk allele [T] correlates largely with increased expression of WRB, LCA5L, SH3BGR, and BRWD1 in multiple tissues (Supplementary Table S7). We specifically investigated if rs2244352 [T] affects expression of WRB or other genes within regions of the brain for which GTEx has data. Similar to findings in nonbrain tissues, rs2244352 [T] increases the expression of WRB, LCA5L, and SH3BGR in all brain regions examined with the exception of the frontal cortex, where WRB expression is noted to be decreased. By contrast, HMGN1 and PSMG1 expression is decreased (Fig. 3; Supplementary Tables S7, S8). 

Finally, the nongenic variant rs912759, which is associated with accommodative ET in the meta-analysis, is predicted by HaploReg v4.1⁴¹ to alter the regulatory motifs for Nkx3_5, Nkx2_11, Hbp1, and HMG-IY_2, but there are currently no data to support its action as an eQTL. 

Testing the Parent-of-Origin Effect for rs2244352

rs2244352 is located 73 base pairs 3′ to CGI-2, one of multiple 5^mCpG islands within WRB (CGI-2; Chr21: 40757603-40757900). Although CGI-2, -3, and -5 are differentially methylated, only CGI-2 is located within the differentially methylated region (DMR) for which WRB was classified as a candidate maternally imprinted gene.^47–49 Several studies have reported a correlation between increased maternal methylation of WRB CpG islands and decreased WRB expression in select tissues,^48,50 while others have found biallelic expression of WRB and its neighboring genes, suggesting that the maternal imprinting of WRB is isoform- and/or tissue-specific.^48,49 Thus, we hypothesized that paternal inheritance of the unmethylated at-risk rs2244352 [T] allele increases expression of WRB and/or other target genes in the TAD and increases the risk of nonaccommodative ET. To test this hypothesis, we asked if there was a parent-of-origin effect of the rs2244352 [T] allele among the nonaccommodative ET discovery cohort. We identified probands from both the discovery and replication cohorts who were heterozygous for the [T] allele and for whom parental DNA was available, and genotyped the parents at rs2244352 using a Taqman assay (Life Technologies). Of 202 informative trios (114 ET, 42 infantile ET, 44 intermittent ET, and 2 acquired ET), 120 (59.4%) probands inherited the [T] allele paternally and 82 (40.6%) inherited it maternally. This paternal skew was significant (P = 0.023, χ²) when compared to 408 informative trios from unrelated birth defect cohorts, of which 193 probands (47.3%) inherited the [T] allele paternally, and 215 probands (52.7%) inherited it maternally. Thus, the OR is 1.63 if the [T] allele is paternally inherited, and is 0.6 if maternally inherited. 

Discussion

Using discovery cohorts from the United States and replication cohorts from the United Kingdom and Australia, we report two significantly associated SNPs on chromosomes 21 and 1 for nonaccommodative and accommodative ET, respectively. We chose the discovery, replication, and control cohorts from Caucasian populations of European ancestry and conducted principal component analysis to confirm absence of population stratification. Because we used controls from the database of Genotypes and Phenotypes (dbGaP), rigorous QC metrics were taken to ensure that allele frequency differences were reflective of true locus-specific associations and not the result of differential genotype calling (batch effects) or population stratification. The two significant SNPs (rs2244352 and rs912759) were genotyped in both discovery and replication cases and controls with no missing calls. Together, these measures provide greater confidence in the discovery and replication results.⁵¹ 

Our most significant association was with SNP rs2244352, which we found to be associated with nonaccommodative ET. rs2244352 falls at the midpoint of intron 1 of the WRB (tryptophan rich basic protein) gene on chromosome 21. WRB is the homolog of yeast Get1; it is located on the surface of the endoplasmic reticulum (ER) and, together with calcium-modulating ligand (CAML), forms a transmembrane receptor complex that mediates the posttranslational TRC40-mediated insertion of tail-anchored (TA) proteins into the ER membrane.^52–54 TA proteins have a single C-terminal transmembrane domain that anchors them to a variety of organelles in both secretory and endocytic pathways; they constitute 5% of total membrane proteins and play critical roles in cell function. Consistent with this, WRB is widely expressed in all fetal and adult tissues examined.^37,55 Wrb loss-of-function zebrafish mutants are blind and deaf secondary to a reduction in synaptic ribbons and surrounding vesicles in the sensory retina and hair cells, and die by 10 days post fertilization.^56,57 Mice with conditional Wrb deletion in inner-ear sensory hair cells also have fewer ribbon-associated vesicles, impaired hair cell exocytosis, and hearing loss.⁵⁸ By contrast to these Wrb loss-of-function phenotypes, our data support the contribution of increased WRB expression to susceptibility to nonaccommodative ET in humans. Importantly, however, rs2244352 falls within the boundaries of a TAD that, in addition to WRB, includes LCA5L, SH3BGR, HMGN1, BRWD1, and PSMG1 (Fig. 2B; Supplementary Table S9). Thus, it remains possible that the susceptibility to ET is mediated through the action of rs2244352 on the regulation of one or more of these other genes. 

The associated SNP rs2244352 maps to chromosome 21q22.2 within the Down syndrome critical region,^59–61 and individuals with Down syndrome have a much higher (∼20%) incidence of ET compared to the general population.^62,63 Thus, it is possible that inheritance of a third copy of chromosome 21 increases expression of WRB and/or nearby genes and contributes to the ET found in this patient population. 

WRB is a maternally imprinted gene; the maternal allele is methylated and silenced while the paternal allele is unmethylated and expressed.^47–49 It is known that the effects of methylation on gene expression can vary with tissue and developmental stage,^64,65 that imprinted genes often function as developmental regulators,⁶⁶ and that genomic imprinting can contribute to complex traits.^67–70 Consistent with this, we found a paternal skew in the inheritance of the rs2244352 [T] allele in cases of nonaccommodative ET. These findings suggest that paternal inheritance of the at-risk allele increases WRB expression in the critical tissue(s) and at the critical developmental timepoint(s) to enhance susceptibility to nonaccommodative ET. Moreover, WRB is the only gene in the TAD demonstrated to be maternally imprinted,^48,49 adding strength to the hypothesis that WRB is the gene providing susceptibility to nonaccommodative ET. While our data support a protective effect of the maternal allele, additional studies would be required to confirm this finding. 

Methylation status is known to be influenced by the environment; and the specific environmental factors that correlate with increased risk of developmental strabismus, including maternal smoking, advanced maternal age, and premature birth,^{22,23,71–73} have been found to also alter the offspring's methylation status.^74–76 Among these, a meta-analysis has specifically identified reduced methylation of WRB in offspring of mothers who smoked during pregnancy,⁷⁷ raising the possibility that genetic and epigenetic influences are working through a common pathway to increase the risk of developmental strabismus. 

Future studies of additional developmental strabismus cohorts are necessary to further replicate the significance of rs2244352 among White European individuals with nonaccommodative ET, determine if it contributes to ET among other ethnicities, and confirm its potential imprinting. These studies will help to establish whether a specific ET subtype primarily drives the signal, and determine if this SNP may also provide susceptibility to XT. Additional studies are also necessary to determine if rs2244352 or a different SNP in high LD with rs2244352 is the causative SNP, and whether altered expression of WRB or another gene in the region is responsible for the increased risk of nonaccommodative ET. Finally, larger studies are necessary to determine if rs912759 replicates in accommodative ET, and to determine if any of the additional SNPs identified in either the nonaccommodative and accommodative cohorts that were suggestive of association are truly significant. 

Figure 3

 

Box plots of WRB, LCA5L, and PSMG1 gene expression for the cis-eQTL rs2244352 in brain cortex. Gene expression levels are stratified by genotype at the SNP rs2244352 cis-eQTL in brain cortex, based on GTEx analysis release V7 (provided in the public domain by The Broad Institute of MIT and Harvard, http://gtexportal.org/). The reference and minor allele [T] is associated with increased risk for nonaccommodative esotropia (P = 2.84 × 10⁻⁰⁹), increased expression of WRB and LCA5L (P = 0.0015, 7.5 × 10⁻¹⁰, respectively), and decreased expression of PSMG1 (P = 0.0099). Homo Ref = [TT]; “Het” = [TG]; “Homo Alt” = [GG]. eQTL, expression of a quantitative trait locus; SNP, single nucleotide polymorphism.

Acknowledgments

The authors thank the families for their participation; Melissa Hoche Culver, Christina Fitzgerald, Michelle Cirioni Janes, Karli Kibby, Kelsey McEntee, Ogo Ndugba, Andrew Perechocky, Allison Rose, Chris Wirth, and Brian Young for research assistance; Kailee Algee, Kamil Gorski, Fiona Harrigan, Jessica Kane, Katelynn MacNeill, Lorraine O'Keefe, Frances Pantano-Abele, Caitlin Rook, Jonathan Russell, and Mariette Tyedmers for orthoptic support; Joel Hirschhorn, S. Alessandro Di Gioia, and members of the Engle lab for enlightening discussions; and the HMS Ocular Genomics Institute and the Broad Institute for technical expertise and support. 

The datasets used in this manuscript were obtained from dbGaP National Institutes of Health (NIH) data repository (provided in the public domain by NIH, http://www.ncbi.nlm.nih.gov/gap). The authors thank the following investigators who contributed the phenotype data and DNA samples from their original studies and the funding organizations that supported these studies: dbGaP accession numbers phs000429.v1.p1 (Dwight Stambolian, PI; NIH contracts HHSN268200782096C and HHSN268201100011I), phs000381.v1.p1 (David Carey, PI; Geisinger Clinic funding), phs000421.v1.p1 (Natalie Afshari, John Gottsch, Sudha K. Iyengar, Nicholas Katsanis, Gordon Klintworth, and Jonathan Lass, PIs; National Eye Institute (NEI) Grants R01EY016482, R01EY016514, and R01EY016835, and NIH contracts HHSN268200782096C and HHSN268201100011I), and phs000428.v1.p1 (David Weir, PI; RC2 AG036495 grant). The authors also acknowledge the NIH, Pediatric Cardiac Genomics Consortium, which is supported by Grants U01-HL098188, U01-HL098147, U01-HL098153, U01-HL098163, U01-HL098123, and U01-HL098162. 

Supported by NIH Grants R01EY015298 and R01EY027421 (ECE), Boston Children's Hospital IDDRC (NIH 1U54HD090255), the NIH Common Fund's Gabriella Miller Kids First Pediatric Research Program (ECE, SUM), Boston Children's Hospital Ophthalmology Department Foundation Discovery Award (ECE), Al-Habtoor Fellowship-Dubai Harvard Foundation for Medical Research (SS), NIH Grant K08EY027850 (MCW), the Australian National Health and Medical Research Council Project Grant APP1031920 (DAM) and Early Career Fellowship support (SY), the Australian Research Council Future Fellowship Scheme (SM), the Ophthalmic Research Institute of Australia (DAM), Medical Research Council Grant MR/N004566/1, and the Ulverscroft Foundation (IG). Centre for Eye Research Australia receives Operational Infrastructure Support from the Victorian Government. ECE is an investigator of the Howard Hughes Medical Institute. 

Disclosure: S. Shaaban, None; S. MacKinnon, None; C. Andrews, None; S.E. Staffieri, None; G.D.E. Maconachie, None; W.-M. Chan, None; M.C. Whitman, None; S.U. Morton, None; S. Yazar, None; S. MacGregor, None; J.E. Elder, None; E.I. Traboulsi, None; I. Gottlob, None; A.W. Hewitt, None; D.G. Hunter, None; D.A. Mackey, None; E.C. Engle, None 

Appendix

Members of the Strabismus Genetics Research Consortium

United States: Anna Baglieri, Brenda Barry, Sarah Bekele, Sarah E. Breau, Kimberley Chan, Frances Corkin, Linda R. Dagi, Alexandra Elliott, Janet Esligar, Caroline Fang, Anne B. Fulton, Gena Heidary, Suzanne Johnston, Melanie Kazlas, Danielle M. Ledoux, Richard L. Levy, Iason S. Mantagos, Kathryn B. Miller, Monte Mills, Darren Oystreck, Christina S. Petersen, Robert A. Petersen, Carrie E. Pierce, Aparna Raghuram, Richard Robb, Josephine C. Sandoval, Sonia Sethee, Ankoor S. Shah, Lois E. H. Smith, Melissa Toffoloni, Deborah K. VanderVeen, Sarah Whitecross, Rupa K. Wong, Carolyn Wu 

Australia: Julie Barbour, Linda Clarke, Joanne C. Dondey, Maree Flaherty, John Grigg, Kate Hanman, Michael Haybittel, Robyn V. Jamieson, Lisa S. Kearns, Lionel Kowal, Geoffrey C. Lam, Troy Lim Joon, John McKenzie, Loren Rose, Jonathan B. Ruddle, Lindsey Scotter, Neil E. Sinclair, Colleen Wilkinson, Robin Wilkinson 

United Kingdom: Viral Sheth, Mervyn G. Thomas