Introduction

Primary open-angle glaucoma (POAG), a leading cause of blindness, is characterized by visual field (VF) loss due to retinal ganglion cell death.1 It is believed that there are multifactorial contributors to POAG aetiology in addition to elevated intraocular pressure (IOP), an established risk factor. One major contributor may be regulation of vascular tone.2, 3, 4, 5, 6 For example, glaucoma patients have impaired blood flow and vascular autoregulation in eye tissues,7, 8, 9, 10, 11 and exhibit impaired flow- 12 or acetylcholine-mediated13 vasodilation in the brachial artery and reduced blood flow in finger nailfold capillaries,14 suggesting systemic vascular abnormalities. Of note, POAG patients with initial paracentral VF loss are particularly more likely to have systemic vascular risk factors, such as hypotension and migraines.15, 16, 17

Many studies support the role of vascular factors in glaucoma. Large epidemiologic studies have observed adverse associations with low ocular perfusion pressure18 or low systemic blood pressure,19 and open-angle glaucoma frequency or risk.20, 21 In addition, the Early Manifest Glaucoma Trial22 identified predictors of POAG disease progression that included lower systolic blood and perfusion pressure.

Recent large POAG genome-wide association studies (GWAS) have identified significant associations for single-nucleotide polymorphism (SNP) in the genomic region containing CAV1,23, 24, 25 which codes for caveolins that are involved in vascular regulation, thereby supporting the role of vascular dysregulation in POAG. Because biologic pathway analyses (where SNP sets from functionally related genes are collectively evaluated in relation to a disease of interest) can enhance the power to discover aetiologically relevant pathways and networks, we tested the hypothesis of whether a collection of SNPs in genes functionally involved in vascular tone regulation and identified as significant hits in GWAS of blood pressure are possibly associated with POAG.

In this study, we used high-throughput genotyping SNP data from the Glaucoma Genes and Environment (GLAUGEN) study and the National Eye Institute Glaucoma Human Genetics Collaboration (NEIGHBOR) consortium to assemble a vascular tone SNP panel (8671 SNPs in 186 genes) and used the Pathway Analysis by Randomization Incorporating Structure (PARIS) analysis software package (freeware developed by Dr Marilyn Ritchie, Pennsylvania State University, University Park, State College, PA, USA)26 for pathway analysis. In exploratory analyses, we assessed whether collections of SNPs within genes that comprise the pathway are associated with POAG.

Material and methods

Study population

We used meta-analysed GWAS data from two studies in the United States: GLAUGEN and NEIGHBOR. The GLAUGEN study included nested case–control samples (976 POAG cases and 1140 controls) from the population-based cohorts of the Nurses’ Health Study (NHS) and Health Professionals Follow-up Study (HPFS), as well as a clinical case–control sample from the Massachusetts Eye and Ear (MEE). The data from GLAUGEN have been evaluated in investigations of genetic factors and gene–environment interactions in POAG.27, 28 A recent publication describes the criteria for inclusion and exclusion for GLAUGEN,23 and the details are also available at dbGaP, accession number phs000308.v1.p1. The NEIGHBOR study included 2132 POAG cases and 2290 controls from 12 institutions and also provided data for genetic investigations of central corneal thickness,29, 30 IOP,31 and POAG.28 The descriptive information about the design, inclusion/exclusion criteria, study sites, and collected data in the NEIGHBOR study was recently published.32

POAG case and control definitions

All POAG cases were required to have at least one eye with either (1) a cup-to-disc ratio (CDR)>0.7 and VF loss consistent with glaucoma on one reliable VF test or (2) the same VF loss evident on at least two reliable VF tests (reliability of a VF test was defined as fixation loss ≤33%, false positives ≤20%, and false negatives ≤20%). The reproducible VF defect(s) evident on reliable tests were compatible with glaucomatous patterns of nerve fibre layer loss. The type of VF test was not specified, although most cases had threshold VF tests on the Humphrey perimeter. On the basis of slit-lamp biomicroscopy, none of the cases had evidence of exfoliation syndrome, pigment dispersion syndrome, or other glaucoma diagnoses that could produce secondary elevation of IOP. In addition, all cases had documentation of open iridocorneal filtration angles in both eyes. POAG cases were included at the full spectrum of IOP at diagnosis, including both normal tension glaucoma cases with IOP ≤21 mm Hg and high tension glaucoma cases with IOP >21 mm Hg. Controls identified in clinical settings underwent eye exams that showed IOP ≤21 mm Hg and CDR ≤0.6 in both eyes. For nested case–control study controls (from NHS and HPFS), controls reported having had an eye exam at the same time period as the index cases and did not report a glaucoma diagnosis.

Each case was categorized into POAG subtypes defined by type of VF loss: early paracentral VF loss (with or without peripheral VF loss) or peripheral VF loss only (without any paracentral loss). The earliest VF report for each case was systematically reviewed such that pattern deviation plots were subdivided into paracentral, Bjerrum, nasal step, and temporal wedge zones above and below the horizontal meridian (Figure 1). Two reviewers who were masked to genotype status assessed VFs, and any differences in VF loss-type designation were adjudicated for establishing a consensus designation. The presence of clusters of three or more contiguous points with reduced retinal sensitivity of one-half log unit (−5 dB) relative to age-matched controls were considered to show VF loss; the location of such loss was determined for each of the four regions. Cases were defined as having POAG with ‘early paracentral loss’ in the two studies differently: in the GLAUGEN data set, where cases were incident cases from population-based studies with early VF loss, cases had VF loss that occurred in the paracentral zone with or without loss in other peripheral zones, whereas in the NEIGHBOR data set, where cases were prevalent cases from tertiary clinic-based case–control studies with advanced disease, to capture ‘early paracentral loss’, cases had to have VF loss that occurred in the paracentral zone only without any loss in other peripheral zones. Cases were defined as having POAG with ‘peripheral loss only’ in both GLAUGEN and NEIGHBOR data sets if cases had VF loss in any of the other peripheral zones without loss in the paracentral zone.

Figure 1
figure 1

Paracentral and peripheral VF loss definitions. Representative VF pattern deviation plot used for defining the type of VF loss. The box with a dashed line indicates the superior paracentral region. The boxes with solid lines represent the superior nasal step and temporal wedge regions, as well as the inferior Bjerrum region, which are all peripheral regions. For a case to be defined as having early paracentral loss, VF loss in the superior or inferior paracentral zones had to be present; in the GLAUGEN data set, VF loss in other peripheral zones were allowed, whereas in the NEIGHBOR data set, VF loss in other peripheral zones were not allowed. For a case to be defined as having peripheral loss, VF loss in the superior and/or inferior nasal step, temporal wedge, or Bjerrum regions had to be present, with no loss in the paracentral regions.

Genotyping and analyses of genotype associations

Genotyping with the Illumina 660W-Quad array (Illumina, San Diego, CA, USA) was conducted at the Broad Institute (Cambridge, MA, USA) for GLAUGEN and at the Center for Inherited Disease Research (Baltimore, MD, USA) for NEIGHBOR. The details of the methods used for DNA collection, extraction and plating, the algorithm for calling genotypes, the SNP quality control (QC) filters used, and the preliminary analyses undertaken to discover independent predictors of genotype calling rates adjusted for in this paper have been provided in previous publications for both GLAUGEN and NEIGHBOR.23, 28 For analyses of genotype associations with POAG, PLINK v.1.0733 was used with SNPs that passed the QC filters (484 419 SNPs in GLAUGEN and 521 683 SNPs in NEIGHBOR). Analyses were conducted among self-reported Caucasians and were adjusted for population structure, age, gender, and study site in both study groups, and analyses in GLAUGEN were additionally adjusted for the three different DNA extraction methods used and DNA source (blood or cheek cell sample).23, 28 Using the METAL software,34 the GLAUGEN and NEIGHBOR data were meta-analysed, and pooled P-values were obtained for the 480 335 overlapping SNPs in the combined data set.

Defining the vascular tone pathway and analysis of the whole pathway and of specific genes

In this paper, we used the P-values for the association with POAG for 8671 unique SNPs in the vascular tone pathway for the analysis that is discussed below. We generated a custom list comprising the vascular tone pathway (see Figure 2) of 8671 SNPs in 186 genes (Supplementary Table 1) across all autosomes. We used the Kyoto Encyclopedia of Genes and Genomes (KEGG)35 online database and other academic sources,36, 37 as well as a recent review38 of the GWAS findings for blood pressure.

Figure 2
figure 2

Vascular tone pathway. This pathway was custom-generated by using sources from the Kyoto Encyclopedia of Genes and Genomes35 online database and other academic sources,36, 37 as well as a recent review38 of the GWAS findings for blood pressure. Proteins coded by genes included in the analysis are in grey. Substrates and reaction products are in green. Overall, the vascular tone pathway was not significantly associated with POAG permuted P=0.87. In exploratory analyses, of the eight proteins of genes with permuted P-values <0.001, seven are indicated with red squares, with the remaining being the HFE gene, which was included in the analysis as part of genes identified from a review38 of known blood pressure SNPs.

We uploaded the meta-analysed P-values of the selected SNPs within these genes residing in the EMSEMBL genomic interval±a 50 kB window to PARIS software from the combined GLAUGEN–NEIGHBOR data set. The analytic approach with the PARIS software has been described.26 Briefly, for the assembled pathway, PARIS first determined important parameters such as the number of linkage disequilibrium (LD) blocks, the number of SNPs per block, and the number of SNPs not in any LD block. The vascular tone pathway contained 1405 complex features (LD blocks with two or more SNPs) and 1779 simple features (single SNPs not in any LD block). Then, PARIS conducts a comparison of the overall association of the vascular tone pathway with POAG to the associations of 1000 randomly generated pathways with genetic architectures similar to the vascular tone pathway with POAG; in this comparison, PARIS provides a permuted P-value, which helps determine whether the assembly of SNPs in the pathway of interest as a whole may be significantly associated with POAG. We chose to produce 1000 random pathways, as generating another order of magnitude more pathways with similar genetic architecture would have been prohibitive computationally. To illustrate the methodology of how the random SNP sets are used to assess statistical significance for SNP sets of interest, we, as described previously,39 evaluated an oestrogen pathway among women using the same data set. The results indicated 31 significant features in this pathway, and PARIS generated a significant permuted P-value=0.007. This indicated that only 7 of 1000 random pathways with a genetic architecture similar to the oestrogen pathway had more features that were statistically significant at P<0.05 level in relation to POAG than the oestrogen pathway (i.e., >31 significant features with P<0.05). This overall significant permuted P-value of 0.007 indicated that the assembly of SNPs in the oestrogen pathway was significantly associated with POAG, supporting the role of oestrogen in the aetiology of POAG.39 We used a similar approach in this study to evaluate the vascular pathway in relation to POAG. Because POAG characterized by VF loss in the paracentral region is considered to be more strongly associated with vascular abnormalities,15, 16, 17 in addition to overall POAG, we evaluated subtypes of POAG defined by the type of VF loss separately (POAG with early VF loss in the paracentral region and POAG with VF loss in the peripheral region only), to investigate differences in association.

In subsequent analyses, to determine whether some genes or SNPs in the vascular tone pathway were significantly associated with POAG, we used the ‘-I’ (Investigate) option within PARIS, which reports the P-values of genes and SNPs within the pathway. In an approach that was similar to the overall pathway approach, we assessed the association with POAG of the aggregate SNPs within specific genes in the vascular tone pathway to gain insight into which aspects of the pathway may be aetiologically important. Because we evaluated 186 genes in this pathway in relation to three outcomes, there was the potential for significant findings arising owing to chance because of multiple testing. In the PARIS methodology, the genome is permuted rather than applying a Bonferroni correction for multiple testing (to avoid the risk of penalizing larger genes as might happen with Bonferroni correction).26 To minimize false-positive findings, we focused on genes that showed permuted P<0.001 (which represented the lowest possible permuted P-value, as such P-values are based on counts of 1000 random pathways per genes; that is, a count of 1 randomly generated pathway that shows more significant features than the vascular tone pathway out of 1000 random pathways yields a permuted P=0.001, and 0 out of 1000 random pathways leads to a permuted P-value of P<0.001).

We certify that all applicable institutional and governmental regulations concerning the ethical use of data from human volunteers were followed during this research, and the Institutional Review Boards of the participating institutions of GLAUGEN and NEIGHBOR approved this study.

Results

A total of 3430 controls (2290 from NEIGHBOR and 1140 from GLAUGEN) and 3108 cases (2132 from NEIGHBOR and 976 from GLAUGEN) and were analysed in the study (Table 1). Of these 3108 cases, 224 had POAG with early paracentral VF loss, whereas 993 had POAG with peripheral VF loss only.

Table 1 Pathway analysis results using data from the combined GLAUGEN and NEIGHBOR: significance of the association of the overall vascular tone pathway in relation to POAG, POAG with early paracentral visual field loss, and POAG with peripheral visual field loss only

The vascular tone pathway was not significantly associated with POAG overall (permuted P=0.87); the results were similarly null (permuted P≥0.81) in men and women analysed separately (data not shown). In addition, the vascular tone pathway was not associated with POAG subtypes defined by VF loss: the permuted P-value in relation to early paracentral VF loss was 0.20, whereas the corresponding P-value for peripheral VF loss only was 0.39.

In gene-based analyses, eight genes were associated with overall POAG at permuted P-value <0.001 (Table 2), and the red-boxed items in Figure 2 are the corresponding proteins for seven of the eight (HFE is not included; HFE was a gene with loci that were significantly associated with blood pressure in a recent review38 of the GWAS studies of blood pressure). Of these eight genes, six were those coding for proteins/enzymes that interact with endothelial nitric oxide synthase (eNOS): PRKAA1, CAV1, ITPR3, EDNRB, GNB2, and DNM2. Supplementary Table 2 shows the association between individual SNPs within these eight genes (residing in the EMSEMBL genomic interval±a 50 kB window) and POAG.

Table 2 Gene analysis results for the vascular tone pathway using data from the combined GLAUGEN and NEIGHBOR

In addition, three of these (CAV1, ITPR3, and EDNRB), along with HFE, were also associated with the subtype of POAG with early paracentral VF loss with the permuted P-value of <0.001, and PRKAA1 showed borderline significance with the permuted P-value of 0.05. In contrast, none of these reached permuted P-value <0.001 for the subtype of POAG with peripheral VF loss only. However, as observed for POAG with early paracentral VF loss, nominally significant permuted P-values were observed for CAV1 (P-value=0.004) and PRKAA1 (P-value=0.05) in relation to POAG with peripheral VF loss only.

For the subtype of POAG with early paracentral VF loss, one other gene, ADCY7, which was not associated with POAG overall, was found to be associated with this subtype with permuted P<0.001 (data not shown). Notably, ADCY7 codes for a protein that influences eNOS activity. For the subtype of POAG with peripheral VF loss only and not overall POAG or the subtype with early paracentral VF loss (data not shown), there were eight different genes that were associated only with this subtype with permuted P<0.001: PRKAA2, PRKAB2, GNB5, PIK3CG, WNK1, PRKCA, ADORA2A, and AGTR1. Of note, four of these 8 (PRKAA2, PRKAB2, GNB5, and PIK3CG) also coded for proteins that interact with eNOS.

Discussion

In this first study of the collective SNPs in genes involved in the vascular tone pathway and POAG, we observed that the assembled vascular tone pathway SNP panel was not related to the main outcome of overall POAG diagnosis or the subtypes of POAG with early VF loss in the paracentral region or POAG with VF loss in the peripheral region only. However, in gene-based analysis, eight genes were associated with POAG with permuted P<0.001, and of these eight genes, six coded for proteins that interact with eNOS. In particular, three of these six were also significantly related (permuted P<0.001) to the POAG subtype that is characterized by early paracentral VF loss. Thus, although this study did not find associations with the vascular tone panel as a whole, it generated further support for the more specific role of eNOS-related activity in POAG.

Because a myriad of extrinsic (factors that originate from outside the local tissue that affect the arterial blood pressure) and intrinsic factors (local factors that regulate local blood flow) affect vascular tone,3 the panel tested included a comprehensive set of 186 genes (8671 unique SNPs) based on the KEGG pathways and GWAS significant associations for blood pressure based on a recent review.38 The overall null association may have been due to the inclusion of SNPs from gene-regulating factors that may affect blood pressure or regulate vascular tone systemically but are not relevant for blood flow regulation, especially to the optic nerve head (ONH). For example, hormones such as epinephrine, vasopressin, and natriuretic peptides40, 41, 42 or neurogenic factors43 are unlikely to affect the regulation of blood flow to the retina and the ONH (e.g., blood vessels in the retina and the ONH are not innervated),40 whereas intrinsic factors such as vascular endothelial vasoactive agents have an important role in influencing blood flow autoregulation and vascular tone in the ONH.3, 44, 45 Of these agents, NO, generated constitutively by eNOS, has a key role in maintaining the basal vasodilator tone in ocular arteries.40, 46, 47, 48, 49, 50, 51, 52 The major contribution of NO to vascular tone in human ocular vasculature is demonstrated by the observation that systemic administration of an NOS inhibitor reduced blood flow in the ONH53, 54, 55 and by the involvement of NO in the POAG-related vascular conditions, such as hypertension, hypotension, and vasospasms.56

Evidence of abnormal NO production in glaucoma is also provided by genetic studies. Although genetic studies have evaluated SNPs of factors related to the renin–angiotensin system57, 58, 59 and the β-adrenoreceptors60, 61 in relation to POAG, they have been small and conflicting. However, relatively more consistent associations are observed with genetic studies of eNOS-related factors and POAG.62, 63, 64, 65, 66, 67 Although not all studies have found associations with variants in the NOS3 gene that codes for eNOS,68, 69 variants or haplotypes of SNPs in NOS3 have been reported to be associated with POAG overall,62, 65 POAG with migraines,66 POAG or high tension POAG among women only,63, 64 and familial POAG.67 Furthermore, the CAV1/CAV2 gene region, which codes for caveolins that control eNOS activity, contain variants associated with POAG in one GWAS,24 and this finding was replicated in another large candidate gene study.23 The caveolin/eNOS interaction is critical to prevent inadequate NO production under basal conditions and to the translation of extracellular stimuli to intracellular NO signals.70

Consistent with this previous literature, it is suggestive that six of the eight most significant associations in the gene-based analyses were genes encoding proteins that interact with eNOS; based on information from large databases of gene expression in various human tissues, all eight genes showed expression in human retinal tissue.71, 72 ITPR3 codes for inositol triphosphate receptor (IP3R); when IP3 binds to IP3R, calcium levels increase, activating calmodulin, and causing dissociation of caveolin from eNOS. The increased calcium levels directly or indirectly via dynamin 2 (Dyn2), coded by DNM2, activate eNOS. EDNRB codes for endothelin receptor type B (ETB, which is particularly highly expressed in retina71, 72) and PRKAA1 codes for AMP kinase (AMPK); both are additional activators of eNOS. Finally, GNB2 codes for G-protein subunit-β (G beta) that activates PI3K, which in turn activates Akt; Akt then activates eNOS. A limitation of this study was that variants in NOS3, which codes for eNOS, were not significantly associated with POAG, although this could be due to the low coverage (only two non-functional SNPs rs7830 and rs3918188 were included) of the NOS3 SNPs on the Illumina 660W quad platform. This might also be observed if NOS3 effects derive primarily from multiple different rare variants rather than from effects of the kinds of common variants used in this study.

Out of the eight genes showing significant association with POAG, only MYL9, coding for myosin regulatory light polypeptide 9 (MLC) and HFE, which codes for hemochromatosis protein involved in iron metabolism, are not directly involved in eNOS activation. Interestingly, HFE was associated with both POAG overall and the POAG subtype with early paracentral VF loss. Although there are no previous reports of an association between hemochromatosis and glaucoma, iron supplementation has been linked to POAG;73 thus, the association between POAG and iron levels or hemochromatosis warrants further study.

Of the eight genes with P<0.001 associated with POAG overall, three (CAV1, ITPR3, and EDNRB) genes involved in eNOS activity were also associated with the subtype of POAG with early paracentral VF loss, even though we included only 224 such cases. In addition, ADCY7, which codes for adenylate cyclase that activates eNOS, was also significantly related to this subtype. In contrast, none of the genes significantly associated with POAG overall or the POAG subtype with early paracentral VF loss were associated with the POAG subtype with only peripheral VF loss. This is consistent with previous findings showing associations between the POAG subtype showing early paracentral VF loss with systemic vascular dysregulation.15, 17

Among the limitations of this study, first, we lacked replication of the gene-based analyses as there is currently no other data set of comparable size for replication. Thus, although suggestive, these gene-based exploratory analyses, especially in relation to the subtype of POAG with early paracentral VF loss, should be interpreted with caution and replicated. Second, this study was restricted to Caucasians and, therefore, the generalizability of these results to other populations, especially East Asian populations with generally lower IOP,15, 17, 57 is limited.

In this large study to evaluate the collective influence of genes coding for proteins influencing vascular tone, we did not observe associations with POAG. However, in gene-based analysis, we observed confirmatory results supporting the role of the dysregulation of eNOS activity in POAG aetiology.