Genetic association between CDKN2B/CDKN2B-AS1 gene polymorphisms with primary glaucoma in a North Indian cohort: an original study and an updated meta-analysis

Background Variants in CDKN2B/CDKN2B-AS1 have been reported to modulate glaucoma risk in several GWAS across different populations. CDKN2B/CDKN2A encodes tumor suppressor proteins p16INK4A/p15INK4B which influences cell proliferation/senescence in RGCs, the degeneration of which is a risk factor for glaucoma. CDKN2B-AS1 codes a long non-coding RNA in antisense direction and is involved in influencing nearby CDKN2A/CDKN2B via regulatory mechanisms. Methods Current study investigated four SNPs (rs2157719, rs3217992, rs4977756, rs1063192) of aforementioned genes in a case–control study in a North Indian cohort. Genotyping was done with Taqman chemistry. In addition, an updated meta-analysis was performed. Results Two SNPs, rs3217992 and rs2157719 were found to be significantly associated with the disease. The frequency of ‘T’ allele of rs3217992 was significantly lower in cases (POAG/PACG) [p = 0.045; OR = 0.80(CI = 0.65–0.99) and p = 0.024; OR = 0.73(CI = 0.55–0.96)], respectively than in controls. Genetic model analysis revealed that TT + CT genotype confers 0.73-fold protection against POAG [p = 0.047; OR = 0.73(CI = 0.54–0.99)] and trend assumed additive model gives 0.53 times higher protection against PACG progression. However the association of rs3217992 with POAG and PACG did not remain significant after Bonferroni correction. For rs2157719, the ‘C’ allele was found to be less prevalent among cases (POAG/PACG) with respect to controls. Cochran Armitage trend test assuming additive model revealed 0.77 and 0.64-fold protection against POAG and PACG respectively. Bonferroni correction (pcorr = 0.003) was applied and the association of rs2157719 remained significant in PACG cases but not among POAG cases (p = 0.024). The ‘CC’ genotype also confers protection against primary glaucoma (POAG/PACG) among males and female subjects. The frequency rs1063192 and rs4977756 did not vary significantly among subjects, however the haplotype ‘CATA’ was found to be associated with increased glaucoma risk. An updated meta-analysis conducted on pooled studies on POAG cases and controls revealed significant association between rs1063192, rs2157719, rs4977756 and POAG except rs3217992. Conclusion The study concludes significant association between INK4 variants and primary glaucoma in the targeted North Indian Punjabi cohort. We believe that deep-sequencing of INK4 locus may help in identifying novel variants modifying susceptibility to glaucoma. Functional studies can further delineate the role of CDKN2B and CDKN2B-AS1 in primary glaucoma for therapeutic intervention.


Background
Glaucoma, a group of optic neuropathies is an outcome of degenerating retinal ganglion cells (RGCs), resulting in optic disc cupping and visual loss [1]. It is the second leading cause of irreversible blindness worldwide [2]. The disease has a complex etiology and a wide clinical spectrum. Amongst different clinical forms, primary open-angle glaucoma (POAG) and primary angle closure glaucoma (PACG) are the most common. Multiple factors contribute to RGCs degeneration including elevated intra-ocular pressure (IOP), oxidative stress in retinal microenvironment and abnormal glial activation [3]. Linkage scans till date have identified 27 genetic loci which may harbor glaucoma related genes [4], yet there exists a wide gap between the heritability estimates in glaucoma and causative genes. Genome-wide association studies (GWAS) have been instrumental in bridging this gap and have led to the identification of several novel genetic loci which might affect genetic risk to POAG and PACG [5]. One such region is INK4 locus at chromosome 9p21.3. Cyclin-dependent kinase inhibitor 2B (CDKN2B) and Cyclin-dependent kinase inhibitor 2A (CDKN2A) are two genes located adjacent to each other at INK4 locus in a stretch of about 80 kb. They encode tumor suppressor proteins p16 INK4A and p15 INK4B , respectively that inhibit cell cycle progression by forming complexes with cyclindependent kinase CDK4 or CDK6 [6]. CDKN2A also encodes for p14 ARF which is a splice-variant produced through an alternate reading frame [7]. Another gene located at INK4 locus is CDKN2B-AS1, also known as ANRIL (antisense non-coding RNA), encodes a long noncoding RNA in the antisense direction and is involved in influencing the nearby CDKN2A and CDKN2B genes via regulatory mechanisms [8]. Genetic variants in three of the genes at INK4 locus have been reported in several GWAS to be associated with glaucoma risk [8][9][10]. Initially thought to be associated with vertical cup disc ratio (VCDR) in a European cohort [11], and later on confirmed by Fan et al. [9], the association of CDKN2B/ CDKN2B-AS1 region with glaucoma was revalidated by several GWAS [8,[12][13][14]. Highly significant association of this region was observed with either the disease or its endophenotypes in many populations [8,[15][16][17][18][19][20][21][22] specifically with IOP. A study conducted by Gao and Jackobs in 2016, demonstrated increased vulnerability of RGCs in response to elevated IOP in mice homozygous for deletion in INK4 locus stressing that the altered expression of these genes might modulate apoptosis of RGCs, eventually contributing to glaucomatous visual field defect [23]. The respective variants in these genes might also affect their levels in eye via miRNA mediated regulation as demonstrated by Ghanbari and co-workers [24]. Two 3′UTR variants (rs3217992 and rs1063192) of CDKN2B enhance affinity for the binding of miR-138-3p and miR-323-5p, respectively to CDKN2B and thereby affect miRNA mediated regulation of this gene [24]. Further strengthening the candidacy of INK4 locus is the observation that p15 INK4B is a potential effector molecule for TGF-beta mediated cell-cycle arrest and thereby triggering axonal damage of optic nerve head (ONH) [6,25]. Replication studies from Indian subcontinent have so far failed to detect significant association between INK4 locus and glaucoma [26,27], yet GWAS data from across the world provide a strong rationale to further validate the region. Therefore, in the present study, we evaluate the genetic association of four SNPs in CDKN2B/CDKN2B-AS1 in a North Indian glaucomatous population.

Study participants
A hospital-based case-control association study was carried out; 461 unrelated cases with POAG and PACG were recruited along with 449 gender matched control individuals from Sardar Bahadur Sohan Singh Eye Hospital, Amritsar, Punjab after a complete ophthalmic examination. Approval for all research procedure was obtained from the Research Ethical Committee of Guru Nanak Dev University, Amritsar, India and the study protocols were according to the tenets of the Declaration of Helsinki. A written informed consent was obtained from all the study participants. A case was defined as having an open angle glaucoma if the patient had (1) an intraocular pressure (IOP) greater than 21 mm Hg in either of the eyes tested using Goldmann Applanation Tonometry and (2) glaucomatous optic nerve head damage defined as a vertical cup-disc ratio (CDR) of 0.7 or greater as adjudged clinically on slit lamp biomicroscopy using hand held + 90 D, this was confirmed using contrast enhanced fundus photograph on optical coherence tomography (OCT) as well as optic disc analysis or glaucomatous visual field defect as detected on Automated perimeter using Humphery's variants modifying susceptibility to glaucoma. Functional studies can further delineate the role of CDKN2B and CDKN2B-AS1 in primary glaucoma for therapeutic intervention.
Keywords: CDKN2B/CDKN2B-AS1, Primary glaucoma, North Indian population, Genetic association study, INK4 locus, SNPs Visual Field Analyser using Swedish Interactive Thresholding Algorithm (SITA) standard protocols. Individuals with known chronic systemic inflammatory, autoimmune or immunosuppressive disease as well as a pre-existing ocular disease (diabetic retinopathy, age-related macular degeneration) were excluded from the study. Individuals with history of corticosteroids, non-steroidal anti-inflammatory drugs, and topical use of steroids or prostaglandin analogues were also not recruited. Study controls, as examined by tonometry, slit lamp examination, CDR measurement and visual field assessment had normal IOP, optic disc and visual field. It was also ensured that the control samples did not have any family history of glaucoma.

Sample collection, DNA isolation and genotyping
Venous blood was collected in EDTA vacutainers for genotyping experiments and stored at -80 ºC till further use. Genomic DNA was extracted using standard phenol chloroform method [28]. Quantification of extracted DNA was done using NanoDrop ND-2000 spectrophotometer (NanoDrop Technologies, Wilmington, DE).
was performed using predesigned TaqMan real time PCR genotyping assay (Applied Biosystem, Foster city, CA; Catalogue no. C_2618046_10, C_341975_10, C_2618013_10 and C_11841829_10 respectively). Reactions were carried out in 48 well plates, in a 10 µL reaction volume using genomic DNA, TaqMan genotyping master mix (2X) (Applied Biosystems) and TaqMan Genotyping Assay (40X). The allelic discrimination (AD) assay was performed using StepOne Real Time PCR System (Applied Biosystems, Foster city, CA). Genotypes were scored by StepOne software v2.3 and manually confirmed by looking at the amplification plots. As a quality control measure, genotypes of few samples were retested by Sanger sequencing. The accession ID of sequences from which the primers were designed are NM_004936.3 (rs1063192), NM_004936.4 (rs3217992), NC_000009.12 (rs2157719 and rs4977756).

Statistical analysis
Descriptive statistical analysis for demographic and clinical characteristics of the study participants was performed using SPSS software. Student's t test was used to compare the baseline data among cases and controls. The results are represented as mean ± standard deviation (SD). Hardy Weinberg Equilibrium (HWE) and genetic association analysis was performed using PLINK software (v1.07). The distributions of genotype and allele frequencies among cases and controls were assessed by using 3 × 2 and 2 × 2 chi-square contingency tables and overall effect was determined in the form of Odd's ratio (OR) under different genetic models; dominant, recessive, co-dominant and additive at 95% confidence interval (CI). Haplotype software was used to determine haplotype frequencies. Pairwise linkage disequilibrium (LD) for the selected gene variants was calculated using the Expectation-Maximization (EM) algorithm as implemented in PLINK software. p values of less than 0.05 were considered to be statistically significant. Bonferroni correction was applied to reduce type I error (4 SNPs × 4 models); thus, p value was set at 0.05/16 = 0.003 for genetic model analysis.

Literature survey strategy
In the present study, a comprehensive literature survey was also performed, to find all the studies investigating the association among CDKN2B/CDKN2B-AS1 gene polymorphisms (rs3217992, rs1063192, rs2157719 and rs4977756), with primary open angle glaucoma. Online search using some popular websites including "Pubmed", "Medline", "Google scholar" and "Embase" was conducted upto 20-08-2019. The following keywords were used in the search: Glaucoma, POAG, CDKN2A/ CDKN2B, rs3217992, rs1063192, rs2157719 and rs4977756 and glaucoma, ocular hypertension and glaucoma. While searching literature, articles published only in English language were included in the study. Reference lists of original articles were screened to include maximum number of papers for the present meta-analysis.

Inclusion/exclusion criteria for meta-analysis
Only case-control studies conducted to examine the association between POAG and INK4 genetic variants were included. The inclusion of the case-control subjects in these studies was done according to the standard diagnosis criteria of glaucoma. Reported odds ratio (OR) or relative risk (RR) estimates with 95% confidence interval (CI), raw p values were incorporated in the present meta-analysis. Studies on secondary glaucoma or other types of primary glaucoma (PCG, PACG) and those for which odds ratio could not be calculated were excluded. Reviews, short/brief reports, abstracts which did not include details about genotypic frequencies, allelic frequencies, p values, odds ratios etc. were excluded from current analysis as shown in Fig. 1.

Data extraction and quality assessment
Two authors (Nanamika Thakur, Manu Kupani) extracted all the required and relevant information regarding the study design, first author, year of publication, place (country and ethnicity) where the study was conducted, definition of glaucoma, sample size, information about genotype number in both cases and controls, outcomes measured with 95% CI and p value adjusted variables. Assessment of the quality of each study was done by two authors using the Newcastle-Ottawa quality assessment scale (NOS) tool criteria [29]. NOS criteria comprises three parameters of quality assurance which include: selection (Case defined adequately, representativeness of cases, selection of controls and definition of controls), comparability (controls matched for age, gender and other confounders either at design stage or during analysis) and exposure (ascertainment of exposure, same selection method for ascertainment of cases/controls and same non response rate) [29]. The total score was either 7 or 8 which indicated that included studies had highquality scores. The individual score for all the studies has been given as Additional file 1: Table S1 (The Newcastle-Ottawa Scale for the assessment of case-control studies included in the meta-analysis). Any kind of conflicts regarding quality assurance was resolved by group discussions.

Statistical analysis
Present meta-analysis was done by software R (Version 3.6.1). Statistical analyses were made on estimates obtained from case control studies. Pooled adjusted odds ratio (OR) and its 95% CI was used as a measure to assess the association between CDKN2B-AS1 gene variants and glaucoma among cases and controls. Pooled OR was calculated for dominant model, heterozygous model, and recessive model. Since all the included studies were conducted on different ethnic groups with different sample sizes, therefore heterozygosity was tested using I 2 statistic. In heterozygosity testing, p value was < 0.05 and I 2 value was > 50% which indicated large difference among studied groups, hence data was combined using random effect models. The fixed effects model was used for studies with low heterogeneity. The I 2 value of > 75%, < 75%, < 50%, and < 25% represent considerable, substantial, moderate, and low heterogeneity, respectively. Sensitivity analysis was also performed by eliminating single study at a time to assess whether a single study could affect the overall results. Publication bias was investigated manually by observing funnel plots visually. Normal distribution (scatter plots) of studies in funnel plot indicated no publication biasing. These results were further verified by using Begg and Egger's regression tests that were applied only if the number of included studies in metaanalysis was ≥ 10. p value < 0.05 was considered as statistically significant for the overall effect.

Results
The association of CDKN2B (rs3217992, rs1063192) and CDKN2B-AS1 (rs2157719, rs4977756) gene polymorphisms with primary glaucoma cases was investigated in a total of 910 samples: 449 controls having cataract, 313 POAG cases and 148 PACG cases. All cases and controls were age and sex matched. Demographic and clinical parameters of the study participants were compared using t test, values represented as mean ± standard deviation (SD) in (Table 1). The success rate of genotyping for four variants was different. The total sample size for each of the variant is mentioned in the respective result tables. Prevalence of POAG was more among males (64.24%) as compared to PACG cases (35.75%), whereas the frequency of PACG was more in females (64.86%) as compared to POAG (35.13%) as given in Fig. 2. Overall there was no significant difference in the number of males and females in both cases and controls (p = 0.234). The genotype counts in controls for three SNPs (rs3217992, rs1063192 and rs4977756) followed HWE frequencies but rs2157719 showed significant deviation in controls (p value = 0.0003). To rule out genotyping errors, the amplification plots were rescored and genotyping was repeated on randomly selected samples. No genotyping discrepancy was observed. Since the genotypes for other three SNPs in the present study were in accordance with HWE frequency in controls, it is highly unlikely that the deviation is due to selection or nonrandom mating. However, to avoid any false positive association, Cochran Armitage (CA) trend test assuming additive model was applied for analyzing rs2157719 in addition to dominant, recessive and codominant models. The advantage of applying CA test is that it does not assume HWE [30].

Frequency distribution of rs3217992 C > T and rs1063192 A > G 3′ UTR polymorphisms of CDKN2B among POAG and PACG cases with respect to control subjects
For genetic association analysis, data was stratified into two groups, POAG cases and PACG cases and comparison was done with controls. Genotype and allele frequencies of the SNPs among cases and controls are given in Table 2. The genotype counts for both the variants were found to be consistent with Hardy Weinberg frequencies. The frequency of minor allele (T) for rs3217992 was significantly lower in cases (   (CI = 0.55-0.96)] respectively than control subjects. The frequency of risk allele i.e. C allele of rs3217992 C > T polymorphism was found to be significantly higher in POAG (60.23%) and PACG cases (62.59%) as compared to controls (55.07%) as mentioned in Table 2. Further, genetic model analysis revealed that TT + CT genotype conferred 0.73 [p = 0.047; OR = 0.73(CI = 0.54-0.99)] and 0.69-fold [p = 0.061; OR = 0.69(CI = 0.46-1.01)] protection against POAG and PACG development and in additive model, TT genotype gave 0.53-fold protection against PACG progression (p = 0.03). However, the association did not remain significant after Bonferroni correction. The distribution of allele frequencies between cases (POAG/PACG) for rs1063192 was not found to be statistically significant as given in Table 2. No significant difference in genotype distribution was obtained under genetic model analysis for any of the groups.

Mutant genotype 'TT' of rs3217992 confers marginal protection against PACG among male subjects
Stratification of POAG and PACG cases on the basis of gender revealed no significant difference in frequencies of CC and CT genotypes among male as well as female subjects (Table 3). However, TT genotype was overrepresented among control males with respect to PACG cases (p = 0.058, OR = 0.36; CI = 0.12-1.03). When the comparison in context to mutant genotype (TT) was done among POAG/PACG females with control female subjects, the study failed to obtain any significant difference in the genotype frequency. Genetic model analysis did not reveal a statistically significant difference in genotypic distribution among primary glaucoma (POAG/ PACG) males and females as compared to control subjects (Table 3). Table 4 shows the comparison of allele and genotype frequencies of CDKN2B-AS1 rs2157719 and rs4977756 intronic polymorphisms. Categorization of two groups was done in a similar fashion as in case of rs3217992 and rs1063192 variants. The genotype count for rs4977756 in controls followed HWE frequency but showed marked deviation for rs2157719. We presume that this deviation might be due to chance factors as other variants were in HWE in controls [31]. The minor allele (C) frequency of rs2157719 was found be significantly higher among controls as compared to cases (POAG/PACG). For  rs4977756 minor allele G was almost equally prevalent among POAG cases and controls but it was found to be higher among controls than PACG cases. For rs2157719, significant difference was observed for both allele as well as genotype frequency distribution on comparing POAG and PACG cases with controls as given in Table 4. Since the distribution of genotypes in the population deviates from HWE, the frequency of genotypes was compared by the Cochran-Armitage test for trend assuming additive model [32]. However since our combined sample is in HWE, the allelic and trend statistic were found to be almost equivalent for the combined dataset [32]. After segregating the samples, the trend test revealed association in POAG and PACG (p values = 0.028 and 0.004 respectively) but the results were non-significant after Bonferroni correction. Genetic data analysis of rs4977756 revealed G allele as the risk allele among POAG cases while A allele was found to confer risk towards PACG, though the distribution of allele frequencies between cases (POAG/PACG) was not found to be statistically significant as given in Table 4.

CC genotype of rs2157719 confers protection against primary glaucoma (POAG/PACG) among male and female subjects
Subsequent segregation of POAG and PACG cases based on gender (Table 5) revealed a different frequency distribution of genotypes among males and females of two groups i.e. POAG and PACG. A lesser number of CC mutant homozygotes were observed among patients (POAG/PACG) of both genders (Table 5). CC genotype conferred 0.33-and 0.20-fold protection among POAG (p = 0.0041) and PACG (p = 0.035) males with respect to control males. In the genetic model analysis, recessive  model provided 0.29-and 0.22-fold protection against glaucoma (POAG/PACG) progression among male subjects. When the comparison of POAG females was done with healthy female subjects, the study failed to obtain any significant difference in the genotypic frequencies.
Genetic model investigation also did not reveal any statistically significant difference in the genotypic distribution among POAG females and control female subjects. In contrast to POAG females, the frequency of CC genotype was significantly different among females having PACG (p = 0.029) with respect to control females. The recessive model unveiled 0.34-fold protection against PACG development in females (p = 0.033, OR = 0.34; 0.12-0.91) but the results were not significant after Bonferroni correction as shown in Table 5.

Meta-analysis Study Characteristics (rs1063192 A > G)
Total 125 articles were retrieved by literature survey.
After examining all publications, 18 related articles were included along with the current study. Therefore, the current meta-analysis was performed with 19 total studies. This meta-analysis comprised all articles representing total cases (n = 13,253) and controls (n = 34,615) from 2011 to 2019. All included studies dealt with association between 4 sequence variants of CDKN2B/AS1 genes and POAG.

Study characteristics (rs3217992 C > T)
Meta-analysis of rs3217992 comprised of 4 studies including total 2088 affected and 10,266 unaffected subjects. It included 2 studies from Indian population (including the current study). 1 report represented Caucasians and 1 from Japan. All the essential information of included studies is given in Table 7. Diverse genotyping techniques used in selected studies were TaqMan, Seq-illumina and PCR-RFLP method. The genotypes followed HWE in controls in all studies.

Meta-analysis of rs1063192 and rs3217992
Main features of results of meta-analysis for rs1063192 are shown in Table 8. Under random effect model, overall a statistically significant association was found between this sequence variant and POAG risk in current metaanalysis. Dominant model (GG + GA vs AA) revealed significant association with OR = 0.76; p < 0.0001, with substantial heterogeneity, I 2 = 64.6%. In recessive model (GG vs GA + AA), highly statistically significant association was observed with overall protective effect of GG genotype from disease progression, OR = 0.64; p < 0.00001, with no heterogeneity, I 2 = 0% as shown in Fig. 4. When we assessed this SNP under heterozygous model (GA vs GG + AA), moderate heterogeneity (I 2 = 43.6%) with OR = 0.84 (95% CI = 0.80-0.88) was observed ( Table 8). The sensitivity analysis for recessive model was done, which investigated the influence of a single study on the overall risk estimate by omitting one study at each turn. This yielded a range of OR from 0.6362 (95% CI 0.5818-0.6956); p < 0.0001 to 0.6625 (95% CI 0.5955-0.7369); p < 0.0001. This suggested that exclusion of any single study did not altered the overall combined OR.

Publication bias for rs3217992
Visual inspection of the funnel plot did not identify obvious asymmetry. The Begg (p = 0.383) and Egger's test (p = 0.348) for funnel plot did not reveal any asymmetry in funnel plot (Fig. 5). For rs3217992, the pooled OR ranged between 1.03 and 1.23, under various genetic model (Dominant, Recessive (Fig. 6), Heterozygous) analysis as given in Table 9 with 95% CI ranging from 0.88 to 1.59 with overall non-significant p value. The sensitivity analysis for recessive model was done, which investigated the overall risk estimate by omitting the present study. This yielded a range of OR: 1.3698 (95% CI 1.1611-1.6160); p < 0.0001 to 0.6625 (95% CI 0.5955-0.7369); p = 0.0002 with 19.2% heterogeneity which differed from pooled OR: 1.2393 (95% CI 0.9612-1.5977), p = 0.0979 with 67.8% heterogeneity (Fig. 6). This suggested that exclusion of present single study obviously altered the overall combined significant association and gave maximum heterogeneity.

Publication bias for rs1063192
Visual inspection of the funnel plot did not identify obvious asymmetry. The Egger test for funnel plot asymmetry was not performed as the power of this test is too low to distinguish chance from real asymmetry when the metaanalysis included less than 10 studies.

Study characteristics (rs2157719A > G)
Samples from 5 studies (cases = 2819 and controls = 3283) were pooled together to calculate the association of rs2157719 with POAG. It included 1 Caucasian and 4 Asian studies. Out of 4 Asian studies, 2 studies were from Indian population including the current study, 1 report was from China and 1 from Japan as given in Table 7. TaqMan, Sequenom, Japonica, and Illumina assays were used to genotype samples and genotypic distribution in the studied controls were consistent with HWE.

Study characteristics (rs4977756 A > G)
Total 12 studies were short listed for the meta-analysis of rs4977756 to find out overall genotypic association between this variant and POAG. Total 8739 POAG patients and 21,174 unaffected healthy age and gender matched controls were involved. Meta-analysis contained 3 studies on Caucasians, 4 reports on African ancestry, 5 studies from Asians (out of which 1 was conducted in Japan, 1 in China, 1 report from Pakistan and 2 from India, including current study) as represented in Table 7. Genotyping techniques included TaqMan chemistry, Seq-Illumina, Sequenom assay, and PCR-RFLP methods. The genotypic distribution in all studied controls were consistent with HWE.

Meta-analysis of rs2157719 and rs4977756
Main features of results of meta-analysis for rs2157719 are given in Table 10. Under random effect model, overall a statistically significant association was found between this sequence variant and POAG risk in current metaanalysis under various genetic models (dominant, recessive (Fig. 7) and heterozygous model) (p = 0.01, 0.007, < 0.0001 respectively) with pooled OR ranging from 0.46 to 0.77), showing overall protection against POAG as given in Table 10 and Fig. 7. The sensitivity analysis for recessive model was done, which investigated the huge impact of a single study conducted by Vishal et al. in 2014 [26] on the overall risk estimate (OR = 0.4048 (95% CI 0.2758-0.5943); p < 0.0001) and hence conferring maximum heterogeneity to the data. After omitting this single study there was no heterogeneity (I 2 = 0%), suggesting that exclusion of this single study altered the overall combined OR.

Publication bias
Visual inspection of the funnel plot did not identify obvious asymmetry. The Egger test for funnel plot asymmetry was not performed as the power of this test was too low to distinguish chance from real asymmetry when the meta-analysis included less than 10 studies. For rs4977756, the pooled OR was 0.85, under various genetic model (Dominant, Recessive) analysis as given Fig. 4 Forest plot depicting association between CDKN2B variant (rs1063192) under recessive model (GG vs AA + GA) and POAG. The black horizontal lines correspond to the 95% CI and the position of black small squares to the OR per study. The size of the square is proportional to the relative weight of that study in % to compute the overall OR (black diamond). The width of the diamond represents the 95% CI of the overall OR. This study has revealed statistically significant association between this SNP and POAG. This meta-analysis illustrates protection conferred by GG genotype against glaucoma progression with overall OR = 0.64 (p ≤ 0.0001) with no heterogeneity (I 2 = 0%; p = 0.82) under fixed effect model  Table 11 with 95% CI ranging from 0.75 to 0.97 with overall highly statistically significant p value (Table 11). Comparison of heterozygous GA with other two genotypes (GG + AA) among total cases and controls also unveiled protection against POAG with OR = 0.90 (95% CI 0.85-0.95) and highly significant p = 0.0005 with very low heterogeneity (I 2 = 27.9%) among all included studies in meta-analysis as shown in Fig. 8. The sensitivity analysis for the overall risk estimate by omitting one study at a time, yielded OR: 0.9033 (95% CI 0.8532-0.9564) with I 2 = 27.9%. This suggested that exclusion of single study obviously did not alter the overall combined significant association.
Publication bias Visual inspection of the funnel plot did not identify obvious asymmetry. The Begg (p = 0.217) and Egger's test (p = 0.095) for funnel plot did not reveal any asymmetry in funnel plot as given in Fig. 9.  Fig. 4. This meta-analysis did not unveil any association between this SNP and POAG with overall OR = 1.24 (p = 0.097) with substantial heterogeneity (I 2 = 67.8%; p = 0.03) under random effect model Table 9 Meta-analysis of association between CDKN2B rs3217992 polymorphism and the risk of POAG   Fig. 4. This meta-analysis unveiled highly significant association between this polymorphism and POAG with overall OR = 0.53 (p = 0.007) with substantial heterogeneity (I 2 = 52%; p = 0.08) under random effect model

Discussion
INK4 locus at chromosome 9p21 has been implicated in affecting genetic susceptibility to glaucoma in different populations. The present study is first of its kind in North India to investigate the role of variants in CDKN2B (rs1063192, rs3217992) and CDKN2B-AS1 (rs2157719, rs4977756) genes located at the INK4 locus as genetic risk factors for POAG and PACG. CDKN2B is cyclin dependent kinase inhibitor gene that influences cell proliferation and senescence in all cell types [33,34] including retinal ganglion cells (RGCs). It is reported to be upregulated in response to elevated IOP in animal model of glaucoma suggesting apoptosis of RGCs may be promoted due to altered expression of CDKN2B and cause glaucomatous visual field loss [8]. Various studies reveal that single nucleotide polymorphisms in CDKN2B can alter its function [16]. In one such study conducted by Pasquale et al. [17] 10 SNPs at CDKN2B were screened in a retrospective observational case series and rs3217992 in CDKN2B was found to be significantly associated with glaucoma (p = 3.34 × 10 −15 ). The study included Caucasian POAG patients that presented increased disease risk with larger CDR (p = 4.74 × 10 −4 ) despite lower IOP in presence of the minor allele (T). Earlier GWA studies have also reported significant association of rs3217992 with POAG in the Japanese and Caucasians [12,35,36]. In initial analysis our findings are in line with the previous findings, however in our population T allele confers protection as POAG cases were found to have significantly lower frequency of the minor allele (T) of rs3217992 variant (39.77%) as compared to controls (44.93%) (p = 0.045). Same relation was observed in PACG, with T allele being significantly lower among cases as compared to controls (p = 0.024). This was in contradiction to the previous finding [17,35] and points towards different genetic structure among various ethnic groups. Further stratification of the data unveiled the ability of TT genotype to confer 0.36-fold protection towards progression of PACG in males, although we could not achieve same significant association after applying correction for multiple testing. No significant difference was observed in the female population, although the overall risk of developing glaucoma is reported to be higher in older women experiencing menopause compared to males which has been attributed to estrogen specific protection in younger females [37]. Another variant, rs1063192 (A > G) resides in the 3′UTR of CDKN2B gene and is strongly correlated with (CDR in a Japanese group [38]. A previous Japanese GWAS reported a strong association of the variant with POAG [(p = 5.2 × 10 −11 ); (OR = 0.75)] [14]. The association reports of this SNP with POAG has been quite consistent in various Caucasian studies [9,12,36,39,40]. Analysis in Afro-Caribbean population including a total of 437 unrelated subjects from Barbados family study of open angle glaucoma also found rs1063192 to be significantly associated with POAG (p = 0.0008) [16]. Yet in the present study, rs1063192 failed to show any significant Fig. 9 Funnel plot of included case-control studies under fixed effect model assessing no publication bias and hence symmetry association with either disease form. A study from African group also suggested no association of the variant with glaucoma contrary to the reports from the Caucasian groups. In that study, the researchers included two groups, African-Americans and Ghanaians and analyzed 57 SNPs in five loci including CDKN2B/CDKN2B-AS1, TMCO1, CAV1/CAV2, chromosome 8q22 intergenic region, and SIX1/SIX6. rs10120688 in CDKN2B-AS1 region was found to be associated with POAG in African Americans (p = 0.0020) [41]. According to the Hap-Map data, this variant is in strong LD with rs2157719 (D′ = 0.99, r 2 = 0.70), rs1063192 (D′ = 0.97, r 2 = 0.68) and rs4977756 (D′ = 0.49, r 2 = 0.67) in Gujarati Indian population. Hence, it is difficult to predict the true functional variant. Nonetheless, in the South Indian population, rs1063192 was associated with one of the endophenotypes of POAG. The researchers studied 97 POAG cases and 371 controls from South India and found C allele of rs1063192 to be associated with decreased axial length in controls suggesting a decreased risk for POAG [27]. Another gene of interest located at INK4 locus is CDKN2B-AS1 which encodes for an antisense RNA that regulates expression of CDKN2A and CDKN2B via forming transcript complexes with polycomb proteins. The antisense region controls the expression of inhibitor genes negatively. Depletion of ANRIL was observed to cause increased expression of CDKN2B gene confirming that ANRIL binds SUZ12 (one of polycomb protein complexes 2) and regulates CDKN2B [40]. The corresponding change in cyclin dependent kinase could affect ganglion cell apoptosis. Due to its suggestive role in ganglion cell apoptosis, the gene has been of special interest in GWASs of various ethnicities including US-Caucasians, Asians, Africans and Europeans [13-15, 19, 35]. rs2157719 (T > C), located in CDKN2B-AS1 has been reported to be associated with exfoliation glaucoma in Caucasian population [15]. A Japanese group recently did genotype-phenotype analysis and revealed significant correlations between the variant and decreased IOP (β = − 6.89 mmHg, dominant model p = 1.70E−04) [42]. Ozel and co-workers also reported a significant association of rs4977756 and rs2157719 at CDKN2B-AS1 with IOP in a meta-analysis in > 6000 subjects of European ancestry collected in three datasets: NEIGH-BOR, GLAUGEN and AMD-MMAP [42,43]. Case-control replication analyses by Yoshikawa and co-workers also yielded strong association of the rs4977756 with POAG in a Japanese population [(p = 3.2 × 10 −4 ); OR = 0.77(CI = 0.68-0.86)] [22]. The present study also found significant association of rs2157719 with primary glaucoma on using CA trend test for association. The advantage of the CA trend test is that it is not dependent on the HWE assumption. Moreover, since three SNPs followed HWE, and genotyping errors were ruled out to the best of our ability, we can assume that the controls were representative of the population targeted in the study. Sex-based data stratification further revealed protection towards progression of both disease forms conferred by CC genotype among males. Similar results were obtained in PACG females where CC genotype provided protection against developing the disease but not in POAG females. Even though the observed association did not survive Bonferroni correction, it should not be rejected out rightly. Since the bonferroni correction is overly conservative as it assumes independence among the tests considered, we would be precautious to reject the positive association with rs321792 and rs2157719.
More so because the SNPs tested in the current study are in strong LD and this may lead to correlation among the tests [44]. rs4977756 (A > G) failed to show association in our study. The results for this variant have been inconsistent; a Brazilian study denied significant association with glaucoma [45]. Cao et al. also didn't find any significant association of this SNP with POAG (p = 0.4507) in Afro-Caribbean population [16].
Nominal association of rs4977756 with NTG was observed in African Americans (p = 0.04) but not in West African population [41]. Ng et al. [40] also reported strong association of the variant with POAG in Australians, especially in the advanced cases. They further established a stronger relation of SNP with POAG in females with OR difference between male and female being statistically significant at p = 1.73 × 10 −4 [40]. However, no such predisposition in the females was seen in our study group for this variant.
Previous findings in India and Pakistan have been consistent in terms that no association was observed between glaucoma and INK4 locus [26,46]. A study in the East Indians included all 4 SNPs targeted in this study and found only nominal association of rs1011970 (p = 0.048) with POAG and rs10120688 (p = 0.048) in NTG patients (IOP < 21 mm of Hg) thereby suggesting lack of significant genetic association of 9p21 locus with POAG [26]. This is in contrast to our findings and further reiterates the heterogeneous nature of populations in the Indian subcontinent due to genetic admixture.
Since there is a significant difference in the results of various studied populations worldwide, we did a meta-analysis to measure the overall association between these variants and risk of primary open angle progression. Pooled analysis on POAG cases and controls revealed a significant association between rs1063192, rs2157719, rs4977756 and POAG except rs3217992. Another meta-analysis, conducted by Hu et al. [47] in 2017 on 11,316 cases and 24,055 controls evaluated CDKN2B rs1063192 and it was found to be associated with decreased risk of glaucoma. Our findings were also in the concordance with previously done meta-analysis [47], rs1063192 was observed to be associated with POAG under fixed and random effect model with pooled OR ranging from 0.64 to 0.84 under dominant, heterozygous, recessive genetic model analysis with overall p ≤ 0.0001. Meta-analysis of rs2147719 and rs4977756 also showed statistically significant association of these SNPs with decreased risk of POAG progression with pooled OR ranging from 0.46 to 0.77 and 0.85 to 0.90 respectively, under various genetic models. rs3217992 did not give any association in meta-analysis although found to be associated with primary glaucoma in genetic association analysis in the current study. Whether the association obtained with rs3217992 in the present study represents real cause-effect model relationship as its not supported by meta-analysis should be further investigated. Significant association in some studies could indicate these sequence variants to be a risk factor in certain ethnic groups, which requires more data from different populations to do a meta-analysis and correlate it with ethnicity.

Conclusion
In conclusion, the present study points towards significant association between rs3217992 and rs2157719 with primary glaucoma. The association with rs2157719 with PACG was strong enough to survive Bonferroni correction for multiple testing. In an updated meta-analysis for INK4 variants significant association was observed with POAG with rs1063192, rs2157719 and rs4977756. Large scale sequencing of the INK4 locus may reasonably detect true functional and rare variants while in vitro and in vivo studies may further assess the functional relevance of these variants in pathogenesis of primary glaucoma.