Introduction

Primary open-angle glaucoma (POAG) refers to a heterogeneous disease resulting in a characteristic excavation of the optic nerve head and progressive optic neuropathy. POAG is a leading cause of irreversible visual impairment and blindness worldwide. Risk of developing POAG is influenced by a number of genetic and environmental factors. To date, several large-scale genome-wide association studies (GWASs) have identified a number of genetic loci linked to POAG susceptibility [1,2,3]. However, these loci only contribute a minor effect on POAG risk with odds ratios (ORs) for individual risk loci not exceeding 1.5. Thus, there remains a large amount of unexplained genetic risk (heritability) for POAG.

In an effort to gain more knowledge of the genetic basis of the disorder researchers are also focusing on quantitative traits associated with POAG [4]. Quantitative endophenotypes that allow individuals to be scored along a continuum can offer increased statistical power for gene mapping studies. Well-established POAG endophenotypes such as intraocular pressure (IOP) and central corneal thickness (CCT) have received considerable attention with large-scale genome-wide association studies identifying many potential risk variants underlying POAG [5, 6].

Founder-effect genetic isolates may offer some advantages for gene mapping due to features such as reduced genetic heterogeneity and extensive familial structure. Linkage and association analyses of endophenotypes for POAG in a Dutch isolate revealed new linkage regions for optic disc morphology at 20p13 (SIRPA and RNF24/PANK2 loci) and 2q33–q34 (IGFBP2 locus). Significant and suggestive genome-wide association studies signals at the previously identified RERE, LRP1B, CDC7, TGFBR3 and ATOH7 loci were replicated in this isolate [7].

Despite the identification of many genetic risk loci contributing to POAG and its underlying endophenotypes, the functional relevance of most genetic variants remains unknown. Expression quantitative trait loci (eQTLs) represent genomic regions harbouring SNPs that have a strong influence on transcript expression level in either a cis (proximal) or trans (distal) location to the gene. Genetic association studies that target known eQTLs and their respective transcripts have the advantage of allowing researchers to restrict their analysis to loci most likely to have a functional effect on the trait of interest.

The Norfolk Island (NI) population, located off the east coast of Australia, is a genetic isolate with well defined familial lineages dating back 12 generations to the originating founders—12 Tahitian women and 6 European men. A previous study of the NI cohort successfully mapped 197 cis-eQTLs in peripheral blood mononuclearcytes. These eQTLs are associated with immune system processes and have been linked to a number cardiovascular disease risk factor traits in NI [8].

As part of the Norfolk Island Eye Study (NIES) the prevalence of chronic ocular diseases has previously been reported [9]. In general, the prevalence of blindness and visual impairment was observed to be lower than mainland Australia, while the prevalence of POAG in adults ≥40 years was 6%, higher than other Australian studies. Individuals with glaucoma were more likely to have a positive family history of disease than those without. Observed differences between NI pedigree and non-pedigree groups for a number of POAG risk factor traits suggests a possible genetic influence [9].

In an effort to better understand the genetic basis of POAG we have conducted an eQTL-centric association analysis focusing on POAG risk factor traits in the NI isolate.

Materials and methods

NI cohort

The NIES has been well described in previous research [9, 10]. In this study we focused on a group of 330 individuals selected from the “core” NI pedigree, who relate back to the original NI founders. For all of these individuals we had ocular phenotype data, SNP genotype data and gene transcript data.

Quantitative POAG endophenotypes

Being a heterogeneous disease, an understanding of the genetic architecture of POAG can be achieved by examining individual traits (endophenotypes) underlying the disease, and this approach can also improve statistical power in genetic analyses. In this study we analysed five quantitative traits representing POAG endophenotypes in the NI cohort. All measurements, ACD = anterior chamber depth (mm), CCT = central corneal thickness (μm), vertical CDR = cup-to-disc ratio, disc size, and IOP = intraocular pressure (mmHg) were obtained as part of the NIES [9, 10].

A larger CDR is associated with increased structural damage of the optic nerve head in POAG and it correlates strongly with IOP, the strongest modifiable POAG risk factor and major therapeutic target [11]. While size of the optic disc does not directly correlate with POAG, its knowledge is important as it influences interpretation of the degree of glaucomatous optic neuropathy [12], and genetic factors are known to influence disc size. CCT is not a true POAG endophenotype, as it has not been shown to have significant association with POAG risk [13]. Nonetheless, it represents an important component of POAG evaluation. CCT influences IOP readings, and it represents an independent and strong risk factor for development of POAG [14]. CCT also represents one of the most heritable traits in humans, reflected by the identification of SNPs associated with this trait [15]. ACD correlates strongly with the size of the eye, and in general, higher values of ACD are associated with POAG with lower values associated with closed-angle (angle-closure) glaucoma.

Phenotypic examination

ACD was measured with the intraocular lens (IOL) master (Carl Zeiss Meditek, Dublin, CA, USA). CCT was measured using the Pachmate ultrasonic pachymeter (DGH, Exton, USA). IOP was measured using a Tonopen handheld tonometer (Reichert, Tustin, CA, USA). Dilated posterior segment examination was performed using the slit lamp with condensing lens (Volk 78D).

eQTL SNP and transcript data

This study is an eQTL-centric design whereby we analysed genotype data for SNPs that have previously been robustly associated with gene transcripts in circulating peripheral blood mononuclearcytes in the NI pedigree [8]. Specifically, we selected the single most significant SNP (as index) for 197 heritable transcripts spanning the genome (see Supplementary Table 1). We chose to focus our analysis on cis-acting eQTL SNPs only given these are most likely to have a direct functional influence on the trait. Briefly, genotype data were generated using the Illumina 660 Quad array and QC was performed using methods described in Benton et al. (2013) [8]. For the current study we included corresponding gene transcript level data for the same individuals obtained using Illumina HT12 expression arrays.

Statistical analysis

Trait heritability estimation

For each POAG endophenotype, we estimated heritability using a variance decomposition based on genotypic similarities among subjects as implemented in Genome-wide Complex Trait Analysis (GCTA). These analyses included age, sex and the top two principal components as covariates. In our GCTA analyses, we employed a genetic relationship matrix (GRM) based on the full set of the autosomal SNPs. The significance threshold was set at 0.01.

eQTL SNP association testing

We performed association analyses of index SNPs from 197 cis-eQTLs previously identified in NI [8]. We performed the analysis in GCTA using linear mixed models fitting a genetic relationship matrix of random effect and age, sex and top two principal components as fixed effect covariates. We then utilised expression data representing the eQTL gene to assess whether the eQTLs could be acting on POAG trait via gene expression. To confirm this hypothesis, we ran associations using linear mixed models in GenABEL with the kinship matrix to account for relatedness and included sex, age and top two principal components as covariates. For each association we inspected the effect of gene expression on phenotype and the effect of the eQTL on trait, conditioned on gene transcript.

Statistical significance and multiple testing considerations

For eQTL–centric association testing, each trait was considered as a separate experiment even though there are intercorrelations among some traits. Within each experiment we used 0.05 as the nominal significance level. Given the relatively modest size of the NI cohort and the over inflation of Type II error under strict Bonferoni adjustment, we chose to slightly relax the significance levels to attain a reasonable balance between Type I and Type II error rates. For the association analysis we tested 197 cis-eQTLs (SNPs). Considering a balanced adjustment we set the significance threshold at 0.005 (Bonferoni = 0.05/197 = 0.0003).

Results

Genetic heritability was estimated for five key POAG endophenotypes in the NI cohort (Table 1). All traits exhibited statistically significant heritability scores (H 2) with estimates ranging from 0.35 for IOP to 0.82 for CCT (P < 0.01). Given the statistical evidence of genetic inheritance, all traits underwent subsequent eQTL SNP and transcript association analyses. To explore the relationship between heritable ocular traits and transcripts, we first performed regression analysis between all 197 eQTL-related transcripts (exclusive of SNP) and the five heritable traits (see Supplementary Table 2). Table 2 shows the most significant transcript associations for each of the traits. The most significant association was between IOP and DDX17 (P = 0.000099). DDX17 is a Dead-box Helicase, which has been implicated in a number of cellular processes involving alteration of RNA secondary structure and linked to pathways associated with antiviral immunity [16].

Table 1 Heritability estimates of ocular traits
Full size table
Table 2 eQTL transcript by trait association
Full size table

We then tested all 197 cis-eQTL SNPs for association with each of the five POAG traits (exclusive of transcript). The top 10 eQTL-by-trait association results are shown in Supplementary Table 3. Table 3 shows the eQTL-by-trait associations that were statistically significant according to our threshold of 0.005 (see Materials and methods). Significant associations were detected for four out of five heritable traits. No associations were seen for CDR. The most significant association observed was between Disc size and LPCAT2 (P = 0.0004). The eQTL for BTN3A2 was associated with both Disc size and CCT. Finally, all eQTL SNPs associated with an ocular endophenotype were tested using conditional regression by considering the relationship of both eQTL SNP and its transcript on the variation in ocular endophenotype. The only model whereby eQTL SNP and transcript were both statistically significant was for BTN3A2 and Disc size. This model included rs853676 (β = 0.23, P = 0.008) and ILMN-1700067 (β = 0.23, P = 0.03). These results indicate that the SNP-influenced transcript for BTN3A2 is associated with variation in Disc size in this cohort (model R 2 = 0.05, P = 0.024).

Table 3 eQTL SNP by trait association
Full size table

Discussion

In this study we performed heritability analysis of five key endophenotypes for POAG in the NI genetic isolate. All traits showed statistically significant evidence of heritability and thus are likely to be influenced, at least in part, by genetic factors. The trait with the highest heritability estimate was CCT (H 2 = 0.82), which is consistent with other studies [17]. Despite the substantial heritability for POAG and risk traits in the general population the identification of genetic variants has been limited to a relatively small number of loci that collectively do not explain a large proportion of the heritability and/or do not have a direct causal relationship to POAG. Alternate research designs that specifically focus on eQTLs and that utilise unique founder-effect genetic isolates may add new insights into the genetics of POAG risk and pathology.

In this study we focused on 197 eQTLs previously identified in the NI isolate [8]. These eQTLs represent heritable gene expression levels derived from peripheral blood mononucleocytes. Given the target cell types in this study it is not surprising that this panel of genes is associated with biological processes of the immune system. Therefore, by design this study was focused on blood cell specific immune-related genes. The relationship between immunity and POAG has received relatively little research attention. However, the evidence is mounting that molecular pathways of inflammatory response mechanisms are associated with POAG pathologies and it is plausible that blood-based eQTLs may influence ocular tissue directly and/or act as biomarkers of disease risk.

IOP is a primary risk factor for POAG and can cause insult to the optic nerve head leading to retinal ganglion cell loss and POAG. It has been hypothesised IOP-related factors can trigger an immune response [18]. Bakalash et al. [19] showed that in rats, resistance of retinal ganglion cells to an increase in IOP is dependent on a T-cell mediated immune response. Ahmed et al. [20] performed a microarray experiment in a rat glaucoma model and found a set of neuroinflammatory gene transcripts were altered in retinal tissue after elevation of IOP. Johnson et al. [21] also performed microarray analysis of the rat glaucoma model and identified 38 genes associated with the biological process of immune response, including complement 1 complex, FC receptor, beta-2 microglobulin, and several histocompatibility class II antigens—all genes expressed by activated microglia in the nervous system.

Our primary finding in this study was the association of a genetically influenced transcript in BTN3A2 with Disc size. Disc size is associated with variation of specific anatomical structures of the optic nerve head and the retinal nerve fibre layer. Thus variation in disc size may influence risk of glaucoma [22]. BTN3A2 (Butyrophilin Subfamily 3 Member A2) is a member of the immunoglobulin superfamily and is part of the major histocompatability class 1 locus on chr6 and may play a role in T-cell response and adaptive immunity. To date, BTN3A2 has not been previously associated with disc size or POAG. However, studies have shown a link between environmental exposures and other immunoglobulins in POAG. Yuki et al. [23] reported elevated serum immunoglobulin G titre against Chlamydia pneumoniae in patients with POAG. Furthermore, Tseng et al. [24] reported an association between glaucoma and immunoglobulin E antibody response to indoor allergens. These findings generally suggest a role for immunoglobulins and the genes that regulate them in POAG risk.

Disc size was also associated with an eQTL for LPCAT2. There was no evidence that this eQTL acted through its transcript in this cohort but this could be a consequence of modest sample size. LPCAT2 encodes an enzyme that is a member of the lysophospholipid acyltransferase family and may have a function in membrane biogenesis and production of platelet-activating factor in inflammatory cells [25]. Although not directly relevant to POAG, Singh et al. [26] used transcriptomics analysis of peripheral blood cells to show that expression of LPCAT2, among other genes, was associated with allergen response in asthmatics and suggested a lipid imbalance resulted in an asthma inducing inflammatory response in this cohort.

In this study we have used an eQTL-centric approach to implicate a number of immune-related genes with POAG risk traits. Most notably, we provide compelling evidence that an eQTL for BTN3A2 is associated with variation in optic disc size. Considered together our results provide another layer of evidence implicating the genetic factors of the systemic immune system in POAG. This study also illustrates the value of isolated founder-effect populations for studying the genetics of complex diseases like POAG.