Genome-Wide Meta-Analysis of Five Asian Cohorts Identifies as a Susceptibility Locus for Corneal Astigmatism

Summary

Corneal astigmatism refers to refractive abnormalities and irregularities in the curvature of the cornea, and this interferes with light being accurately focused at a single point in the eye. This ametropic condition is highly prevalent, influences visual acuity, and is a highly heritable trait. There is currently a paucity of research in the genetic etiology of corneal astigmatism. Here we report the results from five genome-wide association studies of corneal astigmatism across three Asian populations, with an initial discovery set of 4,254 Chinese and Malay individuals consisting of 2,249 cases and 2,005 controls. Replication was obtained from three surveys comprising of 2,139 Indians, an additional 929 Chinese children, and an independent 397 Chinese family trios. Variants in PDGFRA on chromosome 4q12 (lead SNP: rs7677751, allelic odds ratio = 1.26 (95% CI: 1.16–1.36), P_meta = 7.87×10⁻⁹) were identified to be significantly associated with corneal astigmatism, exhibiting consistent effect sizes across all five cohorts. This highlights the potential role of variants in PDGFRA in the genetic etiology of corneal astigmatism across diverse Asian populations.

Introduction

Astigmatism is a condition where light rays are prevented from focusing at a single point in the eye, resulting in blurred vision at any near or far distance. While astigmatism comprises cornea and non-corneal components, it typically results from the unequal curvature of two principle meridians in the anterior surface of the cornea known as corneal astigmatism [1], [2]. The presence of a high degree of astigmatism during early development is believed to be associated with refractive amblyopia [3], [4], [5], as evidenced by decreased best-corrected visual acuity which cannot be remedied by external corrective lenses. Early abnormal visual input caused by uncorrected astigmatism can lead to orientation-dependent visual deficits, despite optical correction of visual acuity later in life [6]. In addition, it has been suggested that optical blurring by astigmatism may predispose the development of myopia, commonly known as nearsightedness [7], [8], [9], [10].

Astigmatism is highly prevalent across most populations and poses a significant burden to global public health with at least 1 in 3 adults above 30 years of age suffering from astigmatism of −0.5 diopters (D) or less [11]. The reported age-adjusted prevalence of astigmatism was 37.8% for Chinese adults [12], 54.8% in rural Asian Indians [13], 37% (≤−0.75D) for Caucasian in Australia [14] and 36.2% in the US [15]. The prevalence of astigmatism in children varies considerably across different studies and ethnic groups. For instance, the prevalence of astigmatism (≤−0.75D) in school-children ranges from 13.6% in Australia [16], 20% in Northern Ireland [17], 28.4% for Singapore school children [8], to 42.7% for Chinese children in urban China [18].

Although the precise cause of astigmatism is unknown, genetic factors have been implicated in the etiology of corneal astigmatism. Studies have reported a higher risk of developing astigmatism in individuals whose sibling or parents have astigmatism [11]. Evidence from twin studies suggests a genetic etiology in astigmatism development, with the estimated heritability ranging from 30% to 60% [19], [20], [21], [22], [23]. For instance, Hammond and colleagues [21] investigated the inheritance of astigmatism for 226 monzygotic (MZ) and 280 dizygotic (DZ) twins in the United Kingdom and found genetic effects accounted for 42% to 61% of the variation in corneal astigmatism. While most of the twin studies have been conducted in Caucasian populations, a study on Chinese twins in Taiwan reported a heritability estimate of 46% for corneal astigmatism, suggesting that genetic factors account for a similar extent in the etiology of astigmatism for Asian populations [22]. However, no genetic loci have been systematically and consistently identified to be implicated in the development of corneal astigmatism.

Here we report the findings from the meta-analyses of five genome-wide association studies (GWAS) performed in 8,513 individuals from three Asian populations. The discovery phase of our study comprises 4,254 individuals from two population-based GWAS performed in adults of Chinese and Malay ethnicities from the Singapore Prospective Study Program (SP2) and the Singapore Malay Eye Study (SiMES) respectively. The replication phase comprises of data from three other genome-wide surveys of: (i) 2,139 Indian adults from the Singapore Indian Eye Study (SINDI); (ii) 929 Chinese school children from the Singapore Cohort Study of the Risk Factors for Myopia (SCORM); and (iii) 397 Chinese trios of parents and astigmatic offsprings from the Singaporean Chinese in the Strabismus, Amblyopia and Refractive Error Study (STARS).

Results

The characteristics of the post-QC samples from the five studies are summarized in Table 1. The post-QC SP2 dataset comprised 2,016 adults, of which 1,231 individuals had corneal astigmatism (≤−0.75 D) and 785 subjects were defined as non-astigmatic controls. The post-QC SiMES dataset comprised 2,238 adults (1,018 cases and 1,220 controls). In total, 462,518 and 515,712 autosomal genotyped SNPs passed stringent quality control criteria for SP2 and SiMES respectively and the genome-wide meta-analysis was conducted on 460,528 SNPs present in both studies.

There was no evidence of over-inflation of statistical significances due to population structure in either of the discovery cohorts (SP2 λ_GC = 1.006, SiMES λ_GC = 1.007) or in the meta-analysis of both studies (overall λ_GC = 1.007). Suggestive evidence of association (defined as 10⁻⁶<P-value<10⁻⁵) were seen in each of SP2 and SiMES (Figure S1A and S1B), as well as in the meta-analysis of SP2 and SiMES where a collection of SNPs deviated from their expected distributions in the quantile-quantile plots of the P-values (Figure S1C).

None of the SNPs in the discovery meta-analysis attained genome-wide significance of P-value<5×10⁻⁸. Seven SNPs exhibited evidence stronger than P-value<1.0×10⁻⁵ and these were found to cluster in the platelet-derived growth factor receptor alpha (PDGFRA) gene on chromosome 4q12 (lowest P = 9.44×10⁻⁷ at rs7677751; Table 2; Figure S2). Interestingly, these SNPs are located within the MYP9 region identified previously as a candidate locus for myopia through linkage scans [24].

In the replication phase with the three additional GWAS cohorts, three SNPs in PDGFRA (rs7677751, rs2307049 and rs7660560) attained genome-wide significance in the combined analysis (Table 2) with the lead SNP rs7677751 from the discovery phase remaining as the strongest signal in the combined analysis (P = 7.87×10⁻⁹; Figure 1). All seven SNPs from the discovery phase exhibited P-values<0.05 in SINDI but not in SCORM or STARS. However the direction and magnitude of the effect sizes at these seven SNPs in all three replication cohorts were highly similar to those seen in the discovery populations of SP2 and SiMES (Table 2, Figure 2). No significant evidence of effect size heterogeneity was detected across the SNPs (heterogeneity I² P-value≥0.246), and the minor allele frequencies of these SNPs are consistently similar across all five studies (Table S1). A genome-wide meta-analysis of the combined five cohorts did not yield any additional locus with genome-wide significance (see Figure S3 for QQ and Manhattan plots, λ_GC = 1.002; Table S2).

At the lead SNP rs7677751 in PDGFRA, the frequency of the risk T-allele ranged from 0.19 to 0.26 in the five cohorts and conferred a 26% higher risk of corneal astigmatism than the C allele (OR = 1.26, 95% CI = 1.16–1.36) in the meta-analysis across all five studies. This SNP alone explains 0.41% of the variation in corneal cylinder power. In addition, a general genetic model identified that the 5.5% of the individuals in the combined cohorts that carry the TT genotype at rs7677751 had a 1.65-fold (95%CI = 1.33–2.06, P-value = 6.23×10⁻⁶) increase in the risk of developing corneal astigmatism compared to those that are not carrying any copies of the risk allele (Figure S4). All of the associated SNPs spanned 10 kb within PDGFRA at 4q12 (Figure 2), and a high degree of linkage disequilibrium is seen at this locus in all three Asian populations (Chinese, Malays and Asian Indians; Figure S5).

Discussion

We have performed a genome-wide survey for corneal astigmatism across 8,513 individuals, where the discovery phase combines the data from two GWAS performed in Chinese and Malay adults, and the replication phase included Asian Indian adults, Chinese children and family trios. We observed a strong and consistent association with the onset of corneal astigmatism at the PDGFRA gene locus on chromosome 4q12 across all five Asian cohorts, with three SNPs in this locus exhibiting evidence stronger than genome-wide significance in the meta-analysis. To the best of our knowledge, this is the first GWAS to investigate the genetic etiology of corneal astigmatism in a genome-wide fashion.

The PDGFRA gene spans 69 kb with 23 coding exons and resides on chromosome 4q12. The receptor for platelet-derived growth factor (PDGF) contains two types of subunit, a- and β- PDGFRA, which are differentially expressed on the cell surface [25]. PDGFR-a binds to three forms of PDGF (PDGF-AA, AB and BB) and mediates many biological process including embryonic development, angiogenesis, cell proliferation and differentiation. The role of PDGFRA in cellular growth and proliferation is underlined by its contribution to the pathogenesis of gastrointestinal stromal tumours [26]. A large body of evidence has shown that both PDGF and its receptors are expressed in corneal epithelium, stromal fibroblasts and endothelium [27], [28] as well as proliferative retinal tissues in eyes [29], [30], [31]. Along with other cytokines (epidermal growth factor, transforming growth factor-a,-β etc), studies have further suggested that PDGF and its receptors can mediate corneal fibroblast migration, matrix remodeling and play an important role in corneal wound healing [28], [32], [33], [34]. The corneal stroma comprises a large portion of the cornea; the sensitivity of stromal tissue to various growth factors is well described [35]. The administration of PDGF resulted in keratinocyte elongation using rabbit corneal stroma tissue [36]. In light of this, a role for PDGFRA in the regulation of ocular development and parameters cannot be excluded, and our study suggests that genetic polymorphisms within PDGFRA may be involved in the regulation of corneal biometrics resulting in the occurrence of corneal astigmatism.

In addition, Hammond et al. reported that 4q12 (MYP9; LOD 3.3) was significant linked with myopia from a genome-wide linkage study of 221 dizygotic twin pairs [24], and subsequent replication revealed nominal significance of 4q12 (P = 0.065) for refractive error in African-American families [37]. We thus undertook a candidate SNP approach with the identified SNPs to investigate the possible association between PDGFRA and (i) the onset of high myopia; (ii) the refractive error as a quantitative trait. We did not observe any striking association between the identified variants with either outcomes, suggesting that the association of PDGFRA with corneal astigmatism is probably not through any shared etiology with myopia.

The lead SNP in our analyses rs7677751 is located in the intro 1 of PDGFRA. Interestingly, among the SNPs identified, rs2228230 is coding-synonymous (valine:GTC>valine:GTT) and resides in exon 18, while rs3690 is within the untranslated-3′ region. These three SNPs (rs7677751, rs2228230 and rs3690) are strongly correlated with each other (pair-wise Pearson correlation coefficient r ranging from 0.77 to 0.81), although the association evidence at the latter two SNPs did not reach genome-wide significance. As the next closest gene (GSX2) from the 5′ end of PDGFRA is 127 kb away and is not within the LD block with our identified SNPs (Figure 1), it is unlikely that the signals observed in our study are attributed to functional variants located beyond PDGFRA.

Our group recently reported a strong association between variants in PDGFRA with corneal curvature [38]. Corneal curvature is an ocular dimension defined as the average of the radius of corneal curvature at the horizontal and vertical meridians. Myopic eyes have been found to have steeper corneas (reduced radius of curvature), but the significant correlation between corneal curvature and refractive error was not consistently observed [39], [40], [41]. Excessively flatter cornea is associated with cornea plana, producing high hyprotropia and likely resulting in angle-closure glaucoma [42], [43]. On the other hand, corneal astigmatism is an eye-disorder, where the cornea is more curved in one meridian direction compared to the other. This fragmentizes the light rays entering the eye, leading to the inability to focus onto a single point in the eye [1]. It is thus interesting that the same PDGFRA gene has been identified in two ocular outcomes that are biologically different, given the presence of a weak correlation between corneal astigmatism and corneal curvature (Spearman correlation coefficient r between 0.088 and 0.192 in our cohorts; Figure S6), pointing to a possible pleiotropic contribution of PDGFRA.

Our study has adopted a binary definition of corneal astigmatism (affected and unaffected) that is commonly adopted in clinical practice and eye-trait epidemiology [16], [44]. One caveat of this definition is the potential for misclassifying the affected status, particularly for samples with cylinder power around the cutoff threshold of −0.75D. To evaluate the robustness of our findings to the choice of threshold used, we additionally performed the association analysis at the identified SNPs with four different combinations of the thresholds used to define cases and controls. We observed that the odds ratios were highly similar across all four scenarios (Table S3), with the combined evidence at rs7677751 ranging from P_meta of 1.5×10⁻⁴ to 6.7×10⁻⁸. Unsurprisingly, the association evidence was weakest in the scenario with the most stringent thresholding (≤−1.5D for cases and >−0.5 for controls), given this stringency comes at the expense of decreasing the number of individuals in each study. We additionally performed a secondary analysis treating corneal cylinder power as a quantitative trait. Strong statistical evidence was consistently observed at the three leading SNPs (rs7677751, P = 1.76×10⁻⁷; rs2307049, P = 3.41×10⁻⁷ and rs7660560, P = 4.41×10⁻⁷; Table S4), indicating that our findings are robust to the definition of the phenotype.

Owing to the relatively small sample sizes within each of the five GWAS studies, we have chosen to prioritize our survey to identify genetic variants that contribute to the etiology of corneal astigmatism in multiple Asian populations. While Malays have been observed to be genetically closer to the Chinese, the Asian Indians tend to be genetically closer to the Caucasians [45]. Our discovery at PDGFRA thus suggests that part of the underlying biological pathway responsible for astigmatism development is common to multiple populations, although there may be population-specific genetic variants that our current study is not sufficiently powered to identify.

Our study has included two pediatric Chinese populations (SCORM and STARS) with school or pre-school children who are still progressing to their final phenotype. It was documented that a high degree of astigmatism occurs during infancy and the early childhood [46]. The prevalence rates remain stable during young adulthood (18 to 40 years), but increase consistently during late adulthood at aged 40 years or older [1], [12]. Studies have also indicated that the age-related change in astigmatism is associated with meridians changes in the cornea [11]. Children and adolescents have a predominance of “within-the-rule” corneal astigmatism in general, where the vertical curve is greater than the horizontal (axis of 1° to 15°); while in older adults, it shifts to “against-the-rule” astigmatism (axis of 75°–105°) [47], [48]. However, our study considers corneal astigmatism without reference to the axis nor the longitudinal changes from children to adults. Whether PDGFRA plays the same role in pediatric and adults populations will however need further investigation.

Materials and Methods

Ethics statement

This study adhere to the Declaration of Helsinki. Ethics approvals have been obtained from the Institutional Review Boards of the Singapore Eye Research Institute, Singapore General hospital, National University of Singapore and National Healthcare Group, Singapore. In all cohorts, participants provided written, informed consent at the recruitment into the studies. For studies involving children who were still minors (SCORM and STARS), written informed consent was obtained from the children's parents.

Discovery cohorts

Singapore Prospective Study Program (SP2)

Participants included in SP2 were recruited from a revisit of four previously conducted population-based surveys in Singapore: Thyroid and Heart Study 1982–1984 [49], National Health Survey 1992 [50], National University of Singapore Heart Study 1993–1995 [51] and the National Health Survey 1998 [50], [51], [52]. These studies comprise random samplings of individuals stratified by ethnicity from the entire Singapore population. From 2003 to 2007, a total of 10,747 adults aged 40 years or older were invited in this follow-up survey of which 5,157 underwent a health examination and blood samples were drawn. The present genome-wide genotyping involved individuals of Chinese descent only (n = 2,867). Complete post-filtering data on corneal astigmatism for GWAS were available for 2,016 individuals.

Singapore Malay Eye Study (SiMES)

SiMES is a population-based cross-sectional survey on eye diseases for Malay adults aged 40 to 80 years living in Singapore. It was conducted between August of 2004 and June of 2006. The details of the study design and methods have been previously described [53]. In brief, a total of 4,168 Malay residents in the Southwestern area of Singapore were identified and invited for a detailed ocular examination which 3,280 (78.7%) participated. Genome-wide genotyping was performed in 3,072 individuals [54], [55]. Individuals with cataract surgery and missing corneal astigmatism measurements were excluded from the study. Complete post-filtering data on corneal astigmatism for GWAS were available for 2,238 individuals.

Replication cohorts

Singapore Indian Eye Study (SINDI)

SINDI is a population-based survey of major eye diseases [56] in ethnic Indians aged 40 to 80 years living in the South-Western part of Singapore and was conducted from August 2007 to December 2009. In brief, 4,497 Indian adults were eligible and 3,400 participated. Genome-wide genotyping was performed in 2,953 individuals [54], [55]. As in the discovery cohorts, participants were excluded from the study if they had cataract surgery and missing corneal astigmatism data. Complete post-filtering data on corneal astigmatism were available for 2,134 participants.

Singapore Cohort Study of the Risk Factors for Myopia (SCORM)

A total of 1,979 children in grades 1, 2, and 3 from three schools were recruited from 1999 to 2001 with detailed information described elsewhere [57]. The children were examined on the school premises annually by a team of eye care professionals. The GWAS was conducted in a subset of Chinese children of 1,116 subjects [58]. The phenotype used in this study was based on the corneal astigmatism measured on the 4^th annual examination of the study (children at age 10 to 12 years). Complete post-filtering data on corneal astigmatism measurements and SNP data were available in 929 SCORM children.

Singaporean Chinese in the Strabismus, Amblyopia, and Refractive Error Study (STARS)

STARS is a population-based survey of Chinese families with children aged 6 to 72 months residing in the south-western and western region of Singapore. Disproportionate random sampling by 6-month age groups resulted in the recruitment and subsequent eye examination of 3,009 Chinese children between May 2006 and November 2008. Details of the study design and methodology have been previously described [59]. A total of 1,451 samples from 440 nuclear families were included for genome-wide genotyping. In all, 397 trio-sets of parents and their offsprings with corneal astigmatism had complete post-filtered genotype data.

Measurements and definition of corneal astigmatism

All studies used a similar protocol to measure ocular phenotypes including corneal curvature, autorefraction and cylinder power by a team of eye care professionals. Participants in SP2, SIMES and SINDI underwent non-cycloplegic automated refractive assessments using the autorefractor (Canon RK-5, Tokyo, Japan). For SCORM and STARS, cycloplegic measurements (Canon RK-F1, Tokyo, Japan) were performed 30 minutes after three drops of 1% cyclopentolate which were administered 5 minutes apart.

Corneal curvature radii in the horizontal and vertical meridian were determined with keratometry in millimeters [60]. The keratometer measured the anterior corneal surface and used a refractive index of 1.3375 to account for the contribution from the posterior corneal surface to derive the corneal refractive power in diopters. Corneal cylinder power was calculated as the difference between the flattest and steepest meridian of the keratometry readings in diopters of power.

As the corneal cylinder power between the right and left eyes are strongly correlated across all five cohorts (Pearson's correlation coefficient r ranging from 0.51 to 0.79; P<2.2×10⁻¹⁶), the mean corneal cylinder power over both eyes was used to define corneal astigmatism. Averaging ocular measurements between two eyes in genetic studies has been suggested to be statistically more powerful than using the information from only one eye [61], and this approach has been consistently adopted in genome-wide studies of myopia [62], [63]. As with previous studies [16], [44], we have defined individuals with average corneal cylinder power ≤−0.75D as cases, and those with average corneal power between −0.75D and 0D as controls.

Genotyping and data quality control

For SP2, a total of 2,867 blood-derived DNA samples were genotyped using the Illumina Human610 Quad and 1Mduov3 Beadchips. For the samples that were genotyped on the two platforms, the genotypes from the denser platform were used in our study. For SiMES (n = 3,072), SINDI (n = 2,593) and STARS (n = 1,451), the Illumina Human610 Quad Beadchips was used for genotyping all DNA samples. For the 1,116 SCORM children, DNA samples were genotyped on the Illumina HumanHap 550 Duo Beadchips.

Detailed data quality control (QC) procedures for each study were provided in the supplementary information (Text S1). In brief, for case-control study design, QC criteria included a first round for autosomes SNP QC to obtain a cleaned set of genotypes for sample QC, by excluding SNPs with: (i) missingness (per-SNP call rate) >5%; (ii) minor allele frequency (MAF) <1%; and (iii) HWE p-value<10⁻⁷. Using the subset of SNPs passing the first round QC, samples were then excluded based on the following conditions: (i) per-sample call rates of less than 95%; (ii) excessive heterozygosity (defined as the sample heterozygosity to be beyond 3 standard deviations from the mean sample heterozygosity); (iii) cryptic relatedness; (iv) gender discrepancies; and (v) deviation in population membership from population structure analysis. A second round of SNP QC was then applied to the remaining samples passing quality checks. We excluded SNPs with missingness >5%, gross departure from HWE (P value<10⁻⁶), MAF<1% and low concordance between duplicate samples on different genotype platforms (relevant to SP2 samples only).

Population structure was ascertained using principal components analyses (PCA) with the EIGENSTRAT program [64]. Population substructure of SP2 and SiMES was examined by PCA with respect to three population panels in the HapMap samples (Figure S7). Due to the presence of population structure within the Malay and Indian samples (Figures S8 and S9 respectively), we adjusted for the top 5 principal components in the association analyses for the SiMES and SINDI datasets.

For the STARS trios, we additionally excluded samples and trio-sets on the basis of excessive Mendelian inconsistencies defined as having >1% of the post-QC SNPs exhibiting Mendelian errors. SNPs with more than 10% Mendelian errors are excluded from the association analyses, and the genotypes leading to Mendelian errors in all other remaining SNPs are coded as missing. As family trios are more robust to the presence of population structure, we did not exclude any samples due to population structure.

Statistical analysis

The genome-wide association tests were performed using PLlNK (version 1.07; http://pngu.mgh.harvard.edu/~purcell/plink/) as the primary analysis tool. A logistic regression adjusted for age and gender is used to model the association of genetic markers with corneal astigmatism. For each of SiMES and SINDI, the top 5 principal components of genetic ancestry from the EIGENSTART PCA were also included as covariates to adjust for population stratification in these populations. We assumed an additive genetic model where the genotypes of each SNP is coded as 0, 1, and 2 for the number of minor alleles carried, in keeping with increments in allelic dosage. For family GWAS association tests in STARS, a transmission disequilibrium test (TDT) is used to measure significant distortions in transmission of an allele from heterozygous parents to the affected offspring under the condition of Mendel's law [65].

We also performed a quantitative trait analysis with the average corneal cylinder power as the outcome. This is performed in PLINK for the unrelated samples, and in FBAT (http://www.hsph.harvard.edu/~fbat/fbat.htm) for the family trios. As the distribution of the quantitative trait of corneal cylinder power is skewed (Figure S10), we performed a normal quantile transformation [66] prior to the association analysis for unrelated samples. For family-based data, no transformation was conducted since the FBAT does not require normal trait [67]

Meta-analyses are performed using weighted-inverse variance estimated from fixed-effect modeling in METAL (http://www.sph.umich.edu/csg/abecasis/metal/). We adopt the method by Kazeem and Farrall [65] to pool the evidence from the case-control analyses and the family trio TDT. For the quantitative trait analysis, the overall z statistics is calculated as a weighted sum of the z-statistics from the linear regressions in the non-familial data and FBAT analysis for the family-based data, weighted by number of unrelated individuals or trios in the respective studies [68].

Results from a genome-wide meta-analysis of the SNPs common to SP2 and SiMES are used in the discovery phase to identify putative variants that are associated with the onset of corneal astigmatism, defined as a P-value<10⁻⁵. The remaining three cohorts (SINDI, SCORM and STARS) are used to validate the putative findings. In addition, a genome-wide meta-analysis of all five datasets is also performed. Genotyping quality of all reported SNPs in this paper have been visually assessed by checking the intensity clusterplots.

Supporting Information

Attachment 1

Attachment 2

Attachment 3

Attachment 4

Attachment 5

Attachment 6

Attachment 7

Attachment 8

Attachment 9

Attachment 10

Attachment 11

Attachment 12

Attachment 13

Attachment 14

Attachment 15