An Application of NGS for WDR36 Gene in Taiwanese Patients with Juvenile-Onset Open-Angle Glaucoma.

Primary open-angle glaucoma (POAG) is one of the most important disease in ophthalmology with high prevalence and risk of irreversible blindness. If diagnosed before the age of 35, it is usually categorized as juvenile open-angle glaucoma (JOAG). The WDR36 gene is reckoned as one of the major causative genes of POAG, and had been studied to be related to the pathogenesis of POAG in the literature. We have selected 61 JOAG patients and 61 JOAG-free individuals, and by next-generation sequencing method, the WDR36 gene of the subjects were analyzed. We identified 26 variations exclusively in JOAG group. Among these 26 variations, there were 3 noteworthy variations. First, a novel variation c.460-650A>G was found in our study which might cause premature termination of splicing of the conserved domain in WDR36; second, c.1494+1111G>T (rs13178997) had significantly different frequency in our JOAG patients compared to the reference frequency on NCBI; third, a variation c.710+30C>T (rs10038177) was found in our study, which had already been reported to be related to high-pressure glaucoma. We offer the profile of WDR36 in JOAG in Taiwan population, and we suggest that WDR36 gene is involved in the pathogenesis of JOAG as a subordinate modifier gene.

Introduction

Primary open-angle glaucoma (POAG) is a major type of glaucoma worldwide, a progressive optic neuropathy typically featuring increased intraocular pressure (IOP) and degeneration of retinal ganglion cells and fibers, with high risk of irreversible blindness 1, 2. In patients younger than 35 years old, open-angle glaucoma is often diagnosed as juvenile open-angle glaucoma (JOAG) 3. The etiology of JOAG remains unknown, but it is believed that the pathogenesis of glaucoma is multifactorial in genetic, environmental, cardiovascular, intraocular pressure-related and aging aspects, with the genetic etiology of glaucoma mostly studied 4, 5.

POAG is a critical ophthalmic issue worldwide, with over a half of all POAG patients in Asia. By 2020, it was estimated that over 73 million people would be affected by POAG 6. In a cross-sectional survey of Chinese population over 40 years old in Singapore, 3.2% prevalence of POAG was observed and POAG was reported to be the leading cause of blindness 7. In another study, glaucoma is also shown to be the leading cause of low vision and blindness in Taiwan 8. Thus, POAG is increasingly concerned as a critical issue in ophthalmology due to its high prevalence and severe optic sequelae.

POAG was heterogeneously associated with at least 20 genetic loci, and among these loci, 11 loci were denominated as GLC1A to GLC1K, among which GLC1A, GLC1E, GLC1G were further identified as myocilin (MYOC), optineurin (OPTN), WD repeat domain 36 (WDR36) genes, respectively 2. Three loci, WDR36, MYOC and OPTN, have been identified as main causative genes, contributing less than 10% of all POAG etiology with WDR36 accounting for 5% alone 2, 5. WDR36 gene was identified within GLC1G locus on 5q22.1, associated with optic neuropathies in glaucoma. Although detailed pathogenesis of mutation of WDR36 remains unclear, genetic screening of WDR36 gene is believed to be useful for early diagnosis of presymptomatic open-angle glaucoma in high risk individuals 9.

In this study, we focus on WDR36 gene, investigating its role in JOAG, although other 6 genes have been reported to be associated with glaucoma, including myocilin (MYOC), optineurin (OPTN), neurotrophin 4 (NTF4), optic atrophy 1 (OPA1), cytochrome P450, family 1, subfamily B (CYP1B1), and latent transforming growth factor-beta binding protein 2 (LTBP2) 10.

The Next-Generation Sequencing (NGS), also named Massive Parallel Sequencing, is a high-throughput DNA sequencing technology that enables scientists to analyze whole genomes or transcriptomes in a more precise but rapid manner, making great advance in genetics or in other related fields 11. It is believed that more detailed genomic alteration could be detected precisely by NGS than by other conventional technologies, improving and renewing our understanding of genomic biology 12.

To the best of our knowledge, the causative genes OPTN and MYOC, bar WDR36, have been studied already in Taiwan 3, 13. In this study, we aimed to compare variants in WDR36 genes of JOAG patients with those of unrelated normal individuals in Taiwan, using NGS technology and comparative genetic analysis methods. Understanding the genetic correlation between JOAG and WDR36 may have benefits to early diagnosis or treatment of this disease.

Materials and Methods

Subject Selection

This part of subject selection is partially modified out of a previous study investigating the same set of 1210 patients for juvenile-onset open-angle glaucoma research 13. In addition to the original subjects, 10 new patients who had undergone identical clinical examinations, including IOP measurement, visual-field test, slit-lamp examination and fundus examination, were added to this study.

Identification of JOAG can be established as clinical presentations of intraocular pressure (IOP) > 22mmHg, a cup-to-disc ratio > 0.5 or one with an asymmetric appearance, a visual field loss characteristic of glaucomatous change, and an open angle width ranging from Shaffer grade II to IV without any other apparent secondary cause like trauma or surgery. In total, 61 unrelated individuals screened were diagnosed as JOAG.

Randomly selected, another 61 individuals over 50 years old comprised the control group. Subjects of the control group also received aforementioned examinations, and glaucoma-free outcomes of the subjects were confirmed. All of 122 individuals, including 61 patients and 61 controls, belong to Han ethic origin in the study. The study protocol was approved by the Institutional Review Board of the Kuo General Hospital, conforming to the World Medical Association's Declaration of Helsinki (2000). All of the participants provided signed informed consent to participate in the study after the details of the study had been explained to them.

Detection of the mutations of the WDR36 gene

Genomic DNA samples were extracted from 10ml peripheral blood collected from every participating individual and were purified by using a Gentra DNA Blood Kit (Gentra Systems, Inc. Minneapolis, USA) according to the manufacturer's instructions. Gel electrophoresis and spectrophotometry were applied to determine the quality and quantity of the purified DNA, respectively. Mutations in the 23 exons and flanking introns of the WDR36 gene were screened by SOLiDTM NGS sequencing (Ensembl Transcript ID. ENST00000323652).

Long range polymerase chain reaction (PCR) was applied for enrichment of WDR36 gene, intragenic primer sequences shown in Table 1. Cycling profile of long range PCR was conducted as follows: one cycle at 94℃ for 2 min, followed by 10 cycles at 94℃ for 10 s, at 55℃ for 30 s, and at 68℃ for 14 min, followed by 30 cycles at 94℃ for 10 s, at 55℃ for 30 s, and at 68℃ for 14 min (+20 s every cycle), followed by a final extension step at 72℃ for 5 min, according to SequalPrepTM Long PCR Kit protocol (Invitrogen, USA). According to the manufacturer's instructions, PCR products were cleaned up using AMPure Reagent beads (Agencourt, Beverly, MA). Amplicons were pooled at an equimolar ratio, and a sequencing library was made up of 5 micrograms of DNA.

Table 1

Oligonucleotide Primer Pairs for PCR

Name	Oligonucleotide
WDR36 E1F	GGGTCTGTGAGTGAATCCCTCTGTC
WDR36 E5R	CCACAACTGCAGGCTTCCTTG
WDR36 E5F	CAAGGAAGCCTGCAGTTGTGG
WDR36 I10R	GCATTGATGAAACTTCCTCCAGTTG
WDR36 I10F	CAACTGGAGGAAGTTTCATCAATGC
WDR36 I15R	TGATCGCATCAACTCCCTGAAA
WDR36 I15F	TTTCAGGGAGTTGATGCGATCA
WDR36 I16R	CAACAGGGAGGAAACAGGAGGA
WDR36 I16F	TCCTCCTGTTTCCTCCCTGTTG
WDR36 I19R	AGCACCCTTGCCGATAAGGC
WDR36 I19F	GCCTTATCGGCAAGGGTGCT
WDR36 E23R	CTCCCAGAATGTCAAAAAGTG

Generally, SOLiD Sequencing Analysis was then performed with the amplified genomic DNA, with the first step being preparation of fragment library, followed orderly by preparation of templated bead, SOLiD analyzer operation (sequencing) and data analysis.

Splicing site prediction by neural network

Both normal and variant WDR36 genomic sequence of the intron 3 were analyzed by using a neural network prediction system of Berkeley Drosophila Genome Project to predict whether there were novel donor or acceptor sites 14.

Result

We screened the WDR36 (GenBank NC_000005) of 61 JOAG patients and 61 normal Taiwanese individuals, including intronic sequence, promoter region and coding region, by PCR amplification and targeting genome sequencing analysis. In general, 27 variants were identified from 61 JOAG patients, whereas 6 variants were identified from the control group. Only one variant c.1964+2910T>C was found in both JOAG patients and normal controls, but no significant difference in frequency was shown between the two groups, with the p value calculated by χ2 test being 0.862, suggesting that this variant has no direct relation to JOAG pathogenesis (Table 2).

Table 2

c.1964+2910T>C Variant Shared by JOAG Patients and Normal Individuals

	Reference Allele Count	Novel Allele Count	Novel Allele Frequency
JOAG	2	3	0.600
Control	27	37	0.578
p value	0.862

Among the 26 variants only found in the JOAG patients, 20 of them were established SNPs, including rs1971050, rs1993465, rs13153937, rs199945131, rs201148106, rs6859041, rs1379298, rs10038177, rs6865932, rs10045255, rs10043631, rs10038058, rs13178997, rs4957924, rs43203, rs6594498, rs7722241, rs7702774, rs4530809 and rs12520738 (Table 3). Among those 20 SNPs, we found that only the rs13178997 (c.1494+1111G>T) had a significant difference in frequency between JOAG patients and reference from NCBI database (p<0.02), while there was another SNP, rs10038177 (c.710+30C>T), which had been reported to be a risk factor for high-tension glaucoma 15. The remaining 18 SNPs had no remarkable findings by comparing to the corresponding frequencies on NCBI database.

Table 3

Variants Found Only in Normal Individuals

HGVS Name	Allele Count		Coordinate	Chromosome Position(GRCh38.p2)	SNP
c.1262-223 G>A	G:0.71	1589	13664	111105834	rs2034896
	A:0.29	647
c.1964+2843 A>G	A:0.34	11	23825	111115996	-
	G:0.66	21
c.1964+2865 G>C	G:0.71	15	23847	111116018	-
	C:0.29	6
c.1965-2867 A>G	A:0.80	156	23875	111116146	-
	G:0.20	40
c.1964+2918 C>T	C:0.16	9	23900	111116071	-
	T:0.84	46

As for the other 6 novel variants, including c.460-650A>G, c.1495-807T>G, c.1965-2849A>G, c.1965-2833C>T, c.1965-2812A>G and c.1965-2786C>G were only found in the JOAG patients, none of them were established SNPs (Table 4). In order to study the influence of the 6 novel variants on mRNA transcription, we utilized a neural network prediction system to predict whether there is any splicing alteration, using a cutoff of minimal score by 0.4. The variant c.460-650A>G was the only one of the six that creates a novel acceptor site AG at c.460-644_643AG, located within intron 3 (Figure 1). The neural network value of the novel acceptor splicing site score was 0.5, crossing the threshold cutoff. The coding DNA sequence (CDS) of WDR36 was translated into 952 amino acids under normal condition, whereas in c.460-650A>G variants, the mRNA were added with 642 nucleotides within intron 3, leading to a premature termination of translation at residue 198 (TGA) (Figure 1). Therefore, we suggest that this variant c.460-650A>G can cause the production of truncated proteins and further impair the function of the proteins by altering the splicing pattern, which might be one of the etiologies of JOAG.

Table 4

Variants Found in JOAG Patients

HGVS Name	Coordinate	Position	Chromosome Position(GRCh38.p2)	SNP	Allele Count and Frequency		Reference Allele Frequency	p value
c.330+925C>T	1372	intron1	111093543	rs1971050	11547	C: 0.55	0.389	0.208
					9448	T: 0.45	0.611
c.459+221A>G	5229	intron3	111097400	rs1993465	8871	A: 0.69	0.667	0.861
					4023	G: 0.31	0.333
c.460-113G>A	6438	intron3	111098609	rs13153937	14732	G: 0.75	0.686	0.578
					4907	A: 0.25	0.314
c.577+624G>T	7292	intron4	111099463	rs199945131	15	G: 0.01	-	-
					1112	T: 0.99	-
c.577+637T>G	7305	intron4	111099476	rs201148106	983	T: 0.55	-	-
					793	G: 0.45	-
c.577+694G>A	7362	intron4	111099533	rs6859041	8768	G: 0.72	0.667	0.670
					3455	A: 0.28	0.333
c.578-561T>C	7857	intron4	111100028	rs1379298	11861	T: 0.71	0.617	0.458
					4918	C: 0.29	0.383
c.710+30C>T	8581	intron5	111100751	rs10038177	2932	C: 0.77	0.637	0.268
					892	T: 0.23	0.363
c.710+432G>C	8983	intron5	111101153	rs6865932	12948	G: 0.69	0.663	0.809
					5781	C: 0.31	0.338
c.765+259A>G	10488	intron6	111102658	rs10045255	8709	A: 0.70	0.663	0.765
					3757	G: 0.30	0.337
c.1494+90C>T	15359	intron12	111107529	rs10043631	5403	C: 0.76	0.617	0.221
					1688	T: 0.24	0.383
c.1494+143A>G	15412	intron12	111107582	rs10038058	4633	A: 0.65	0.651	0.976
					2524	G: 0.35	0.349
c.1494+1111G>T	16380	intron12	111108550	rs13178997	6059	G: 0.75	0.434	0.011
					1995	T: 0.25	0.566
c.2316+256C>T	29226	intron19	111121397	rs4957924	4948	C: 0.54	0.384	0.214
					4172	T: 0.46	0.616
c.2316+1144G>A	30114	intron19	111122285	rs43203	5569	G: 0.61	0.622	0.901
					3613	A: 0.39	0.378
c.2317-1152G>A	30482	intron19	111122653	rs6594498	11061	G: 0.78	0.617	0.158
					3069	A: 0.22	0.383
c.2437-559T>A	32878	intron21	111125049	rs7722241	7100	T: 0.73	0.711	0.885
					2657	A: 0.27	0.289
c.2437-455G>T	32982	intron21	111125153	rs7702774	4397	G: 0.70	0.705	0.961
					1891	T: 0.30	0.295
c.2707-202A>G	34361	intron22	111126532	rs4530809	3968	A: 0.71	0.711	0.989
					1626	G: 0.29	0.289
c.*1045T>G	35757	exon23	111127928	rs12520738	585	T: 0.64	-	-
					323	G: 0.36	-
c.460-650A>G	5901	intron3	111098072	-	773	A: 0.36	-	-
					1391	G: 0.64	-
c.1495-807T>G	17211	intron12	111109382	-	815	T: 0.64	-	-
					463	G: 0.36	-
c.1965-2849A>G	23993	intron16	111116164	-	2	A: 0.10	-	-
					18	G: 0.90	-
c.1965-2833C>T	24009	intron16	111116180	-	249	C: 0.54	-	-
					211	T: 0.46	-
c.1965-2812A>G	24030	intron16	111116201	-	234	A: 0.47	-	-
					263	G: 0.53	-
c.1965-2786C>G	24056	intron16	111116227	-	1	C: 0.02	-	-
					50	G: 0.98	-
c.1964+2910T>C	23892	intron16	111116063	-	2	T: 0.40	-	-
					3	C: 0.60	-

Figure 1

Spliced site prediction result of variation c.460-650A>G by a neural network prediction system. A: In WDR36, variant c.460-650A>G compared to normal transcript created a novel acceptor site c.460-644_643AG (star) with the score of neural network value increased from 0.34 to 0.50. B: In variants c.460-650A>G of WDR36, a 642-nucleotide fragment would be spliced out, leaving 38 residual nucleotides between c.459 and c.460. To avoid redundancy, 542 nucleotides within the spliced fragment are abbreviated. C: Comparison of predicted protein sequence of normal and variant of WDR36 showed substantially shortened amino acid sequence with a premature termination.

Another 2 SNPs found in JOAG patients, c.710+30C>T (rs10038177), reported to be a risk factor of high-tension glaucoma, and c.1494+1111G>T (rs13178997), calculated to have significant difference in frequency between JOAG patients and the control group, were also analyzed by the neural network prediction system in the same way, which revealed no difference of splicing result in c.710+30C>T (rs10038177). A minimal alteration of increased score of an acceptor site from 0.81 to 0.91 was detected in c.1494+1111G>T (rs13178997) (Table 5).

Table 5

Splicing Site Prediction Results of WDR36 in Variants and Control Group

Variant	Novel acceptor site prediction	Score	Start	End
c.460-650A>G	ataaaggacctttAtgcacAGgtgtggtcagggttagggga	0.34	348	388
	ataaaggacctttGtgcacAGgtgtggtcagggttagggga	0.50
c.1494+1111G>T	ttctatcttacgctcccGtAGggacatttacctacatagtt	0.81	484	524
	ttctatcttacgctcccTtAGggacatttacctacatagtt	0.91

Discussion

The WDR36 gene was recognized as one of the major causative genes of POAG. Similarly, in a yeast model, WDR36 sequence variant made a functional defect in rRNA processing that may have contributed to the development of POAG 16. However, the role of WDR36 in POAG is open to debate. Gallenberger et al. described in 2014 that, despite the decrease in the expression of Wdr36 mRNA, heterozygote Wdr36-deficient mice showed no significant difference from wild type individuals with regards to structures, intraocular pressure, susceptibility of retinal ganglion cells to excitotoxic damage, susceptibility of optic nerve to high IOP damage, and analysis of rRNA processing 17. Moreover, the role of WDR36 in POAG had been studied in different populations. Mookherjee and colleagues investigated 10 SNPs proposed to be associated with POAG in eastern India, and found that there was little correlation between POAG and SNPs in WDR36, except for rs10038177 (c.710+30C>T) which may be a risk factor for high tension glaucoma 15. In Italian population, it was described that WDR36 sequence variances play a minor role in the pathogenesis of POAG 18. Also, WDR36 was reported to be irrelevant to the development of glaucoma among German population 19.

Mookherjee and colleagues had studied 10 SNPs of WDR36 (rs1971050, rs1993465, rs13153937, rs10038177, rs11241095, rs10043631, rs10038058, rs10491424, rs17553936, and rs13186912) in East Indian population that were suspected to be contributive to POAG 15. In comparison, 4 SNPs (rs13186912, rs17553936, rs10491424 and rs11241095) were not identified in our patients of POAG, suggesting that these 4 variants were merely genetic polymorphisms between Indian and Taiwanese populations, with minimal role in pathogenesis of POAG.

In different ethnic groups, the frequency of WDR36 variant varies greatly, ranging from 2.9% to 5.6% 10. Moreover, different WDR36 mutations were identified among studies. For examples, Huang et al. identified several mutations within WDR36 which were mainly exons and were not found in our patients 10. Likewise, Bao Jian Fan et al. did not find significant mutations in coding exons or splicing junctions of WDR36 associated with POAG 1. We agree with the idea that, in consistency with Jonathan et al. 20, WDR36 may be just a modifier gene of POAG, playing a causative role in the pathogenesis of POAG to an uncertain extent.

WDR36, resembling yeast Utp21p in structure, with 14 WD40 repeats folded into two connected seven-bladed β-propellers, is known to be an essential nuclear protein localized in nucleolus that plays a critical role in the processing of 18S rRNA 16, 21. The variety and importance of function of WDR36 had been studied. Besides heart, brain and other non-ocular organs, expression of WDR36 were found in many of the ocular tissues, including lens, retina, optic nerve, and so on 9. With WDR36 gene deleted or depletion by RNA interference, mouse embryos die before implantation and, in human trabecular meshwork cells, apoptosis is triggered and the formation of small subunit ribosomal RNA delays 21. In zebrafish models, the function of Wdr36 protein had also been studied, and it was reported that Wdr36 interacted with the p53 pathway. In POAG, apoptosis causes the loss of retinal ganglion cells, which may be in accordance with the interaction between WDR36 and p53 pathway. However, without sufficient evidence, the WDR36 variant alone might not be a decisive factor in the pathogenesis of POAG, but its contributory role in the disease had be discussed 20. Moreover, WDR36 was investigated by cell experiments. It had been described that WDR36 acted as a new scaffold protein, interacting with β isoform of thromboxane A2 receptor, and enhancing the interaction between Gαq and PLCβ. Despite of the unclear association between the findings and POAG, understanding the functions of WDR36 can benefit future studies to elucidate its role in POAG 22.

According to NCBI database, within the amino acid sequence translated from WDR36 mRNA, it has been established that there are 2 types of conserved domain, including Utp21 and WD40. Utp21, known as an rRNA-processing protein, is one of the components of U3 snoRNP which plays an important role in 18S rRNA synthesis 16. WD40, typically featuring 11-24 residues of GH dipeptide at its N-terminus and 40 residues of WD dipeptide at its C-terminus, is able to form a special propeller-like platform structure and can be found within many proteins of various functions, including signal transduction, pre-mRNA processing and cytoskeleton system 23. The Utp21 conserved domain was found to be located at residue 728-948, whereas the residues associated with WD40 were identified at residue 155-420, 157-582, 332-650, 567-717 and 615-692, all of which were involved in regions behind 198, suggesting that, in c.460-650A>G variant, the premature termination at 198 residue could lead to the loss of these conserved domains, probably inducing the apoptosis of retinal ganglion cells which is one of the pathological features of POAG.

Conclusions

In summary, first, we have found a mutation c.460-650A>G within WDR36 which is likely to cause abnormal splicing and thus produce truncated proteins without conserved domains in Taiwanese population; further functional analysis of the mutations in this study is necessary in order to determine their exact role in POAG pathogenesis, which is also the limitation of this study. Second, by reviewing the literatures from different countries, we suggest that WDR36 is only a subordinate factor playing a role of modifier in POAG, without consistency among races or populations.