A genome-wide association study of corneal astigmatism: The CREAM Consortium.

Purpose: To identify genes and genetic markers associated with corneal astigmatism.
Methods: A meta-analysis of genome-wide association studies (GWASs) of corneal astigmatism undertaken for 14 European ancestry (n=22,250) and 8 Asian ancestry (n=9,120) cohorts was performed by the Consortium for Refractive Error and Myopia. Cases were defined as having >0.75 diopters of corneal astigmatism. Subsequent gene-based and gene-set analyses of the meta-analyzed results of European ancestry cohorts were performed using VEGAS2 and MAGMA software. Additionally, estimates of single nucleotide polymorphism (SNP)-based heritability for corneal and refractive astigmatism and the spherical equivalent were calculated for Europeans using LD score regression.
Results: The meta-analysis of all cohorts identified a genome-wide significant locus near the platelet-derived growth factor receptor alpha (PDGFRA) gene: top SNP: rs7673984, odds ratio=1.12 (95% CI:1.08-1.16), p=5.55×10-9. No other genome-wide significant loci were identified in the combined analysis or European/Asian ancestry-specific analyses. Gene-based analysis identified three novel candidate genes for corneal astigmatism in Europeans-claudin-7 (CLDN7), acid phosphatase 2, lysosomal (ACP2), and TNF alpha-induced protein 8 like 3 (TNFAIP8L3). Conclusions: In addition to replicating a previously identified genome-wide significant locus for corneal astigmatism near the PDGFRA gene, gene-based analysis identified three novel candidate genes, CLDN7, ACP2, and TNFAIP8L3, that warrant further investigation to understand their role in the pathogenesis of corneal astigmatism. The much lower number of genetic variants and genes demonstrating an association with corneal astigmatism compared to published spherical equivalent GWAS analyses suggest a greater influence of rare genetic variants, non-additive genetic effects, or environmental factors in the development of astigmatism.

Introduction

Astigmatism is a commonly occurring refractive error that leads to impaired visual acuity if uncorrected and is a risk factor for amblyopia [1-4]. The two major sources of refractive astigmatism in the human eye are the cornea and the crystalline lens. In emmetropic eyes, a low degree of with-the-rule (WTR) corneal astigmatism is typically compensated by a low degree of against-the-rule (ATR) lenticular astigmatism [5]. For individuals with higher levels of refractive astigmatism, corneal astigmatism is usually the major contributor, while lenticular astigmatism is within the normal range [6].

Studies in chicks have recently shown that the eye can compensate for experimentally induced astigmatism through the alteration of corneal curvature [7]. This suggests that the reduction in innate astigmatism seen during infancy in children occurs via active emmetropization [8]. Potential reasons why astigmatism still arises despite the presence of an emmetropization system include (a) astigmatism of too high a degree to be compensated within the juvenile period, (b) astigmatism outside the “operating range” of the emmetropization system, for example, producing a retinal image that is not detected as being caused by astigmatism or that arises at an age beyond that at which emmetropization normally acts, and (c) a failure of the emmetropization response [2].

Several lines of evidence support the role of genetics in the etiology of astigmatism. First, epidemiology studies have shown marked differences in the prevalence of astigmatism across ethnic groups, even after accounting for differences in spherical refractive error. For instance, 78% of native American Tohono O’odham children aged 0–8 have at least 1 diopter (D) of corneal astigmatism, and in Australian children aged 12, at least 1 D of corneal astigmatism was found in 19% of European individuals versus 50% of East Asian individuals [9,10]. Second, corneal and refractive astigmatism have been reported as being moderately/highly heritable (heritability of 0.3–0.6) in twin studies [11,12]. Third, a genetic segregation study in families with high-degree astigmatism found evidence of Mendelian inheritance [13]. Finally, genetic association studies have identified specific genetic variants associated with susceptibility to either refractive and/or corneal astigmatism [14-17]. Despite these latter studies, our understanding of the genetic contribution to astigmatism has lagged behind that of spherical refractive errors, for which dozens of genetic variants have been discovered [18-21].

Previously, the Consortium for Refractive Error and Myopia (CREAM) reported a genome-wide association study (GWAS) of refractive astigmatism that examined approximately two million genetic markers in 45,931 individuals [17]. Only a single marker reached genome-wide significance (rs1401327 in the NRXN1 gene, p=3.92E−8). Reasoning that the paucity of genome-wide significant hits in the previous CREAM study may have been due to phenotypic uncertainty when studying refractive astigmatism that arose from the combination of both corneal and lenticular influence, CREAM has now undertaken a GWAS of corneal astigmatism. Fan et al. [16] performed a GWAS of corneal astigmatism using a discovery sample of 4,254 East Asian individuals and identified a genome-wide significant locus near the PDGFRA gene. In view of the success of the Fan et al. [16] study, the current analysis has adopted the same phenotype definition.

Methods

The research study followed an analysis plan that was agreed upon by members of CREAM before starting work. This plan was designed to standardize methods across participating CREAM groups and to set timelines for the completion of specific tasks. All research groups known to CREAM with relevant genotype and phenotype data were invited to contribute to the study. Ethical approval for the study was obtained locally for each CREAM study group, and participants gave informed consent. The research was carried out in accordance with the tenets of the Declaration of Helsinki.

Study sample

The demographics of the participating study groups are shown in Table 1. The participants comprised 22,250 European individuals from 14 studies and 9,120 Asian individuals from 8 studies. There were 5,470 European participants and 947 Asian participants aged <25 years.

Table 1

Subject demographics of participating CREAM study groups.

Study	Ancestry	N (cases/controls)	%Female	Age (years)	Corneal astigmatism (D)
				Mean (SD)	Median (IQR)	Range
ALSPAC	European	2279(985/1294)	53.1	15.5 (0.3)	0.683 (0.469–0.959)	0.000–5.680
BMES	European	1238(720/518)	41.8	73.3 (7.6)	0.863 (0.565–1.295)	0.155–8.615
EPIC	European	857(456/401)	58.5	68.7 (7.4)	0.780 (0.527–1.152)	0.075–3.997
FITSA	European	127(62/65)	100	67.9 (3.1)	0.733 (0.530–1.033)	0.270–2.020
GenerationR	European	2071(981/1090)	49.9	6.09 (0.4)	0.725 (0.480–0.995)	0.000–3.370
GHS 11	European	2398(1003/1395)	48.7	55.9 (10.8)	0.65 (0.400–0.950)	0.050–4.350
GHS 21	European	851(383/468)	50.9	55.1 (10.8)	0.65 (0.450–1.000)	0.050–3.800
RAINE	European	1028(407/621)	50.9	20.0 (0.4)	0.649 (0.445–0.905)	0.280–2.440
Rotterdam-I	European	5537(2064/3473)	59.3	69.5 (9.2)	0.601 (0.334–1.007)	0.000–9.663
Rotterdam-II	European	1982(633/1349)	53.8	64.8 (8.0)	0.539 (0.294–0.884)	0.000–6.789
Rotterdam-III	European	2925(1180/1745)	56.2	57 (6.9)	0.618 (0.356–1.019)	0.000–4.869
OGP-A2	European	92(37/55)	44.6	16.0 (4.5)	0.682 (0.512–0.942)	0.185–3.070
OGP-B2	European	446(181/265)	43.7	50.6 (15.4)	0.650 (0.430–0.970)	0.130–4.240
TwinsUK	European	419(201/218)	92.7	64 (10.5)	0.729 (0.476–1.105)	0.000–5.432
BES-610K3	Asian	553 (240/313)	65.6	62.1 (8.4)	0.666 (0.407–1.056)	0.000–3.620
BES-OmniE3	Asian	469 (208/261)	60.1	64.7 (9.5)	0.676 (0.429–1.016)	0.000–5.082
SCES-610K3	Asian	1745 (787/958)	48.7	57.6 (9.0)	0.703 (0.476–1.060)	0.109–5.868
SCES-OmniE3	Asian	545 (257/288)	48.6	59.2 (8.8)	0.723 (0.470–1.065)	0.117–5.404
SCORM	Asian	947 (768/179)	48.6	10.8 (0.8)	1.205 (0.851–1.624)	0.138–3.911
SIMES	Asian	1778 (750/1028)	51.7	59.5 (10.8)	0.662 (0.432–1.016)	0.078–5.618
SINDI	Asian	2261 (814/1447)	48.6	56.5 (9.1)	0.614 (0.411–0.912)	0.115–4.727
STARS	Asian	822 (525/297)	50.0	38.5 (5.3)	1.000 (0.625–1.380)	0.125–3.875

1Association tests were undertaken separately for samples recruited in different waves. 2Association tests were undertaken separately for different age strata (stratum A, age >3 and <25 years; stratum B, age ≥25 years). 3Association tests were undertaken separately for samples genotyped on different platforms.

Phenotype definition

Following Fan et al. [16], cases were defined as participants with corneal astigmatism >0.75 D, and controls were defined as those with corneal astigmatism ≤0.75 D. Corneal astigmatism was averaged between the two eyes, except for participants with data available for only one eye. For the conversion of keratometry readings in millimeters to diopters, we used a conversion factor of 332 divided by the K-reading in mm [22].

Phenotyping, genotyping, and genetic imputation

Anterior corneal curvature was measured using keratometry (the keratometer used by each CREAM study group is listed in Appendix 1), and corneal astigmatism was calculated as the difference in curvature between the steepest and flattest meridians. Participants known to have keratoconus, corneal scarring, ocular surgery, or any corneal/ocular condition that would impair keratometry were excluded from the analysis. DNA samples were extracted from blood or saliva and genotyped on a high-density single nucleotide polymorphism (SNP) platform, as previously described [17]. Each CREAM study group imputed non-genotyped markers from an ancestry-matched reference panel from the 1000 Genomes Project [17] using IMPUTE2 [23] or Minimac [24]. Quality-control filtering was performed in accordance with standard GWAS practices [25]. In general, markers with per-study missingness <0.95, minor allele frequency (MAF) <0.05, or a Hardy–Weinberg disequilibrium p value <1×10−6 were excluded, along with samples with per-study missingness <0.95, extreme heterozygosity, sex mismatch, unaccounted for relatedness, or outlying ancestry [25]. Poorly imputed markers (IMPUTE2 info ≤0.5 or Minimac Rsq ≤0.5) were also excluded.

Genome-wide association studies and meta-analyses

Tests of association between corneal astigmatism case/control status and SNP genotype were performed genome-wide by each participating CREAM study group. The analysis was performed using PLINK [26] for marker genotypes coded 0, 1, or 2 or using mach2dat [24] or ProbABEL [27] for marker genotypes coded as imputed dosage on the scale 0–2. Age and sex were included as a continuous and a binary covariate, respectively. The first five major principal components were also included as continuous covariates if there was evidence of population stratification from Q-Q plots or the genomic control inflation factor (λGC). For samples of related individuals, the analysis method took account of genetic background by treating this as a random effect in the analysis model. Tests of association were conducted separately for participants of European ancestry and participants of Asian ancestry and for younger (age >3 and <25) and older (age ≥25) participants.

Summary statistics from the participating CREAM study groups were submitted to a central site for meta-analysis. Using the approach implemented in easyQC [28], the summary statistics were evaluated by examining quality control plots and metrics, including effect allele frequency (EAF) plots, p value versus z-score (P-Z) plots, standard error versus sample size (SE-N) plots, effect size (odds ratio) distributions, and genomic control inflation factors. Queries were resolved by discussion with study groups analysts, and, where indicated, imputation or association testing was repeated.

Meta-analyses were performed separately for the four demographic strata—younger/older, European/Asian ancestry individuals. Fixed effects, standard error-weighted meta-analysis [29] was performed initially, followed by a random effects meta-analysis [30] for highly associated markers showing excessive between-study heterogeneity of I2>0.5, where I2 is a measure of heterogeneity derived from Cochran’s Q statistic [31]. A p value of 5×10−8 was adopted for declaring genome-wide significant association in the GWAS meta-analyses [32]. Regional association plots were created using LocusZoom [33]. Conditional analysis was performed on GWAS meta-analysis summary statistics using GCTA-COJO [34].

Gene-based tests and pathway analysis

Two gene-based tests, VEGAS2 [35] and MAGMA [36], were used to explore whether specific genes were enriched with strongly associated variants in the GWAS meta-analysis of older European individuals. Attention was restricted to the older European samples because gene-based testing relies on consistent patterns of linkage disequilibrium (LD) across genes, and the sample size was larger for the European meta-analysis compared to that for the Asian cohorts. Markers within 50 kb upstream and downstream of a gene were included in the gene-based tests, with the aim of detecting variants that altered the expression level of genes. The gene-based testing using MAGMA was repeated using an extended flanking region of 200 kb upstream and downstream of genes.

VEGAS2 [35] uses a fast approximation of a permutation-based test to determine whether genes are enriched for highly associated markers and makes use of LD information from an ancestry-matched reference panel to account for association signals shared by markers in LD. The test was implemented to analyze all markers in each gene. MAGMA [36] overcomes the low statistical power inherent when a gene contains many markers, some of which may be in strong LD, by first carrying out a principal components analysis (PCA) for the markers in each gene and then carrying out a per-gene linear regression analysis using the PCA eigenvectors as predictor variables. High statistical power is attained by limiting the regression to the major eigenvalues. Permutation-based p values are calculated to account for the use of a binary outcome as the dependent variable in the linear regression analysis [36].

Gene-set “pathway analysis” was also performed using MAGMA [36]. This was performed using a competitive approach whereby the test statistics for all genes within a gene set were combined to form a joint association statistic. This statistic was compared against that for all other genes not in that set while accounting for the number of SNPs within each gene, gene density, and differential sample size (unequal sample size contributing to each gene) [36]. Gene sets were defined using the Molecular Signatures Database (MSigDB) [37]. Gene definitions and their respective association signals for genes contributing to gene sets were taken from the MAGMA gene-based analyses with the aim of identifying potential biologic processes that may be influenced by these variants.

Shared genetic contribution to traits

LD score regression [38,39] was used to quantify the degree of shared genetic contribution between corneal astigmatism and two related traits, refractive astigmatism and mean spherical equivalent refractive error. GWAS summary statistics for refractive astigmatism and for spherical equivalent refractive error were obtained from previous CREAM studies [17,18]. LD score regression utilizes LD information from an ancestry-matched reference panel and requires large sample sizes; therefore, analyses were limited to European GWAS samples. Specifically, LD score regression was performed using the LDSC program [38,39] for variants present in the HapMap3 CEU reference panel with MAF ≥0.05. The prevalence of corneal astigmatism (defined as an amount >0.75 D) in the general population was taken as 42%, which was calculated as the average for the European ancestry population-based studies contributing to this meta-analysis. LD score regression remains valid when two traits are measured in overlapping samples [39], which was the case for these CREAM GWAS samples.

Results

Meta-analysis of genome-wide association studies

Meta-analyses were performed using a fixed effects model for approximately six million genetic variants (approximately 5,500,000 SNPs and 380,000 indels) in each of the four ancestry/age strata (younger European, older European, younger Asian, and older Asian). However, none of the markers had a p value below the pre-determined threshold of 5×10−8 used to declare genome-wide significance (Appendix 2, Appendix 3, Appendix 4, Appendix 5, Appendix 6, and Appendix 7). Therefore, to increase power, a meta-analysis was performed using data for all four ancestry/age strata, under the assumption that the genetic determinants of corneal astigmatism are consistent across ancestry groups and lifespan. This yielded 49 markers with p values <5×10−8, all of which were located in a narrow interval on chromosome 4 close to the PDGFRA gene (Figure 1, Figure 2, and Figure 3). This locus has previously been identified in GWAS analyses for corneal astigmatism [16], refractive astigmatism [17], and corneal curvature [15,40,41]. Table 2 lists the most strongly associated marker in each region showing suggestive association, defined as a region with at least one marker with p<1×10−5. Both the European and Asian meta-analyses contributed to the association signal at the PDGFRA locus; the most strongly associated marker, rs7673984, had an effect size (odds ratio) of OR=1.15 (95% CI:1.07–1.24; p=1.76×10−4) in Asians, OR=1.11 (95% CI:1.06–1.16; p=5.64×10−6) in Europeans, and OR=1.12 (95% CI:1.08–1.16; p=5.55×10−9) in the meta-analysis of Asians and Europeans. The association of rs7673984 in the individual cohorts examined is summarized in Figure 4. Conditional analysis using GCTA-COJO yielded no additional association signals at the PDGRFA locus independent of rs7673984.

Figure 1

Manhattan plot showing most strongly associated markers in the GWAS fixed-effects meta-analysis for European and Asian participants of all ages combined (n=31,375). Red line: p=5×10−8, blue line: p=1×10−5.

Figure 2

Q-Q plot for the GWAS fixed-effects meta-analysis for European and Asian participants of all ages combined (n=31,375).

Figure 3

Region plot for the most strongly associated region in the GWAS fixed-effects meta-analysis for European and Asian participants of all ages combined (n=31,375).

Table 2

Most strongly associated marker in each region in the GWAS meta-analysis of all samples (Europeans and Asians of all ages combined).

SNP	Chr	Pos	Effect allele	Other allele	EAF	OR (95%CI)	P value	Nearest gene
rs7673984	4	55,088,761	T	C	0.22	1.12 (1.08–1.16)	5.55×10−9	PDGFRA
rs34751092	4	24,129,037	A	G	0.28	1.09 (1.05–1.13)	6.07×10−7	PPARGC1A
rs630203	5	141,444,269	T	G	0.74	0.92 (0.88–0.95)	8.83×10−7	MRPL11P2
rs75607298	8	128,611,496	A	G	0.72	1.14 (1.08–1.21)	2.28×10−6	CASC11
rs62401199	6	43,813,341	T	C	0.14	1.15 (1.09–1.22)	3.29×10−6	LINC01512
rs753992	11	47,349,846	A	G	0.29	0.91 (0.87–0.95)	3.48×10−6	MADD
rs3214101	11	114,009,408	A	T	0.68	1.08 (1.05–1.12)	3.75×10−6	ZBTB16
rs10985068	9	123,629,724	C	G	0.12	1.13 (1.07–1.19)	5.87×10−6	PHF19
rs62128379	2	26,960,055	T	C	0.85	1.13 (1.07–1.19)	6.00×10−6	KCNK3
rs9939114	16	84,023,972	A	G	0.05	0.54 (0.41–0.70)	6.01×10−6	NECAB2
rs60083876	7	34,228,819	A	T	0.95	1.31 (1.17–1.48)	6.53×10−6	BMPER
rs859362	1	175,495,090	T	C	0.19	1.10 (1.05–1.15)	7.14×10−6	TNR
rs11775037	8	108,317,615	A	G	0.20	1.10 (1.05–1.14)	7.31×10−6	ANGPT1
rs7036824	9	96,149,894	T	C	0.94	0.83 (0.77–0.90)	8.78×10−6	C9orf129
rs142168171	7	71,253,651	I	R	0.09	1.37 (1.19–1.57)	9.14×10−6	CALN1
rs7278671	21	41,047,876	A	G	0.51	0.93 (0.90–0.96)	9.31×10−6	B3GALT5
rs191640722	1	119,264,997	C	G	0.09	0.87 (0.82–0.93)	9.53×10−6	LOC100421281
rs36107906	2	44,162,800	D	R	0.29	1.09 (1.05–1.13)	9.61×10−6	LRPPRC
rs4896367	6	138,807,281	T	C	0.72	1.09 (1.05–1.14)	9.75×10−6	NHSL1
rs35587414	1	153,174,958	T	C	0.15	1.13 (1.07–1.20)	9.79×10−6	LELP1

EAF=effect allele frequency, OR=odds ratio.

Figure 4

Forest plot and summary table for lead variant rs7673984 across all cohorts. Studies listed above the dotted line are new cohorts not included in the only prior GWAS for corneal astigmatism [16]. EAF=effect allele frequency. (Note that rs7673984 was excluded from the Rotterdam-I cohort analysis during quality control filtering).

Gene-based analyses

To explore whether specific genes were enriched for markers with low p values in the GWAS meta-analysis, we performed gene-based tests using VEGAS2 [35] and MAGMA [36]. These programs use different approaches to test for such enrichment (see Methods). Due to the requirement for an ancestry-matched reference panel, analyses were conducted using the results of a GWAS meta-analysis of European samples of all ages (however, similar results were obtained when attention was restricted to the meta-analysis of older Europeans). The 10 most strongly associated genes from the VEGAS2 and MAGMA analyses are shown in Appendix 8 and Appendix 9. There was a high degree of overlap between the results of the two programs, with the genes ACP2, CLDN7, ELP5, and CTDNEP1 showing the strongest association in both analyses (Appendix 8 and Appendix 9). In the MAGMA gene-based test, these four genes and TNFAIP8L3 achieved p<0.05 after stringent Bonferroni correction, whereas this was not the case for VEGAS (Appendix 8 and Appendix 9). A further exploratory gene-based analysis that included markers up to 200 kb upstream or downstream of each gene—an approach that has been successful for certain traits [42]—failed to identify any additional genes associated with corneal astigmatism.

Pathway analysis

As biologic processes tend to involve multiple genes, a gene-set analysis was performed with MAGMA [36] using the gene-based analysis results for the European samples. Gene-set analyses seek to identify potential biologic mechanisms enriched for genes with markers attaining low p values in the GWAS meta-analysis. However, no gene sets were identified as demonstrating a greater level of association with corneal astigmatism than would be expected by chance (when flanking regions of either ±50 kb or ±200 kb upstream or downstream of genes were tested).

SNP heritability and genetic correlation between traits

LD score regression was used to quantify the heritability explained by commonly occurring genetic variants (“SNP heritability”) and the degree of genetic sharing between corneal astigmatism and two related traits, refractive astigmatism and spherical equivalent refractive error (Table 3 and Table 4). The SNP heritability (h2) estimates for corneal and refractive astigmatism (~5% and ~1%, respectively) were lower than for the spherical equivalent (~23%); indeed, the SNP heritability estimates for corneal and refractive astigmatism were not significantly different from zero. The genetic correlation estimates also had high standard errors and therefore yielded very imprecise estimates (Table 4). These hinted at a high genetic correlation between corneal astigmatism and the spherical equivalent; however, in view of the low SNP heritability estimate for the astigmatism traits, these findings imply that much larger sample sizes and/or a more homogeneous population sample is needed to obtain robust findings.

Table 3

SNP-heritability estimated using LD Score Regression (Europeans only).

Trait	No. of Markers	SNP-heritability (SE)	P value
Corneal Astigmatism	1,024,525	0.0555 (0.0381)	0.15
Refractive Astigmatism	1,056,658	0.0136 (0.0218)	0.53
Mean Spherical Equivalent	1,056,658	0.2326 (0.0175)	2.60×10−40

Table 4

Genetic correlations between pairs of refractive error traits in samples of European ancestry from the CREAM consortium (using LD Score Regression).

Trait Pairs	No. of Markers	Genetic Correlation (SE)	P value
RA and CA	934,512	0.2327 (0.703)	0.7406
MSE and CA	1,024,525	−0.0238 (0.1599)	0.8815
RA and MSE	1,056,658	0.7732 (0.6504)	0.2345

RA=refractive astigmatism, CA=corneal astigmatism, MSE=mean spherical equivalent. P values refer to likelihood of non-zero correlation between traits.

Discussion

This GWAS for corneal astigmatism in a combined sample of Europeans and Asians identified a single genome-wide significant locus in the promoter region of the PDGFRA gene, replicating the previous discovery of this corneal astigmatism locus by Fan et al. [16] in a predominantly Asian sample. Therefore, despite a fourfold increase in sample size (n=31,370 versus n=7,719) compared to the only previous GWAS meta-analysis for corneal astigmatism [16], the standard, single-marker GWAS analysis performed here did not identify any new loci. GWAS analyses for spherical equivalent and other morphological traits in equivalently sized samples have identified dozens of independent risk loci [18,19]. This paucity of GWAS loci for corneal astigmatism mirrors that observed in a previous large-scale GWAS for refractive astigmatism [17]. Our LD score regression-based SNP heritability estimates for corneal astigmatism (h2 ~5%) and refractive astigmatism (h2 ~1%)—the first ever estimates for these traits—were also much lower than those for the spherical equivalent (h2 ~23%), suggesting that common, additively acting SNPs make a relatively minor contribution to the development of astigmatism. In the study by Fan et al. [16] that originally identified the association between SNPs close to the PDGFRA gene and corneal astigmatism, the authors speculated that the underlying causal mechanism was common to populations of diverse ancestry and not specifically to those of Asian origin. This was based on the knowledge that their GWAS included individuals of Indian ancestry, who are more closely genetically related to Europeans than East Asians [16]. Our findings support this theory.

The association between PDGFRA SNPs and corneal astigmatism has been replicated in a previous study of Europeans (n=1968) but not in another smaller study (n=1013) [41,43]. Variants at this locus were not associated with refractive astigmatism in GWAS meta-analyses of n=45,287 participants [17] yet were associated with corneal curvature in an Asian sample [15] and with both corneal curvature and axial length (but not refractive error) in a European sample [41]. This complex series of findings suggests a role for PDGRFA in the regulation of eye size and corneal astigmatism; however, the underlying mechanism of action remains uncertain.

In contrast to the single-marker analyses, gene-based analysis did provide new insight into the genetic basis of corneal astigmatism, implicating the genes ACP2, CLDN7, CTDNEP1, ELP5, and TNFAIP8L3. Three of these five genes—CLDN7, CTDNEP1, and ELP5—are tightly clustered on chromosome 11, with their respective (gene-based) association signals sharing many variants in common. Therefore, a parsimonious interpretation is that only one of the genes has a causal association with astigmatism and that the other two genes are false-positive associations detected due to the signal from the causal gene. Of the three genes, CLDN7, which encodes the claudin-7 membrane protein [44], appears to be the most biologically plausible candidate. Claudins are responsible for tight junction formation and function [45], with claudin-7 being the subtype present in human corneal epithelium and endothelium [46]. Currently, how claudin-7 may contribute to the development of corneal astigmatism is unclear. The acid phosphatase 2, lysosomal gene (ACP2) is located on chromosome 17 and codes for the beta subunit of the degradative enzyme, lysosomal acid phosphatase (LAP). Interestingly, LAP activity is enhanced in keratoconic corneas [47,48]. The TNFAIP8L3 gene located on chromosome 15 codes for TNF alpha-induced protein 8 like 3, which is preferentially expressed in secretory epithelial cells [49]. TNFAIP8L3 is implicated as a negative regulator of inflammation (and carcinogenesis) through its role in TNFα and phospholipid signaling. Based on this evidence, the CLDN7, ACP2, and TNFAIP8L3 genes are promising susceptibility genes for corneal astigmatism. It is important to note that while the statistical support for the above three genes was much stronger in the MAGMA analysis than in the VEGAS2 analysis, the two software programs similarly ranked the most strongly associated genes. This commonality between the MAGMA and VEGAS2 results provides greater confidence that the findings are robust than would be the case for findings identified using either software program alone, as the statistical models and hypothesis tests used by the two programs differ, especially regarding the adjustment for variants in LD.

The strengths of this investigation are that data from multiple population samples were combined and meta-analyzed and that gene-based and pathway-based follow-up analyses were undertaken to leverage new biologic insights into the genetics of astigmatism. The weaknesses were that although the samples included both European and Asian ancestry individuals, trans-ethnic meta-analysis [50] was not performed due to the small size of the Asian sample compared to the European sample and that the age spectrum of the participants was very broad. The latter point is an important consideration because astigmatism does not remain constant during life, with changes in both magnitude and orientation occurring with age [1]. For example, in childhood, astigmatism tends to be WTR, whereas in older adults this orientation typically changes to ATR. Our study design sought to overcome some of this variation by considering only the magnitude of corneal astigmatism (i.e., no consideration of astigmatism axis) and by using a case-control classification scheme, with the aim of reducing the impact of the subtle changes in astigmatism that commonly occur with age.

In conclusion, this GWAS meta-analysis for corneal astigmatism replicated the discovery of a genome-wide significant locus near the PDGFRA gene [16] and provided strong evidence that this locus is important in both Asians and Europeans (Figure 4). Three novel candidate genes, CLDN7, ACP2, and TNFAIP8L3, were identified using gene-based analyses that leveraged data from across genomic regions rather than from examining one genetic marker at a time. These novel genes warrant further investigation to understand their role in the pathogenesis of corneal astigmatism. Finally, exploiting the recently introduced LD score regression technique, we estimated the SNP heritability of corneal astigmatism (and refractive astigmatism) to be much lower than that for spherical equivalent refractive error (Table 3) [51]. This implies that astigmatism must be under greater influence of rare genetic variants or environmental risk factors than spherical equivalent or that the common genetic variants that contribute to astigmatism have non-additive effects.