Primary closed angle glaucoma in the Basset Hound: Genetic investigations using genome-wide association and RNA sequencing strategies.

Purpose: To investigate the genetic basis of primary closed angle glaucoma (PCAG) in European Basset Hounds using genome-wide association and RNA sequencing strategies.
Methods: DNA samples from 119 European Basset Hounds were genotyped on the 170 K SNP CanineHD BeadChip array (Illumina) comprising 37 with normal iridocorneal angles (controls), 57 with pectinate ligament abnormality (PLA cases), and 25 with PCAG (PCAG cases). Genome-wide association studies (GWASs) of the PLA and PCAG cases were conducted. Whole transcriptome sequences of iridocorneal angle tissues from five Basset Hounds with PCAG were compared with those from four dogs with normal eyes to investigate differences in gene expression between the affected and unaffected eyes in GWAS-associated loci. A variant in NEB, previously reported to be associated with PCAG in American Basset Hounds, was genotyped in cohorts of European Basset Hounds and non-Basset Hounds.
Results: The GWASs revealed 1.4 and 0.2 Mb regions, on chromosomes 24 and 37, respectively, that are statistically associated with PCAG. The former locus has previously been associated with glaucoma in humans. Whole transcriptome analysis revealed differential gene expression of eight genes within these two loci. The NEB variant was not associated with PLA or PCAG in this set of European Basset Hounds. Conclusions: We identified two novel loci for canine PCAG. Further investigation is required to elucidate candidate variants that underlie canine PCAG.

Introduction

Canine glaucoma is a heterogeneous group of neurodegenerative diseases that lead to progressive retinal ganglion cell death, optic nerve degeneration, and visual field deficits. Intraocular pressure (IOP) is a consistent risk factor for the development of canine glaucoma [1]. Glaucoma is described as primary if it occurs in the absence of an antecedent ocular disease process, such as intraocular inflammation, neoplasia, or lens instability [2]. Primary glaucoma is thought to be caused by unidentifiable and inherent abnormalities in the aqueous humor outflow apparatus of the eye and in adult dogs is further subdivided into primary open angle glaucoma (POAG) and primary closed angle glaucoma (PCAG) based on the appearance of the iridocorneal angle (ICA). In POAG, the IOP is usually pathologically elevated, even in the presence of an open ICA and normal pectinate ligament anatomy in the early disease stages [3]. PCAG, however, is associated with pectinate ligament abnormality (PLA), a form of goniodysgenesis, in several dog breeds [4-10]. PLA describes the broad sheets of tissue that span the ICA and is associated with increasing age [11-14]. Although PLA appears to be required, but not sufficient, for the development of PCAG, only a small fraction of dogs with PLA develop PCAG [15,16].

The increased prevalence of primary glaucoma in certain dog breeds implies a genetic etiology. The Basset Hound (BH) is a dog breed known to be affected by POAG and PCAG [4,17]. In common with three other canine breeds, POAG in the BH is an autosomal recessive trait due to a mutation in ADAMTS17 (Gene ID 170691, OMIM 607511) [17-20]. PCAG, however, is thought to be complex and caused by multiple genetic and environmental factors [21-23]. Two previous studies of the genetics of PCAG have been reported in the American BH. A genome-wide association study (GWAS) of 37 PCAG cases and 41 controls revealed associations at two novel loci that contain the candidate genes COL1A2 (Gene ID 1278, OMIM 120160) and RAB22A (GENE ID 57403, OMIM 612966), but no candidate variants were reported [22]. In a subsequent paper, the same group reported a non-synonymous variant in NEB (Gene ID 4703, OMIM 161650) segregates with PCAG in the BH (BROADD2 chr19:55,885,214 A>G) [21]. COL1A2 and NEB are promising candidate genes based on their functions. COL1A2 encodes the pro-α2 chain of collagen type 1 which is an important component of the trabecular meshwork, and previous studies have implicated the role of collagen genes in the pathogenesis of human primary glaucoma [24-28]. NEB encodes the muscle contractility regulating protein, nebulin, which is expressed in the ciliary body musculature, and the belief that muscle-related mechanisms are involved in the aqueous humor outflow pathways makes NEB a feasible candidate [21,29,30]. RAB22A is known to be an oncogene, and thus, this gene’s potential role in primary glaucoma is more difficult to explain [31].

In this study, we investigated the genetic basis of PLA and PCAG in the European BH. First, we tested for a possible association between the reported NEB variant and PLA and PCAG. Second, we used a GWAS to identify loci associated with PLA or PCAG and then used an additional cohort of American BHs with PCAG to further evaluate statistically significant associations. We used RNA sequencing (RNA-Seq) to investigate differential gene expression and identify candidate genes in the GWAS-derived loci.

Methods

Genotyping of the NEB variant

The NEB variant (BROADD2 chr19:55,885,214 A>G) reported by Ahram et al. [21] was genotyped in a cohort of 158 BHs (10 PCAG cases, 52 PLA cases, and 96 controls) and 83 non-BH dogs, comprising 31 different breeds, with Sanger sequencing using primers detailed in Table 1. The non-BH dogs formed part of a multibreed screening panel that was created before the inception of this study. These dogs had not undergone previous ophthalmological examination.

Table 1

Primer details used for Sanger sequencing of the NEB variant.

Genomic coordinate of NEB variant (CanFam2)	Forward primer sequence Reverse primer sequence	Amplicon size (bp)	Annealing temp (°C)
chr19:55,885,214	ACCAGTAAGGTGAGTGCTTTCC AGGCTATGATCTCAGAACTGATGC	103	57

Genome-wide association study

All DNA samples from privately owned pet dogs were collected following fully informed and written owner consent and with approval from the Animal Health Trust’s Research and Ethical Committee (approval number 36–2016). This study was performed in accordance with the ARVO Statement for Use of Animals in Research. DNA samples from American PCAG dogs were provided by one of the coauthors (MK) and included members of a PCAG colony and client-owned dogs all originating within the USA. Dogs were designated as controls or PLA cases based on the results of gonioscopy performed by the primary author (JO), a board-certified veterinary ophthalmologist, and as previously described [12-14,32]. Dogs were designated as PCAG cases following examination by JO and other board-certified veterinary ophthalmologists. The inclusion criteria for the controls, PLA cases, and PCAG cases were as follows:

1. Controls: Dogs with normal appearing ICAs aged five years and older.

2. PLA cases: Dogs with PLA affecting at least 50% of the ICA of each eye. The controls and the PLA cases had IOP within normal limits and normal optic nerve head anatomy as assessed with direct and indirect ophthalmoscopy.

3. PCAG cases: Dogs with one eye with IOP >50 mmHg (arbitrary cutoff) without any possible cause of secondary glaucoma and the finding of severe PLA (>90% of the ICA affected) in the contralateral eye. No diagnostic techniques were performed on the affected eye to confirm closure of the ICA or to assess the pathology of the optic nerve head.

DNA for the GWAS was extracted from buccal mucosal swabs as previously described [17]. Only dogs clear of the published BH POAG mutation were used [17]. DNA samples were submitted at a concentration of 20 ng/μl and a volume of 25 μl to an external laboratory (Neogen; Lansing, MI) for genotyping on the CanineHD BeadChip (Illumina; San Diego, CA) which contains 172,115 single nucleotide polymorphisms (SNPs) [33]. The GWAS data were analyzed using the freely available software package PLINK [34]. Data were filtered for quality control parameters, including the sample call rate, SNP call rate, and minor allele frequency (MAF). SNPs were excluded from analysis if they had a MAF <0.05 or had a call rate of <97%. Individuals were excluded if >10% SNP genotypes were missing. GWASs were conducted using a standard unadjusted allelic chi-square test for association (1 degree of freedom (df)). Inflation of observed statistics due to relatedness or the population substructure was estimated with genomic control (λ). Data were subsequently corrected for population stratification using a mixed model in the freely available software GEMMA [35]. Five case and control GWASs were performed as detailed in Table 2. The p values for each SNP derived from the GWAS using GEMMA were used to create a Manhattan plot for each analysis. A p value of 0.05 after correction for multiple testing using the Bonferroni correction (0.05/number of tests) was the threshold for statistical association in all studies (p<5.28 × 10−7) and prompted further investigation of the associated loci. Associations between the GWAS SNPs and disease in individual dog data sets were assessed using logistic regression using an additive per allele model to compute an odds ratio and confidence interval for the association of one or more SNPs with disease using STATA 10.0 (College Station, TX).

Table 2

Summary of number of dogs, results of quality control filtering steps and correction of population stratification for the five GWAS analyses.

Analysis	Phenotype of cases	Phenotype of controls	No. cases	No. controls	SNPs with MAF <0.05	SNPs with <0.97 call rate	SNPs remaining for analysis	λ before correction	λ after correction	P value
1	PCAG	Controls	24	37	71,687	11,749	96,259	1.96	1.04	1.7×10−4
2	PCAG + PLA	Controls	81	37	71,633	10,166	97,282	1.34	1.00	7.7×10−5
3	PLA	Controls	57	37	73,253	11,935	94,781	1.28	1.01	1.0×10−4
4	PCAG	PLA	24	57	72,150	9,858	96,868	2.23	1.01	3.8×10-7*
5	PCAG	Controls + PLA	24	94	71,633	10,166	97,282	2.25	0.96	1.4×10-7*

No.=number, λ=lambda (genomic inflation factor), * denotes significant association (Bonferonni)

Definition of associated loci

Associated loci were defined based on pairwise linkage disequilibrium (LD) estimates (conservative r2 ≥0.5) of the SNPs using the Tagger program embedded in Haploview [36,37]. A list of SNPs in LD with the most strongly associated GWAS SNP for each locus was generated which could denote a region tagged by the GWAS association signal. As for the GWAS, SNPs were excluded from analysis if they had a MAF <0.05 or a call rate of <97%.

Gene expression using RNA-Seq

Eyes were removed from five BHs with PCAG on welfare grounds and from four non-BH dogs with normal eyes that had been euthanized via intravenous overdose of pentobarbitone for reasons unrelated to this study. Following enucleation, each globe was transected along the sagittal plane into two equal halves. Each half globe was immersed in 15 ml RNAlater solution (Ambion; Foster City, CA) before being stored at −80 °C until later use as previously described [38,39]. RNA was extracted from selected tissue samples using the Qiagen RNAeasy Midi Kit (Qiagen; Hilden, Germany) according to the manufacturer’s instructions. Following thawing of each half globe, a section of the ICA was dissected under an operating microscope with microsurgical instrumentation aiming for approximately 200–250 mg tissue. The posterior boundary of the dissected tissue was the posterior pars plana, and the anterior boundary was the corneoscleral limbus. No lens tissue was included in the dissection, but the sclera was included. Adherent conjunctiva and episcleral tissues were removed. Total RNA (20 ng–1 µg per sample, at a concentration of ≥20 ng/μl with RNA integrity number (RIN) ≥8.0 and optical density (OD)260/280 ≥1.8) was submitted to Otogenetics Corporation (Atlanta, GA) for library preparation and sequencing: 1–2 μg of cDNA was generated using the Clontech Smart cDNA kit (Clontech Laboratories; Mountain View, CA) from 100 ng of total RNA. cDNA was fragmented using Covaris (Covaris; Woburn, MA) and profiled using an Agilent Bioanalyzer 2100 (Agilent; Santa Clara, CA). Libraries were prepared using NEBNext reagents (Catalog No. 634,925. New England Biolabs; Ipswich, MA): mRNAs were purified using Poly(A) selection from the total RNA sample and then fragmented. cDNA was then synthesized using random priming, followed by end repair, phosphorylation, A-tailing, adaptor ligation, and finally, PCR amplification. The quality, quantity, and size distribution of the Illumina libraries were determined using an Agilent Bioanalyzer 2100. The libraries were then loaded on an Illumina HiSeq2500 for clustering and sequencing according to the standard operation. Paired-end 90–100 nucleotide reads were generated, and data quality was assessed using FASTQC (Babraham Institute; Cambridge, UK). After optimum quality control (QC) results were achieved, the samples were analyzed. Sequence reads (minimum of 50 M per sample) were mapped to the CanFam3.1 reference genome, and the read count values of gene expression (fragments per kilobase million, FPKM) were calculated. The differential expression of the FPKM (the read counts) between the PCAG cases and controls was calculated as the base 2 log of the fold change of controls/PCAG cases (log2 (the FPKM of the PCAG cases/the FPKM of the controls)). The test statistic was used to compute the statistical significance of the observed change in the FPKM. The p value represents the uncorrected p value of the test statistic, and the Q value represents the false discovery rate–adjusted p value of the test statistic. For genes shown to be differentially expressed, their functions and associated phenotypes (for potential candidacy for canine PCAG) were assessed using VarElect GeneCards®.

Results

NEB variant genotyping

One hundred fifty-seven BHs (ten PCAG cases, 52 PLA cases, and 95 controls) were homozygous for the reported NEB variant (G/G), and only one control was heterozygous (A/G). In the 83 non-BH dogs, 59 were homozygous for the variant (G/G), 15 were heterozygous (A/G), and only nine were homozygous for the reference allele (A/A). Thus, as well as not being associated with PLA or PCAG, it appeared that the SNP at this location is a common polymorphism, and the reported variant may actually be more common in European dogs than the reference allele.

Genome-wide association studies

The results of the five GWASs, including case and control designation, filtering steps, the genomic inflation factor, and the p value of the top SNP from each analysis can be found in Table 2. Statistically significant SNP associations were found in GWAS Analyses 4 and 5. In Analysis 4 (Cases: PCAG cases; Controls: PLA cases), a SNP on canine chromosome 24 (BICF2G630505097; CanFam3.1 chr24:17,381,226; p = 3.8 × 10−7) reached the threshold of genome-wide statistical association (Figure 1). In Analysis 5 (Cases: PCAG cases; Controls: PLA cases and control dogs), an additional SNP in the same region increased the strength of the statistical association further (BICF2P544799; CanFam3.1 chr24:18,739,902; p = 1.4 × 10−7; Figure 2).

Figure 1

Manhattan plot for GWAS Analysis 4 (Cases: PCAG cases; Controls: PLA cases) in the European BH. The horizontal red line denotes the threshold for genome-wide statistical association.

Figure 2

Manhattan plot for GWAS Analysis 5 (Cases: PCAG cases; Controls: PLA cases and control dogs) in the European BH. The horizontal red line denotes the threshold for genome-wide statistical association.

GWAS conditioning on the most strongly associated SNP was performed to assess whether any additional loci were associated with PCAG. To achieve this, SNP BICF2P544799 was introduced as a covariate in the mixed model in Analysis 5. This revealed an additional SNP at a distinct chromosomal location on canine chromosome 37 to reach the threshold of genome-wide statistical association (BICF2P928441; CanFam3.1 chr37: 24,747,131; p = 1.9 × 10−7; Figure 3). Conditional analysis using both associated SNPs as covariates in the model revealed no further association signals.

Figure 3

Manhattan plot for GWAS analysis 5 (Cases: PCAG cases; Controls: PLA cases and control dogs) in the European BH using BICF2P544799 as a covariate. The horizontal red line denotes the threshold for genome-wide statistical association.

The interrelationships between the three strongest associated SNPs and PCAG risk were further assessed in the individual SNP data sets using logistic regression. This was done for the original GWAS set and with the addition of PCAG cases from the USA. All three SNPs were strongly associated with PCAG, and the strength of association increased following the addition of the American dogs with PCAG to the analyses (Table 3). A combined analysis of the two SNPs and the risk of PCAG in all 153 dogs (59 cases and 94 controls defined as Table 3) using an additive (per allele) model resulted in an odds ratio of 19.30 (95% confidence interval 7.62–48.85; p = 1.4 × 10−24).

Table 3

Association of GWAS top SNPs with PCAG in the BH (Analysis 5 (see Table 2)

Top SNP from GWAS	USA PCAG cases included?	No. PCAG cases	No. PLA cases and controls	OR	Lower 95% CI	Upper 95% CI	P-value
chr24:17381226 BICF2G630505097	No	24	94	13.73	5.06	37.26	4.2×10−11
	Yes	59	94	18.84	7.97	44.51	1.6×10−21
chr24:18739902 BICF2P544799	No	24	94	28.35	8.50	94.58	1.6×10−11
	Yes	59	94	15.81	6.03	41.45	1.4×10−12
chr37: 24,747,131 BICF2P928441	No	24	94	18.03	4.44	73.28	1.1×10−6
	Yes	59	94	14.76	4.19	51.96	3.6×10−8

No.=number, OR=odds ratio, CI=confidence interval

A critical region spanning the top associated SNPs was defined based on pairwise LD of these top three SNPs with other SNPs on the respective chromosome in the GWAS data set. Although the top two chromosome 24 SNPs BICF2G630505097 and BICF2P544799 were correlated at an r2 of 0.79 in this study set, they were both used as tag SNPs to ensure coverage. From these analyses, the following chromosomal loci were defined for subsequent investigation: CanFam3.1 chr24:17,381,226–18,739,902 and chr37:24,747,131–24,958,250.

Identification of candidate genes from RNA-Seq

Whole transcriptome data derived from ICA tissues were used to compare differential gene expression between BH PCAG cases (n=5) and non-BH controls (n=4) in the extended loci determined from GWAS and LD analysis (Table 4). In the chromosome 24 locus, differential gene expression was present for five genes or loci: SIGLEC1 (Gene ID 6614, OMIM 600751), C24H20orf194, SLC4A11 (Gene ID 83959, OMIM 610206), PROSAPIP1 (Gene ID 9762, OMIM 610484), and OXT (Gene ID 5020, OMIM 167050). In the chromosome 37 locus, differential expression was present for three genes: CXCR2 (Gene ID 3579, OMIM 146928), CXCR1 (Gene ID 3577, OMIM 146929) and ARPC2 (Gene ID 10109, OMIM 604224). A survey of the functions and associated phenotypes of these genes revealed no evidence of reported association with glaucoma. All genes, however, appear to be involved in inflammation or immunity (Table 5).

Table 4

Differential gene expression results for the BH PCAG loci.

Gene	CanFam3.1 coordinate	Status	FPKM Controls	FPKM PCAG cases	log2(fold change)	Test statistic	P-value	Q-value	Significant?
Locus chr24:17,381,226–18,739,902
LOC100855425	24:17377802–17459234	OK	3.66177	4.27275	0.222624	0.105002	0.90045	0.972499	no
RNF24	24:17377802–17459234	OK	13.5062	15.8018	0.226461	0.549112	0.2461	0.564651	no
PAX1	24:1743907–1752578	NOTEST	0.751058	0.206085	−1.86569	0	1	1	no
PANK2	24:17463086–17489560	OK	19.9698	18.2754	−0.127921	−0.289163	0.55265	0.833548	no
MAVS	24:17511842–17529197	OK	29.2054	32.6664	0.161573	0.390852	0.42265	0.743935	no
AP5S1	24:17539786–17543605	OK	8.54101	8.52501	−0.00270593	−0.00393991	0.99425	0.998863	no
CDC25B	24:17552205–17561823	OK	12.1744	14.8	0.28175	0.556191	0.2447	0.563837	no
CENPB	24:17569396–17572311	OK	31.9742	38.916	0.283456	0.799702	0.1328	0.393219	no
SPEF1	24:17572374–17578571	OK	1.5338	1.92433	0.327243	0.400218	0.3999	0.724145	no
C24H20orf27	24:17586019–17600147	OK	22.4481	29.0792	0.373394	0.779776	0.1107	0.352522	no
HSPA12B	24:17600534–17619889	OK	6.57984	9.70452	0.560605	1.14769	0.0296	0.142427	no
SIGLEC1	24:17636164–17667080	OK	2.18752	12.5593	2.52139	4.87528	5.00E-05	0.0008345	yes
ADAM33	24:17668462–17682217	OK	10.7869	8.23265	−0.389845	−0.829875	0.0807	0.285699	no
GFRA4	24:17686816–17690747	OK	3.16594	2.14218	−0.563553	−0.711984	0.2384	0.555935	no
ATRN	24:17698237–17859816	OK	18.7085	14.7908	−0.338994	−0.812765	0.08165	0.288073	no
C24H20orf194	24:17908366–18050002	OK	14.0202	9.25396	−0.599364	−1.46229	0.00165	0.0153647	yes
SLC4A11	24:18057411–18068425	OK	28.8724	6.6876	−2.11013	−4.36281	5.00E-05	0.0008345	yes
ITPA	24:18070979–18082187	OK	11.9247	16.3186	0.452557	0.814939	0.0945	0.316902	no
DDRGK1	24:18085722–18098868	OK	17.4284	15.1717	−0.200057	−0.450069	0.33115	0.659854	no
LOC102157332	24:18111855–18122295	NOTEST	0.531426	0.412288	−0.366218	0	1	1	no
PROSAPIP1	24:18111855–18122295	OK	20.0459	12.6993	−0.658565	−1.63743	0.00105	0.0107012	yes
FASTKD5	24:18124262–18166506	OK	13.9635	14.315	0.0358696	0.0772293	0.86905	0.965583	no
CST11	24:181590–184213	NOTEST	0	0	0	0	1	1	no
AVP	24:18183056–18184827	NOTEST	0.182349	0.0787568	−1.21123	0	1	1	no
OXT	24:18193380–18194236	OK	18.6151	5.43912	−1.77503	−2.94744	0.0006	0.0068551	yes
MRPS26	24:18215054–18217105	OK	24.9994	26.2617	0.0710653	0.176449	0.7648	0.934576	no
LOC102151131	24:18217903–18222755	NOTEST	0.0563875	0.0232463	−1.27838	0	1	1	no
LOC102154818	24:18224354–18395654	NOTEST	0.867838	0.884208	0.0269605	0	1	1	no
PTPRA	24:18224354–18395654	OK	30.9735	28.6433	−0.112839	−0.288467	0.5554	0.835146	no
VPS16	24:18401450–18424606	OK	11.5555	11.6659	0.0137143	0.0286225	0.95285	0.989023	no
PCED1A	24:18424775–18429698	OK	18.7163	16.9792	−0.140525	−0.340822	0.5442	0.828421	no
LOC102155387	24:18442763–18443512	NOTEST	0.12029	0.206892	0.782359	0	1	1	no
TMEM239	24:18443575–18445637	NOTEST	0.308423	0.133061	−1.21283	0	1	1	no
C24H20orf141	24:18445693–18446873	NOTEST	0.202272	0.317907	0.652307	0	1	1	no
CPXM1	24:18453867–18460350	OK	8.40695	13.7502	0.7098	1.29913	0.0077	0.0506978	no
EBF4	24:18470596–18528305	OK	1.69066	1.6627	−0.0240522	−0.0278686	0.96215	0.991001	no
IDH3B	24:18556489–18561401	OK	46.7683	38.3285	−0.287114	−0.659815	0.1486	0.420704	no
NOP56	24:18561409–18566780	OK	28.0915	30.953	0.139945	0.325807	0.48615	0.78904	no
TMC2	24:18577342–18642204	NOTEST	0.308586	0.147869	−1.06136	0	1	1	no
ZNF343	24:18661008–18665685	NOTEST	0.245175	0.133892	−0.872746	0	1	1	no
SNRPB	24:18674022–18683315	OK	62.0096	61.9386	−0.00165212	−0.004023	0.99365	0.998542	no
TGM6	24:18700241–18721827	NOTEST	0.0171946	0	-	0	1	1	no
TGM3	24:18753637–18793790	NOTEST	0.0342808	0.138228	2.01158	0	1	1	no
Locus chr37:24,747,131–24,958,250
RUFY4	37:24809358–24830254	OK	1.71109	1.18547	−0.52945	−0.51443	0.27475	0.599772	no
CXCR2	37:24831509–24846517	OK	3.35508	1.02764	−1.70702	−1.59836	0.00145	0.013826	yes
CXCR1	37:24861913–24866062	OK	1.82409	0.616257	−1.56557	−1.58159	0.0041	0.030963	yes
ARPC2	37:24914388–24936198	OK	84.4912	137.67	0.704343	1.97128	0.00035	0.004435	yes
GPBAR1	37:24940929–24943399	OK	1.26415	2.05007	0.697505	0.919835	0.14675	0.417894	no
AAMP	37:24943657–24948585	OK	40.9075	43.1513	0.077039	0.181187	0.70665	0.913414	no
PNKD	37:24948682–25008542	OK	6.30845	7.78033	0.302547	0.202254	0.70135	0.911242	no

OK=test successful, NOTEST=not enough alignments for testing

Table 5

Details of differentially expressed genes, their functions and associated disorders and phenotypes (assessed using VarElect GeneCards®)

Gene	Gene full name	Function	Associated disorders/phenotypes
SIGLEC1	Sialoadhesin	Endocytic receptor mediating clathrin dependent endocytosis	Glomerulonephritis Sclerosis X-linked intellectual disability
C20orf194	Chromosome 20 open reading frame 194	May act as an effector for ARL3	HIV-1 Hepatitis
SLC4A11	Solute Carrier Family 4 Member 11	Transporter which plays an important role in sodium-mediated fluid transport in different organs.	Corneal dystrophy HIV-1 Hepatitis
PROSAPIP1	Proline Rich Synapse Associated Protein Interacting Protein 1	May be involved in promoting the maturation of dendritic spines	Hepatitis
OXT	Oxytocin	Contraction of smooth muscle of uterus and mammary gland	Persistent genital arousal Endometritis Inhibited male orgasm Epignathus Chorioamnionitis Parturition Lactation
CXCR1	Chemokine (CXC) Receptor 1	Neutrophil activation, neutrophil count	HIV-1 Pyelonephritis & urinary tract infections Idiopathic anterior uveitis
CXCR2	Chemokine (CXC) Receptor 2	Neutrophil activation, neutrophil count	Congenital neutropaenia Neutrophil migration Pyelonephritis Septicaemia Granulocytic anaplasmosis
ARPC2	Actin Related Protein 2/3 Complex Subunit 2	Actin filament assembly. Platelet, reticuloycte and neutrophil count	Platelet, reticuloycte and neutrophil count

Discussion

Canine PCAG is considered complex, likely involving multiple genetic and environmental factors [40]. GWASs have been used extensively in the human field over the last two decades to improve the understanding of the genetic basis of many complex traits by revealing associated susceptibility loci and have provided invaluable insights into the allelic architecture of many multifactorial traits [41]. In humans, GWASs of thousands of individuals are usually required to find significant associations with complex diseases [42,43]. More recently, findings from human GWASs of complex disease have supported the original hypothesis that their genetic architecture is mainly comprised of common variants with modest or intermediate effect, as well as a smaller subset of high-penetrant (familial) low-frequency variants of high effect [44]. In the dog, however, many fewer canine samples are required to detect a statistically significant association compared to human studies, owing to the more limited genomic architecture and haplotype structure of domestic dog breeds [22,23,45-47]. In the present GWAS, a locus on canine chromosome 24 was found to be associated with PCAG relative to PLA in Analysis 4 using a total of only 81 dogs (24 PCAG cases and 57 PLA cases), and the level of statistical association was further increased when the controls were combined with the PLA cases in Analysis 5 (a total of 118 dogs). An additional locus was found on canine chromosome 37 following conditional analysis. Owing to our employment of correction for population stratification and the stringent threshold for a statistically significant association, it is unlikely that these are false positive associations [48]; in addition, the strength of the associations was augmented by the addition of the American BH PCAG cases. This was particularly marked for BICF2G630505097 which was also the most strongly associated SNP in GWAS Analysis 4 (Cases: PCAG cases; Controls: PLA cases). This SNP is located in an intergenic region (chr24:17,381,226) with RNF24 (Gene ID 11237, OMIM 612489; chr24:17,423,245–17,454,478) and PANK2 (Gene ID 80025, OMIM 606157; chr24:17,463,603–17,489,677) the nearest upstream genes. These genes may prove to be candidate genes for PCAG as a previous study of a glaucoma-related trait in humans, optic nerve head morphology, revealed a statistically significant linkage signal at the RNF24/PANK2 locus [49].

The present GWAS results did not corroborate those of previous studies in the BH. A GWAS of PCAG in American BHs reported two regions were associated with PCAG [22]. These susceptibility loci were situated on canine chromosomes 14 and 24 (chr14:19,911,001–20,161,358 and chr24:43,091,222–43,595,979) and contained the collagen gene COL1A2 and oncogene RAB22A, respectively. The inability to reproduce these findings in the present GWAS might be related to population differences in the frequency of the risk loci between European BHs and American BHs. Another study in American BHs by the same investigators used linkage analysis followed by exome-sequencing analysis to identify an SNP in NEB to segregate with PCAG [21]. The present study, however, did not find any association between this variant and PCAG, and the multibreed genotyping suggested that this variant is a common polymorphism. As the resequencing study that identified this NEB variant captured only the exome, it is possible that the variant is tagging a true functional variant underlying the linkage study results [21]. However, although a recessive mode of inheritance was assumed, the logarithm (base 10) of odds (LOD) score was considerably lower than those achieved for other recessive diseases that have been mapped in a similar way [21,50-52]. Therefore, it is possible that this result was a chance finding; the magnitude of associations identified in the present GWAS, which used stringent correction for multiple testing and a set of American BHs that validated the results, does not support PCAG in the BH to be a recessive disease from the sample set we analyzed; PCAG is more likely to be of complex etiology.

It was interesting that a statistically significant association of SNPs with PCAG was found in Analysis 4 (Cases: PCAG cases; Controls: PLA cases) and not in Analysis 2 (Cases: PCAG and PLA cases; Controls: control dogs). PLA is considered a consistent risk factor for, and therefore, likely to be on the causal pathway to, canine PCAG [5,9,15,53]. Thus, we considered it probable that there are shared genetic factors between PCAG and PLA cases, and a statistically significant association with disease (PCAG or PLA or both) would more likely be discovered when PCAG is analyzed combined with PLA cases against controls owing to the greater number of cases, and thus, power, of such analysis. Instead, these studies revealed two loci that are statistically significantly more common in PCAG cases than in PLA cases and controls. This is not a completely unexpected finding, however, as although PLA is a risk factor for PCAG, only a minority of dogs with PLA develop PCAG. This study supports the theory that PLA is required for but is not sufficient for PCAG. Thus, it is likely that the two chromosomal loci identified are involved in the trigger of the progression from PLA to PCAG. This finding will influence our ongoing investigations of the allelic architecture of PCAG in the BH. Future investigations to identify causal variants for PCAG in these loci will likely utilize a next-generation whole genome sequencing (WGS) approach using multiple PCAG and PLA cases and subsequent genotyping of segregating variants in extended case-control sets. For complex disease, such as canine PCAG, a relatively large sample set is likely to be needed. Furthermore, the efficiency of a WGS approach will likely be improved by selecting dogs based on a homozygous genotype for the most strongly associated GWAS SNPs.

The main challenge of GWASs in investigating complex disease is to pinpoint possible causal variants underlying association signals as the majority of GWAS hits are in non-coding or intergenic regions because complex disease is often caused by disturbance to biologic networks, not by isolated genes or proteins [54]. Regulatory SNPs can influence gene expression through several mechanisms that include the three-dimensional organization of the genome, RNA splicing, transcription factor binding, DNA methylation, and long non-coding RNAs [55,56]. In an attempt to overcome these limitations, we used RNA-Seq. We first used the most strongly associated SNPs from the present GWASs to define chromosomal loci based on their ability to be tagged by other nearby SNPs. We then used RNA-Seq to compare gene expression between affected and unaffected dogs within these loci. Eight genes were found to be differentially expressed, and none of these genes are obvious candidates for canine PCAG. Instead, most of the genes identified have functions related to inflammation and immunity. This finding was not completely unexpected as the results of several other studies have implicated the role of inflammation and immunity in the pathogenesis of canine and human glaucoma [57-62]. Thus, the present results appear to corroborate these findings. However, intraocular inflammation may have occurred secondary to the onset of PCAG in the cases we used for RNA-Seq. Inflammation is a common histological finding in eyes affected by primary and secondary glaucoma, and there is no published evidence of inflammation being present in eyes considered at risk of PCAG (i.e., those with severe PLA) [63-65].

There are several study limitations that should be discussed. This GWAS included many fewer individuals compared to studies of complex disease in humans. In the dog, much more modest sample numbers are generally used mainly as a result of limitations in sample availability, robust phenotypic characterization, and funding. The limitations in sample availability in the present study also, unfortunately, meant it was not possible to repeat the GWAS with an additional replication set of samples. There were also limitations in the RNA-Seq investigations. The lack of ocular tissue from BHs with normal eyes meant that tissues from other breeds had to be used that represent sub-optimal control samples. Furthermore, the sequencing of RNA from multiple cell types that make up the ICA and neighboring structures could have led to a misleading pattern of differential gene expression than is actually occurring in the specific region of the ICA that is directly implicated in the pathogenesis of PCAG.

In conclusion, we successfully used GWASs to identify two novel loci associated with canine PCAG. Comparison of the gene expression profiles implicated genes involved in inflammation and immunity as involved in pathogenesis. WGS, involving multiple dogs, is likely required to elucidate candidate causal variants for canine PCAG that underlie the present GWAS findings.