{'fraction_non_alphanumeric': 0.04213084112149533, 'fraction_numerical': 0.04085981308411215, 'mean_word_length': 4.935433769691591, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Variation in PTCHD2, CRISP3, NAP1L4, FSCB, and AP3B2 associated with spherical equivalent.

Purpose Ocular refraction is measured in spherical equivalent as the power of the external lens required to focus images on the retina. Myopia (nearsightedness) and hyperopia (farsightedness) are the most common refractive errors, and the leading causes of visual impairment and blindness in the world. The goal of this study is to identify rare and low-frequency variants that influence spherical equivalent. Methods We conducted variant-level and gene-level quantitative trait association analyses for mean spherical equivalent, using data from 1,560 individuals in the Beaver Dam Eye Study. Genotyping was conducted using the Illumina exome array. We analyzed 34,976 single nucleotide variants and 11,571 autosomal genes across the genome, using single-variant tests as well as gene-based tests. Results Spherical equivalent was significantly associated with five genes in gene-based analysis: PTCHD2 at 1p36.22 (p = 3.6 × 10−7), CRISP3 at 6p12.3 (p = 4.3 × 10−6), NAP1L4 at 11p15.5 (p = 3.6 × 10−6), FSCB at 14q21.2 (p = 1.5 × 10−7), and AP3B2 at 15q25.2 (p = 1.6 × 10−7). The variant-based tests identified evidence suggestive of association with two novel variants in linkage disequilibrium (pairwise r2 = 0.80) in the TCTE1 gene region at 6p21.1 (rs2297336, minor allele frequency (MAF) = 14.1%, β = –0.62 p = 3.7 × 10−6; rs324146, MAF = 16.9%, β = –0.55, p = 1.4 × 10−5). In addition to these novel findings, we successfully replicated a previously reported association with rs634990 near GJD2 at 15q14 (MAF = 47%, β = –0.29, p=1.8 × 10−3). We also found evidence of association with spherical equivalent on 2q37.1 in PRSS56 at rs1550094 (MAF = 31%, β = –0.33, p = 1.7 × 10−3), a region previously associated with myopia. Conclusions We identified several novel candidate genes that may play a role in the control of spherical equivalent. However, further studies are needed to replicate these findings. In addition, our results contribute to the increasing evidence that variation in the GJD2 and PRSS56 genes influence the development of refractive errors. Identifying that variation in these genes is associated with spherical equivalent may provide further insight into the etiology of myopia and consequent vision loss.

basis of these associations have yet to be experimentally validated [36].
Evolutionary theory and empirical data suggest that rare and low-frequency variants could have relatively large effects on the risk of developing complex traits, and collectively rare variants may make significant contributions to disease susceptibility in a given population [37]. However, these variants were poorly characterized by earlier generations of GWAS genotyping arrays [38]. Sequencing of the exome or whole genome directly assays all variants, including rare and low-frequency variants, but it is currently costly to perform on thousands of study subjects. Exome arrays, with enhanced resolution of coding variants, provide a less costly alternative to sequencing to evaluate the role of rare (minor allele frequency (MAF) < 1.0%) and low-frequency (MAF = 1.0 -5.0%) variants in phenotypes.
To examine the impact of rare and low-frequency coding variants on spherical equivalent, we conducted an exome array analysis on participants in the Beaver Dam Eye Study (BDES). Additionally, data permitting, we evaluated known regions of association reported in other populations for spherical equivalent and myopia.

Study participants:
The BDES is a population-based cohort established in 1987 to study age-related eye disorders [39]. A total of 5,924 individuals aged from 43 to 84 years were identified from 3,715 households through a private census of Beaver Dam, Wisconsin. Between 1988 and 1990, 4,926 (83.14%) individuals participated in the baseline examination. Familial relationships of 2,783 participants were confirmed and the data was assembled into pedigrees. To improve the sampling efficiency, we genotyped individuals with intraocular pressure or spherical equivalent values at the baseline, resulting in a nested sub-cohort of 2,032 individuals with the full spectrum of spherical equivalent values. The study adhered to the ARVO statement on human subjects and was approved by the institutional review board at the University of Wisconsin. The study complied with the guidelines of the Declaration of Helsinki. Informed consent was obtained from all participants in the study.
Clinical evaluation: All participants in the BDES received a complete eye exam that included standardized evaluation of refraction using the Humphrey 530 refractor (San Leandro, CA) [39]. When the automated refraction yielded visual acuity of 20/40 or worse, a modification of the Early Treatment Diabetic Retinopathy Study (ETDRS) protocol [40] was followed to obtain the best-corrected visual acuity. The mean spherical equivalent (MSE) in diopters for each eye was calculated using the standard formula (MSE = sphere + 0.5 × cylinder). Nuclear lens opacity was determined by grading slit-lamp lens photographs using a standard protocol, resulting in a five-level scale [41]. Eyes without a lens, with an intraocular lens, or with best-corrected visual acuity of 20/200 or worse were excluded. Age, sex, and education were also obtained at the baseline visit.
Genotyping and quality control: DNA from 2,032 BDES participants, 22 blind duplicate samples, and 44 HapMap samples were genotyped on the Illumina Infinium HumanExome-12 v1.1 (Illumina, Inc., San Diego, CA) array at the Genetic Resources Core Facility (GRCF), the Johns Hopkins School of Medicine, the Institute of Genetic Medicine (Baltimore, MD). Of the 2,032 BDES participants, 1,908 (93.9%) were successfully genotyped and had a call rate >98%. After related individuals (n = 126) and individuals with unresolved gender inconsistencies (n = 15), missing spherical equivalent measurements of either eye (n=118), differences in spherical equivalent > ± 4 D between the eyes (n = 15), and missing values in age, sex, education, or nuclear sclerosis (n = 74) were excluded, 1,560 individuals remained in the final analytical cohort.
Of the 242,901 variants genotyped, we excluded the following: sex chromosomal variants (n = 5,465), insertion/ deletion variants (n = 134), variants with call rates ≤98% (n = 6,027), monomorphic variants (n = 132,213), and duplicate variants (n = 314). The final data set included 98,748 polymorphic variants. We evaluated Hardy-Weinberg equilibrium (HWE) on these variants and reported the p value estimated from the HWE exact test for significant associations.
To evaluate genetic ancestry, we performed principal component analysis (PCA) using SMARTPCA in EIGEN-STRAT (version 4.2) [42]. All BDES participants genotyped were self-reported non-Hispanic white. The first two principal components were included in the analysis to control for potential population stratification, resulting in a genomic inflation factor (λ) of 1.02.
Trait definitions: Age, education, and nuclear sclerosis have known effects on spherical equivalent and were significant covariates in our previous linkage analysis [43]. Using linear models, we regressed age (mean centered), education (mean centered), sex, nuclear sclerosis (sum of the nuclear lens opacity score of both eyes), and the first two principal components on the average of spherical equivalent in the right and left eyes. The residuals from this regression were used as the quantitative outcome in our association tests. For replications of loci associated with myopia in prior GWASs, we performed logistic regression analysis that compared myopia to hyperopia and myopia to non-myopia, adjusting for the covariates.
We defined individuals with spherical equivalent < −1 D in either eye and < −0.5 D in the opposite eye as myopic and those with spherical equivalent > +1 D in either eye and > + 0.5 D in the opposite eye as hyperopic. Emmetropia was defined as spherical equivalent ranged from −1 D to +1 D in both eyes. The non-myopic group included emmetropes and hyperopes. Individuals who fell outside these categories were excluded from the logistic regression analysis.
Statistical analysis: We applied single-variant and gene-based tests to investigate the association of rare and low-frequency genetic variants with spherical equivalent as a quantitative trait. Single-variant association analysis was performed on variants with MAF ≥1.0% (n = 34,976) under an additive model using PLINK (version 1.07) [44]. To determine the number of independent tests, we tabulated the number of single nucleotide polymorphisms (SNPs) with pairwise linkage disequilibrium (LD) <0.2. There were 23,342 independent variants, which translated to a Bonferroni corrected significance threshold of p<2.1 × 10 −6 .
Gene-based analysis can be a powerful alternative approach to single-variant analysis in detecting effects of rare variants [45]. Rare variants in our analysis were defined as those with MAF of 1% or less. To examine the impact of rare variants on spherical equivalent, we used the following tests: (1) the burden test, where all variants were collapsed into a gene-based score [46], and (2) the sequence kernel association test (SKAT) [47]. Analysis was restricted to 11,571 genes with at least two rare variants. To optimize the statistical power of SKAT, we applied the Madsen-Browning weighting which slightly up-weights the rare variants [48]. A p value of <4.3 × 10 −6 , corresponding to a Bonferroni correction for 11,571 independent gene tests, was used to determine significance.
We evaluated gene regions with GWAS significant associations in published studies or the region spanning two genes if the GWAS association was intergenic [23][24][25][26][27][28][29][30][31][32][33][34]. For each known GWAS locus, we performed a single-variant association analysis for quantitative refraction and myopia, targeting (1) the exact variant or (2) variants in strong LD (r 2 >0.8) with the reported variants. Variants were considered as replicated at p<0.01 if the magnitude and direction of their effects in our study were consistent with those reported in previous studies.

RESULTS
Characteristics of the Beaver Dam Eye Study and the genotyped sub-cohort included in these analyses are presented in Table 1. The baseline mean spherical equivalent in the subset was −0.20 D with a standard deviation of 2.9 D (Appendix 1). Compared to the entire BDES cohort, our subset included a higher proportion (38.6% versus 21.5%, respectively) of myopes and a similar proportion of hyperopes (37.0% versus 37.8%, respectively).
A total of 11,571 gene regions were tested in the genebased analysis. In these analyses, five gene regions showed significant (p<4.3 × 10 −6 ) evidence of association (Table 2, Appendix 3) with spherical equivalent. The most significant association was with the FSCB gene (OMIM 611779) on chromosome 14q21.2 (p = 1.5 × 10 −7 ). Six protein-altering (one stop-gain and five missense) variants contributed to the gene-based finding.  We identified four additional gene regions with significant association with spherical equivalent: the AP3B2 gene (OMIM 602166) on chromosome 15q25.2 (p = 1.6 × 10 −7 ), the PTCHD2 gene (OMIM 611251) on chromosome 1p36.22 (p = 3.6 × 10 −7 ), the NAP1L4 gene (OMIM 601651) on chromosome 11p15.5 (p = 3.6 × 10 −6 ), and the CRISP3 gene on chromosome 6p12.3 (p = 4.3 × 10 −6 ). Although the SKAT and the burden test supported a significant association between rare variants in the AP3B2 gene and quantitative refraction, the p value was more significant in the burden test, possibly suggesting that the variants included in the test influence spherical equivalent in the same direction and with similar magnitudes. Although these findings are intriguing, due to the low frequency of rare variants in these four genes (  [30,34], 4q25 [24], and MIPEP (OMIM 602241) [26]. Variants of GJD2, RASGRF1, and PRSS56 had been previously associated with spherical equivalent and myopia. In our analysis, we had evidence of replication for two of these eight regions (Table 3). We replicated the association of rs634990 located 37.2 kb upstream of the GJD2 gene on chromosome 15q14 with spherical equivalent (Appendix 4). In the Beaver Dam Eye study, the minor allele C of rs634990 was associated with a lower spherical equivalent (MAF = 47%, β = -0. 29

DISCUSSION
In this study, we performed an exome array analysis on 1,560 BDES participants to identify rare and low-frequency variants associated with spherical equivalent. We discovered novel associations with rare variants in five gene regions, PTCHD2, CRISP3, NAP1L4, FSCB, and AP3B2, and suggestive associations with common variants in TCTE1. In addition, we replicated the previously reported association of GJD2 with spherical equivalent and myopia and discovered an association of spherical equivalent on PRSS56, a locus that had been previously identified for myopia.
The rare variants in the identified genes may not explain much of the phenotypic variation in quantitative refraction in the general population due to their low frequencies. However, the novel loci associated with spherical equivalent may not only be of great influence in a small number of families that harbor these rare variants but also provide valuable insights into novel genes and pathways involved in the etiology of refractive errors. The AP3B2 gene is highly expressed in adult RPE cells [54,55]. The AP3B2 protein is a subunit of the neuron-specific AP3 complex involved in neuronal membrane trafficking via mediating the direct transport of lysosomal proteins from the trans-Golgi network to the lysosomes and lysosome-related organelles [56,57]. A recent study found binding activity between the AP3B2 protein and a membrane scaffold protein, 4.1G, whose deficiency can lead to mislocalization of photoreceptor terminals and visual acuity impairments in mice [58]. The AP3B2 gene has been implicated in the lysosome pathway [59][60][61][62]. At the significance level of 0.01, we also found evidence of association with two other genes in the lysosome pathway in the gene-based analysis: SLC17A5 (OMIM 604322) on chromosome 6q13 (p = 2.0 × 10 −3 ) and MFSD8 (OMIM 611124) on chromosome 4q28.2 (p = 2.7 × 10 −3 ; Appendix 5). These findings may suggest the potential involvement of one or more genes of the lysosome pathway in the development of refractive errors.
The gene PTCHD2 on chromosome 1p36.22 is located in the myopia locus 14 (MYP14, OMIM 610320), a known linkage region for myopia and spherical equivalent [63,64]. PTCHD2 is also overlapped with a linkage region for the Volkmann type of cataract (CTRCT8, OMIM 115655) and a locus for primary infantile glaucoma (GLC3B, OMIM 600975), suggesting its potential role in the development of several ocular conditions. This gene is highly expressed in retinal ganglion cells [54]. The PTCHD2 protein interacts with several other proteins in the Hedgehog signaling pathway and may be involved in regulation of lipid transport and cholesterol homeostasis [65,66].
The genes FSCB, NAP1L4, and CRISP3 are all expressed in the retina [54] but do not have a clear mechanism by which they may influence spherical equivalent. The product of the FSCB gene is colocalized with a calcium-binding protein, CABYR, in the mouse testis, and was postulated to be involved in the later stages of fibrous sheath assembly [67]. Alterations of the FSCB gene have also been reported in patients with osteosarcoma [68] and in families with cerebral palsy [69]. The NAP1L4 gene encodes a member of the nucleosome assembly protein (NAP) family. This protein is a highly conserved histone chaperone protein that regulates the histone H2A-H2B heterodimer and mediates nucleosome formation [70,71]. This gene is located near the imprinted gene domain of 11p15.5, an important tumor-suppressor gene region that has been linked to growth-related disorders such as Beckwith-Wiedemann syndrome and various human cancers [72]. The CRISP3 gene is widely distributed in tissues with an exocrine function [73] and is suggested to play an role in the pathogenesis of Sjögren's syndrome [74,75] and prostate cancer [76][77][78][79].
To further infer the potential mechanisms of PTCHD2, CRISP3, NAP1L4, FSCB, and AP3B2 gene underlying refractive errors, we identified genes that were predicted with strong evidence (confidence score >0.70) to be functional partners of these five genes from the STRING (Search Tool for the Retrieval of Interacting Genes/Proteins) database [80]. Due to the limits of the exon array, not all pathway members were included among the 11,571 genes in our gene-based analysis. For those tested, we found no association at the significance level of 0.01 (Appendix 6), suggesting that these five genes may contribute to the development of refractive error either independently or through other functional partners that have yet to be identified. As refractive errors often develop due to incorrect morphology (i.e., incorrect axial length), further studies are necessary to understand how variation in these genes may influence ocular growth and development.
Our results contribute to an increasing body of evidence that variants in the GJD2 gene are associated with spherical equivalent. Multiple SNPs in vicinity of the GJD2 gene on chromosome 15q14 have been reported to be associated with refractive error and myopia. rs634990 was first discovered as a susceptibility locus for spherical equivalent in a Dutch population-based study (MAF = 47%, β = -0.25, p = 1.03 × 10 −7 ) [31]. This association was later replicated by the Blue Mountain Eye Study (BMES) [81] and the Consortium for Refractive Error and Myopia (CREAM) [82], as well as in this analysis of the BDES cohort. This variant, rs634990, is also associated with myopia. When we compared the myopes to hyperopes, as was reported by Verhoeven [82], our findings were consistent with this previous study, confirming the association of rs634990 with myopia and indicated a protective effect of its C allele on hyperopia. We did not have data on rs524952, a locus that is strongly correlated with rs634990 (r 2 = 0.90) and associated with spherical equivalent [30,81,82] and myopia [34,83] in European populations. Associations of GJD2 variants with refractive errors and myopia have also been replicated in studies conducted in Asian populations [84,85]. GJD2 encodes a member of the connexin family of proteins (hCx36) that form hexameric intercellular gap junction channels. Expression of the GJD2 gene has been found in the inner and outer plexiform layers of the retina in mice, rats, and humans [86]. Although the function of GJD2 in the human retina remains unknown, animal studies suggest an essential role for gap junction channels containing Cx36 (counterpart protein of hCx36 in mice) in the transmission of visual signals through photoreceptor cells, particularly the rod-cone pathways [87]. Disruption of Cx36 expression in mice can cause visual transmission deficits [88]. However, further investigation is necessary to understand how variation in the GJD2 gene influences spherical equivalent.
Our results support the potential involvement of PRSS56 in spherical equivalent. Associations with spherical equivalent or myopia have been found for two SNPs (rs1656404 and rs1550094) in the PRSS56 gene region. We did not have data on rs1656404, a locus that was first reported for spherical equivalent (MAF = 21.0%, β = -0.151, p = 2.38 × 10 −9 ) in a multiancestry cohort [30] and later for myopia [83]. However, rs1550094 was discovered for myopia development (i.e., age of onset) in the GWAS of 23andMe participants of European ancestry (HR = 1.09 [1.07-1.11], MAF = 30.5%, p = 1.3 × 10 −15 ) [34] and is associated with spherical equivalent and myopia in the present study. Due to the low correlation (r 2 <0.50) between these two SNPs, rs1550094 and rs1656404 were independently associated with a decreased spherical equivalent. Experimental evidence also supports a potential role for PRSS56 in controlling spherical equivalent. The PRSS56 gene encodes a trypsin-like serine peptidase expressed in the human retina, cornea, sclera, and optic nerve [89]. Mutations of the PRSS56 gene cause autosomal recessive posterior microphthalmos and familial nanophthalmos, two rare ocular disorders characterized by abnormally small eyes and severe hyperopia due to reduced axial length [89,90]. In animal studies, mutations of the Mfrp gene caused a progressive accumulation of Prss56 throughout postnatal development and together led to shortening of photoreceptor outer segments and retinal degeneration [91].
We employed an extreme value sampling strategy to maximize the phenotypic differences among individuals. This could lead to the oversampling of variants that are more common among individuals with extreme trait values and thus increase the statistical power of our study. In this study, we directly examined coding variants for association with spherical equivalent. We applied single-variant analysis and gene-based analysis to optimize the power to detect associations for low-frequency and rare variants. However, the exome array was designed using exome sequencing data from ethnically diverse populations. Given that our study population is ethnically homogenous, more than 50% of the genotyped variants were monomorphic, and exclusion of these variants resulted in the mean coverage of six variants in each gene region. This low coverage can be particularly influential when variants in a gene region are mostly singletons or doubletons, resulting in an inflated association driven by a limited number of extreme observations. Therefore, we are cautious with the interpretation of such genes in our genebased analysis. In addition, variants in sequencing data have to reach certain frequencies to be selected into the exome array. Thus, there is a lack of coverage in the exome array of extremely rare variants, which prohibited us from examining these variants for association with quantitative refraction.
Additional studies focused on assessing the role of lowfrequency coding variants across diverse populations are necessary to further validate the reported findings. Such studies would complement the existing GWAS that highlight common variations, and together, these approaches may begin to elucidate the causal alleles and the underlying mechanisms for refractive errors. A better understating of the genetics of spherical equivalent, in combination with the interaction of environmental exposure, may be useful in determining individual risks for myopia and hyperopia and lead to novel or improved measures of prevention, management, and treatment.

APPENDIX 1. SPHERICAL EQUIVALENT MEASUREMENTS IN THE BEAVER DAM EYE STUDY SAMPLES
(A) The average of spherical equivalent in left and right eyes at baseline visit was distributed with a mean of -0.20 D and a standard deviation of 2.9 D (B) The residual spherical equivalent after adjusting for age, sex, years of education, nuclear sclerosis and the first two principal components was distributed with a mean of 0 D and a standard deviation of 2.7 D. (C) Compared to the full BDES cohort, the distribution of spherical equivalent in the analyzed cohort was more heavy-tailed. To access the data, click or select the words "Appendix 1."

APPENDIX 2. BOX PLOTS OF SPHERICAL EQUIVALENT BY GENOTYPES AT RS2297336 AND RS324146 OF TCTE1
The distributions of spherical equivalent residuals were plotted for each observed genotype at rs2297336 (left) and rs324146 (right) of TCTE1 on chromosome 6p21.1. Genotypes are showed on x-axis. The mean spherical equivalent residuals for each genotype is marked by a dark red diamond with the value listed above. The three horizontal bars from the bottom to the top indicate the 25th, 50th and 75th percentile of the spherical equivalent residuals. Observations beyond the 25th and 75th percentile range are represented by black dots. For both SNPs, each copy of its minor allele is associated with decreased spherical equivalent. To access the data, click or select the words "Appendix 2."

APPENDIX 3. DISTRIBUTIONS OF SPHERICAL EQUIVALENT FOR FIVE IDENTIFIED GENES
(Continue) Five gene regions reached the significant threshold in the gene-based analysis (p < 4.3 × 10-6). Spherical equivalent residuals were plotted by the presence of the minor allele of each variant in the gene region. Majority of study samples are homozygotes for major alleles (Red). The remaining groups are individuals with one copy of the minor allele at the variant indicated on the X-axis (Black). The dash line represents the group mean of spherical equivalent residuals. No homozygote or compound heterozygous of the minor allele was seen for these genes. To access the data, click or select the words "Appendix 3."

APPENDIX 4. BOX PLOTS OF SPHERICAL EQUIVALENT BY GENOTYPES AT RS634990 OF GJD2 AND RS1550094 OF PRSS56
The distributions of spherical equivalent residuals were plotted for each observed genotype at rs634990 of GJD2 (left) and rs150094 of PRSS56 (right). After adjusting for age, sex, years of education, nuclear sclerosis and the first two principal components, (A) each copy of the minor allele C of rs634990 is associated with a 0.29 D decrease in spherical equivalent (p = 1.8 × 10 -3 ), and (B) the minor allele G of rs150094 is associated with a 0.33 D decrease in spherical equivalent (p = 1.7 × 10 -3 ). To access the data, click or select the words "Appendix 4."

APPENDIX 5. GENE-BASED ASSOCIATION RESULTS FOR GENES IN LYSOSOME PATHWAYS
a Genes in Lysosome pathway were identified from KEGG (the Kyoto Encyclopedia of Genes and Genomes) database. b Genes reached significant threshold at p < 0.01 were highlighted in bold. To access the data, click or select the words "Appendix 5." a Genes in the STRING database with strong evidence (confidence score > 0.70) of protein-protein interactions with reported genes. The number in the parentheses indicates the proportion of genes tested in our gene-based analysis. b The number in parentheses represents the number of genes out of the total number of genes identified as functional partners that were tested in our gene-based analysis. To access the data, click or select the words "Appendix 6." Articles are provided courtesy of Emory University and the Zhongshan Ophthalmic Center, Sun Yat-sen University, P.R. China. The print version of this article was created on 14 July 2016. This reflects all typographical corrections and errata to the article through that date. Details of any changes may be found in the online version of the article.