Genome-wide association study identifies seven novel susceptibility loci for primary open-angle glaucoma Free

Abstract

Introduction

Glaucoma, one of the leading causes of blindness worldwide, is a neurodegenerative optic neuropathy that leads to irreversible visual field loss (1). Primary open-angle glaucoma (POAG), the most common form of glaucoma, has a complex etiology that includes genetic predisposition. Recent progress using genome-wide association studies (GWASs) has uncovered 15 POAG-susceptibility loci (2–10). However, the majority of these loci were discovered in European (EUR) ancestry populations (2–7). Although previous GWASs in Asian populations have identified five POAG-associated loci, most of these loci were also POAG-associated in EUR populations (8–10). Larger GWASs in Asian populations will be needed to identify new POAG-associated loci that provide new insights into the genetic basis of POAG.

Here, we performed a large-scale GWAS and replication studies for POAG using a total of 43 763 subjects from a Japanese population. Additionally, we examined the association of novel susceptibility variants in a multiethnic population using 46 457 subjects from Singaporean Chinese, EUR and African populations. We also carried out functional annotation of the associated variants, an expression analysis of the positional candidate genes at novel POAG loci in macaque ocular tissues, and a gene-set enrichment analysis. Finally, we investigated genetic correlations between POAG and systemic diseases.

Results

The overall study design is outlined in Supplementary Material, Figure S1. First, we performed a GWAS using 3980 POAG cases and 18 815 controls from a Japanese population. The association analysis was performed using 6 109 376 single nucleotide polymorphisms (SNPs) following whole-genome imputation of 515 326 genotyped SNPs using East Asian (EAS) samples from the 1000 Genomes Project (phase1 version 3) for reference. We observed a slightly high genomic inflation factor (λ_GC = 1.14). However, when we consider the sample size of the GWAS, the observed inflation was within an acceptable range (λ₁₀₀₀ = 1.02). Therefore, we used uncorrected P-values for further analyses (Supplementary Material, Fig. S2). The GWAS identified six loci that were significantly associated with POAG (P_GWAS < 5.0 × 10⁻⁸), including three known (CDKN2B-AS1, ABCA1 and SIX6) and three novel loci for POAG (ANKRD55-MAP3K1, TAP2 and LOXL1). Conditional analysis of the top SNPs demonstrated that none of the variants had significantly associated secondary signals at these loci. When we checked the associations of 15 previously identified SNPs, 9 SNPs had nominal associations (P_GWAS < 0.05) and 6 of them reached the threshold for replication after Bonferroni correction (P < 0.0033 = 0.05/15; Supplementary Material, Table S1).

To confirm the association of GWAS and to search additional susceptibility loci, we selected 33 top SNPs from 33 loci including the 6 loci with P_GWAS < 5.0 x10⁻⁸ that showed association at P_GWAS < 1.0 x10⁻⁵ (Supplementary Material, Table S2). For the replication study, two independent sets comprising a total of 3398 cases and 17 570 controls from Japan were used (see Materials and Methods for detail). We genotyped 33 SNPs from the 33 loci mentioned earlier in these sets and combined the results from these two data sets using the inverse-variance method. Seven SNPs met criteria for replication after Bonferroni correction (P < 0.0015 = 0.05/33; Supplementary Material, Table S2). Then, we further performed a combined analysis by pooling the results from the GWAS and the two replication sets. We identified a total of eleven POAG-associated loci, including four known and seven novel loci at P_combined < 5.0 × 10⁻⁸, increasing the total number of POAG-associated loci from 15 to 22 (Fig. 1 and Table 1). Regional association plots of the newly identified loci are shown in Supplementary Material, Figure S3. To further assess the clinical relevance of 11 POAG-associated loci, we performed the stratified analyses of GWAS data by dividing POAG cases into normal-tension glaucoma and high-tension glaucoma. Although most of the associations were weakened probably due to the smaller sample sizes, all SNPs showed the same direction of effect in both subtypes (Supplementary Material, Table S3), implying that these POAG-associated loci may affect both glaucoma subtypes.

Figure 1.

 Manhattan plot depicting the results of the genome-wide association study (GWAS) and replication studies in the Japanese population. The x-axis shows the chromosomes and their positions. The y-axis shows −log10 P-values of single nucleotide polymorphisms (SNPs). The orange horizontal line is the genome-wide significance level (P=5.0 × 10⁻⁸). Plots of the seven novel associated loci discovered by the Japanese GWAS are highlighted in red. Combined-analyses of the GWAS and Japanese replication sets are shown as diamonds for the loci which have not been reported in East Asians. Previously identified loci that showed a significant genome-wide significant association in the Japanese GWAS are highlighted in blue.

Open in new tabDownload slide

To check the relevance of the novel POAG-associated loci in other ethnic populations, we looked up the association of 7 SNPs in GWAS data from the Singaporean Chinese (9), EUR (5,7) and Africans (11,12) (see Supplementary Material Note). We found the nominal association for 2 SNPs in Singaporean Chinese and 4 SNPs in EUR, however, no variant showed significant association after Bonferroni correction (P < 0.0071 = 0.05/7; Table 2 and Supplementary Material, Table S4).

To investigate the biological roles linked to these loci, we annotated the variants in linkage disequilibrium (LD, r² ≥ 0.8) with the lead SNPs at the seven newly identified loci using HaploReg (ver. 4.1) (13) (Supplementary Material, Table S5). We found that six lead variants were in LD with variants located at conserved regions predicted by SiPhy-omega (14) or GERP (15). Although the lead variant at LOXL1 was a missense variant (rs1048661, R141L) that increases the risk of exfoliation syndrome (16), other variants were not in LD with variants at coding regions. Recent study of exfoliation syndrome identified the association of another missense variant (rs3825942, G153D) at LOXL1 (17). However, the association of rs3825942 was weak (P_GWAS = 0.002), and the association of rs1048661 remained after conditioning on rs3825942 (P_GWAS = 1.13×10⁻⁸ and P_condition = 1.57 × 10⁻⁶). This indicates that these associations may be independent. Since we excluded the exfoliation syndrome from POAG cases in advance, we consider that rs1048661 in LOXL1 is associated with Japanese POAG. Annotations using the 15-state model of the Roadmap Epigenomics project (18) indicated that all lead variants were in LD with promoter or enhancer regions in at least one of the tissues investigated. Next, we evaluated the overlaps between lead SNPs and cis-expression quantitative trait loci (eQTLs) (19,20), and found a link between rs10819187 and eQTLs of LMX1B in adipocytes. However, we could not identify any other loci that overlapped with eQTLs. We defined 8 positional candidate genes, including 2 flanking genes for the ANKRD55-MAP3K1 intergenic loci, from the newly identified loci based on their positions from each lead SNP (see Materials and Methods), and evaluated their transcript levels in macaque ocular tissues. We found that all candidate genes were expressed at key tissues (retina or the trabecular meshwork) linked to POAG-pathogenesis (Fig. 2).

Figure 2.

 Expression of positional candidate genes at novel primary open-angle glaucoma (POAG)-associated loci in ocular tissues. Expression of the positional candidate genes (ANKRD55, MAP3K1, LMX1B, LHPP, HMGA2, MEIS2) at five novel POAG-associated loci in different parts of the monkey eye. Two loci (FNDC3B, LOXL1) were used as positive controls because the expression of these genes was confirmed in the retina and the trabecular meshwork. Lysates from the COS-7 cell line were used as controls. The internal control gene was GAPDH.

Open in new tabDownload slide

Next, we investigated pleiotropy. Of the seven newly identified loci, two were reportedly associated with other traits (r² between lead SNPs ≥ 0.8 in the 1000 Genomes Project; Supplementary Material, Table S5): HMGA2 in type 2 diabetes (21) and LOXL1 in exfoliation syndrome (22). The lead SNP at FNDC3B (rs7636836) was relatively close (227 kilobases) to that of a known intraocular pressure locus (23) (rs6445055); however, these SNPs were not in LD (r² < 0.03 in both EAS and EUR populations). Because rs6445055 was putatively associated with POAG in our GWAS (P = 5.98 × 10⁻⁵, OR = 1.12), we performed a conditional analysis that adjusted for rs7636836 and found that POAG-association with rs6445055 was likely to be independent of that of rs7636836 (P_condition = 5.36 × 10⁻⁴, OR = 1.10).

To identify biological pathways involved in POAG pathogenesis, we performed a pathway analysis using MAGENTA (24). This analysis revealed significant enrichment of the epidermal growth factor receptor (EGFR) signaling pathway (false discovery rate [FDR] q < 5%; Supplementary Material, Table S6). In addition to MAP3K1, a positional candidate gene identified in this study, other EGFR-pathway genes were nominally associated with POAG (Supplementary Material, Table S7).

POAG is an ocular disease and although previous epidemiological studies have suggested links between POAG and systemic diseases, these relationships are not well defined (25). To evaluate the genetic links between POAG and systemic traits, we conducted a cross-trait LD score regression analysis (26) using the Biobank Japan Project GWAS data for four diseases (chronic kidney disease, ischemic stroke, myocardial infarction and type 2 diabetes) and three clinical measurements (body mass index, and systolic and diastolic blood pressure) (Supplementary Material, Table S8). We found a significant positive genetic correlation with type 2 diabetes (q < 0.05, r_g = 0.27 ± 0.07; Table 3) and nominal genetic correlations with myocardial infarction and ischemic stroke (P < 0.05).

Discussion

In this study, we set out to search for novel POAG-associated loci by taking the advantage of analyzing Japanese population, a distinct East Asian population with poorly characterized genetic background with regard to POAG. To accomplish this, we combined the results of the GWAS and a replication study in Japanese and successfully identified 11 POAG-associated loci at the genome-wide significance level, including 7 novel loci (FNDC3B, ANKRD55-MAP3K1, LMX1B, LHPP, HMGA2, MEIS2 and LOXL1). The functional annotation of the seven newly identified loci revealed that all the lead SNPs were in LD with functional SNPs in promoter or enhancer regions in at least one of the tissues investigated. This suggests that the lead SNPs linked to POAG-associated loci could potentially influence the expression of nearby genes. Next, we found a link between rs10819187 and eQTLs of LMX1B in adipocytes. However, we could not identify any other loci that overlapped with eQTLs. One possible explanation for this is that the existing databases do not include ocular tissues. For pleiotropy analysis, we performed a conditional analysis that adjusted for rs7636836 and found that POAG-association with rs6445055 was likely to be independent of that of rs7636836, suggesting that distinct causal variants may exist at this locus. Future fine-mapping or sequencing studies will clarify this.

Two of the positional candidate genes identified in this study, MAP3K1 and LMX1B, have previously been linked with ocular development (27,28). Mutations in LMX1B, which cause nail patella syndrome, are known to cause glaucoma due to abnormal development of the anterior segment of the eye (29). In the LMX1B locus, the top variant showing strongest eQTL effect for LMX1B among analyzed variants in our GWAS was rs10987373, which showed third top association for POAG (P_GWAS = 1.31×10⁻⁶, r² with rs10819187 = 0.85 in EAS samples of 1KG). Risk allele of this SNP (G allele) was associated with the increased expression of LMX1B in GTEx, implying that increased level of this gene product might confer the risk of POAG. Nevertheless, considering that the individuals with loss of function variants in LMX1B (Nail–Patella syndrome) was reported to be at higher risk for POAG, these observations appeared to be inconsistent. Further investigations will be needed to clarify the role of LMX1B in the pathogenesis of POAG. MEIS2 interacts with FOXC1, a gene located at a known POAG-susceptibility locus (7), and is expressed in the trabecular meshwork (30). HMGA2 is also important for the development of the eye (31). These results suggest that genes related to ocular development may play a role in POAG pathogenesis. Further analyses that include assessment of eQTL in these ocular tissues will be required to confirm which genes are associated with POAG.

Although we excluded exfoliation syndrome from POAG cases, we found a significant association of a missense variant in LOXL1 (rs1048661, R141L) with POAG susceptibility. The association between rs1048661 and glaucoma is known to be explained by exfoliation glaucoma in European population (22), however previous meta-analysis of LOXL1 gene polymorphisms (32) suggested the association of rs1048661 with POAG in Asian population (P = 0.03). LOXL1 gene product has been reported to modulate extracellular matrix biogenesis by cross-linking elastin and collagen in connective tissues (33). Elastin is a major component of the elastic fibers of the extracellular matrix of the lamina cribrosa, and deformity of the lamina cribrosa has been postulated to damage retinal ganglion cell axons that penetrate its pores (34). Since rs1048661 is reported to affect the LOXL1 expression (22), LOXL1 dysregulation may have a possibility to alter the structural integrity of retinal ganglion cells and their axons, and resulted in the pathogenesis of POAG.

The results of the pathway analysis suggested the involvement of the EGFR signaling in POAG pathogenesis. The EGFR signaling affects the development and fate of the retinal ganglion cells, e.g. axon growth and death of these cells (35,36). A previous study showed that the EGFR signaling in concert with complement C3 signaling, both of which are active in astrocytes at the optic nerve head (37,38), may act against glaucoma (39). Conversely, intraocular pressure-dependent activation of EGFR may induce NOS-2 expression, which may cause oxidative stress that are deemed damaging for the retinal ganglion cells (38), particularly in glaucoma (40). Furthermore, EGFR appears to regulate various pathways in astrocytes including MAPK signaling and PI3K-AKT signaling that may be important for the biology of retinal ganglion cells in glaucoma (39). For example, increased MAPK signaling activation has been reported in neurons and glia in postmortem eyes from glaucoma patients (41). Furthermore, we have reported that an activation of PI3K-AKT signaling protects retinal ganglion cells in an animal model of glaucoma (42). Taken together, despite some opposing effects the EGFR that remains to be defined in the context of glaucoma neurobiology, it is not surprising that EGFR signaling were computationally predicted as the important pathway potentially contributing to POAG pathogenesis.

We performed an extended analysis by checking the associations of the newly identified POAG-loci in multiethnic populations. Interestingly, we did not find a significant association of the lead SNPs with POAG in any of the populations including the Singaporean Chinese. This could be partly explained by combination of several reasons including the differences in the allele frequency, limited sample size or the difference in LD structures. Therefore, at this point, it is difficult to conclude that the novel POAG-associated loci/genes are unique to Japanese population. Indeed, at least for the Singaporean Chinese, although none reached statistical significance, two (FNDC3B and LMX1B) loci were nominally associated and all SNPs showed the same direction of the effect.

Finally, we evaluate the genetic links between POAG and systemic traits and found a significant positive genetic correlation with type 2 diabetes, and nominal genetic correlations with myocardial infarction and ischemic stroke. These relationships were consistent with the epidemiological findings (43–45). A genetic link between POAG and type 2 diabetes was also suggested by a Mendelian randomization analysis (46). Although our GWAS cases included 956 individuals (24.0%) with type 2 diabetes, 182 (4.6%) with myocardial infarction, 415 (10.4%) with ischemic stroke and 334 (8.4%) with angina as the underlying diseases, overlapping of analyzed samples should not influence the results of genetic correlation analysis because we applied bivariate LD score regression which provides robust and unbiased estimates even if samples were overlapped in analyzed GWASs (26). These results imply that some of the genetic components conferring an increased risk of POAG also increase the risk of type 2 diabetes and cardiovascular diseases, providing some genetic support for epidemiological relationships between POAG and systemic diseases.

In summary, we performed a large-scale association analysis in Japanese and identified seven novel POAG-susceptibility loci, none of which were relevant in other populations in our study. Although this may limit the applicability of our findings to other ethnicities, our results imply that genes related to ocular development and the EGFR signaling pathway may play a role in POAG pathogenesis, and that some of the POAG risk variants are also risk factors in the development of some systemic diseases. These findings provide novel insight into the genetic factors underlying POAG.

Materials and Methods

Subjects

The characteristics of study subjects are shown in Supplementary Material, Table S9 and the Supplementary Material Note. For the GWAS cases, we used 3980 POAG patients recruited by the BioBank Japan (BBJ) project (47,48). We used genome-wide screening data from 18 815 subjects derived from two Japanese population-based cohorts as controls for the GWAS. These cohorts were from the Japan Multi-Institutional Collaborative Cohort Study (J-MICC) (49) and the Japan Public Health Center-based Prospective Study (JPHC).

For the Japanese replication study, we used two independent case–control sets matched according to their locations within Japan. The cases in the first replication set included 1101 POAG patients recruited at Tohoku University and its related hospitals. Control data for the first replication set comprised genome-wide screening data from 7716 subjects collected by the Tohoku Medical Megabank Project (50). The cases in the second replication set included 2297 POAG patients recruited from all areas of Japan (Kyoto University, Fukuiken Saiseikai Hospital, Kobe University, Hiroshima University, Kyushu University and Oka Eye Hospital for the Japan Glaucoma Society Omics Group (JGS-OG)). Control data for the second replication set comprised genome-wide screening data from 9854 subjects collected by the BBJ project after the patients’medical records were used to exclude those with glaucoma.

The diagnoses of POAG cases in the BBJ project and in the replication studies were determined by ophthalmologists in the collaborating hospitals. Detailed criteria are described in our previous report (10). The diagnoses of POAG were determined by expert ophthalmologists from the participating hospitals based on the criterion: glaucomatous optic neuropathy with corresponding glaucomatous visual-field defects measured using different opthalmic perimeters. Combinations of the following examinations were also used to make diagnoses: slit-lamp microscopy, examination of the optic disc and retinal nerve-fiber layer using an opthalmoscope, measurement of the vertical cup-to-disc ratio, tonometry to determine intraocular pressure and gonioscopy to assess the iridocorneal angle (an open angle is required for a diagnosis of POAG). Patients with secondary glaucoma including those caused by exfoliation syndrome, closed-angle glaucoma or congenital glaucoma were excluded from this study. For the subtype analysis of the associated loci, we divided GWAS cases into 1026 normal tension glaucoma (defined by peak IOP of < 22 mmHg) and 808 high tension glaucoma (defined by peak IOP ≥ 22 mmHg). Since peak IOP data was not available in 2146 cases, they were excluded from the subtype analysis.

For the multiethnic validation study, we extracted the association results for 7 SNPs from GWAS data sets of the Singaporean Chinese (1008 cases and 591 controls), the EUR descendants from two independent cohorts (a total of 5008 cases and 35 472 controls) and the African descendants comprising three geographically distinct populations (a total of 2341 cases and 2037 controls). Details from these data sets are summarized in the Supplementary Material Note.

All subjects were enrolled in the study after written informed consent and approval from the relevant national and regional institutional review boards were obtained from the collaborating hospitals in accordance with the Declaration of Helsinki.

Genotyping and quality control

For the GWAS, cases and controls were genotyped using HumanOmniExpressExome BeadChip arrays or a combination of HumanOmniExpress and HumanExome BeadChip arrays (Illumina, Inc., San Diego, CA). To apply stringent quality control (QC), we excluded subjects with a genotype call rate of < 98%, those with close kinship based on PI_HAT values (threshold = 0.1), and outliers from the EAS population estimated by principal component analysis (PCA) using EIGENSTRAT software (ver. 6.0.1). For stringent QC of SNPs, we applied the following filters for autosomal variants: variant call rate ≥ 0.99; P-value for Hardy–Weinberg equilibrium (HWE) in controls ≥ 1.0 × 10⁻⁶ and minor allele frequency (MAF) ≥ 0.5%. We applied QCs for X chromosome variants as follows: variant call rate ≥ 0.99 in both males and females, P-value for HWE ≥ 1.0 × 10⁻⁶ in females and MAF > 0.5%.

For the replication study in the Japanese population, selected SNPs were genotyped in cases from Japanese replication sets using the multiplex polymerase chain reaction (PCR)-based Invader assay (Third Wave Technologies, Madison, WI) (51). The genotypes were called by visual inspection. For QC of the Invader assay, we excluded individuals with a call rate < 95% from further analysis. All SNPs had a call rate ≥ 99%. Details of the genotyping methods and QC are described in a previous report (10). The replication set controls were genotyped using HumanOmniExpressExome BeadChip arrays or a combination of HumanOmniExpress and HumanExome BeadChip arrays (Illumina, Inc.). The same QC filters used for the GWAS were applied for Japanese replication set controls. Genotyping and QC filters used in the multiethnic replication were described in the Supplementary Material Note.

Imputation

EAS samples from the 1000 Genomes Project (phase 1 version 3) were used for imputation reference. As part of QC for the reference panel, we excluded 11 closely related samples based on IBS and used the remaining 275 samples for reference. We also excluded variants with MAFs of < 1.0% or P-values for HWE of ≤ 1.0 × 10⁻⁶. In the GWAS, we phased haplotypes using MACH software. Imputation was carried out by minimac (ver. 0.1.1). For the imputation of X chromosome variants, we phased the haplotypes and performed imputation separately in males and females.

Control sample genotypes from the replication studies were imputed using the same reference panel that was used in the GWAS. The QC procedure for the genotypes used for imputation is identical to that described earlier. We phased the haplotypes using SHAPEITv2 (v2.r778) and performed imputation using Minimac3 (ver. 1.0.13). In the replication studies, best-guess genotypes were used for the association analysis. To evaluate the accuracy of best-guess genotypes, we used the Invader assay to genotype seven lead SNPs at the novel loci using a subsample of replication set 2 controls (3488 out of 9854 controls). The concordance rate of the genotypes was 96.5%.

SNP selection for the replication study

For the Japanese replication study, we initially selected candidate SNPs with P-values of < 1.0 × 10⁻⁵ from the GWAS results. Among these SNPs, we identified 33 lead SNPs with the lowest P-values that were within a 1 Mb genomic region (flanked by 500kb upstream and downstream). Finally, 33 SNPs located at 33 loci were analyzed in the two independent Japanese replication sets.

For the validation study performed in the multiethnic population, we selected 7 SNPs with P-values of < 5.0 × 10⁻⁸ in the combined analysis of the GWAS and the two Japanese replication sets after excluding previously reported loci.

Association analysis

In the GWAS, logistic regression was performed using mach2dat software (version 1.0.19) and assuming an additive genetic model. We used the top 10 principal components as covariates to control for population stratification. To analyze X chromosome variants, we assessed males and females separately and combined the results using the inverse-variance method. In the replication study, we performed a logistic regression analysis using PLINK (ver. 1.9).

The combined analysis of the association results from each population (i.e. Japanese, EAS, EUR and African) was performed using the inverse-variance method and assuming a fixed-effect model. Cochran’s Q test was used to assess heterogeneity across the studies.

General analyses were performed using R software for statistical computing. Regional association plots were generated using LocusZoom (ver. 1.1) (52).

Pathway analysis

MAGENTA software (24) was used for the pathway analysis. To use the LD data in the analysis, we adjusted EAS population variants from the 1000 Genomes Project with PLINK. Genes located at major histocompatibility complex regions were excluded from the analysis. We only considered those pathways with FDR values of < 5% as significant.

Genetic correlation analysis

To evaluate genetic correlations between POAG and systemic traits, we performed cross-trait LD score regression analysis (26) using LDSC software. We used predetermined LD scores provided by the software for the EAS population. GWAS results for these traits were obtained from the BBJ project. Details for each GWAS are shown in Supplementary Material, Table S7. We estimated FDRs using the Benjamini–Hochberg method.

Functional annotation, eQTL and pleiotropy analysis

Variants in LD (r² ≥ 0.8) with the lead SNPs from the newly identified loci were annotated using HaploReg (ver. 4.1) (13). We investigated the overlap between GWAS lead SNPs and eQTLs using two publicly available databases: the Blood eQTL browser (20) and GTEx (release ver. 6) (19). We determined significant overlaps if LD between the top variant of the eQTL and that of GWAS was r² > 0.8 in either EAS or EUR population in the1000 Genomes Project. We only considered eQTLs satisfying FDR < 0.05 in each database. When we determined eQTL top variants, we only considered the variants included in our GWAS.

To analyze pleiotropy, we downloaded the GWAS catalog on 22 May 2017. We defined pleiotropy as when the lead SNP from this study was in LD (r² ≥ 0.8) with those of other traits in either the EAS population or the study population in the report.

Expression analysis in ocular tissues

We defined the positional candidate genes for each locus according to their positions relative to the lead SNPs. For each identified locus, we calculated the range of variants in LD (r² ≥ 0.5) with the lead SNPs using EAS samples from the 1000 Genomes Project, and extending 25 kb upstream and downstream. Positional candidate genes were defined as those which overlapped this predetermined range. For the ANKRD55-MAPK1 locus, the two nearest genes flanking the region of interest were chosen because we could not find any gene that satisfied the criteria described.

Positional candidate gene expression was measured in macaque ocular tissues using qRT-PCR. Ocular tissue from postmortem macaque eyes was obtained from the Tsukuba Primate Research Center at the National Institutes of Biomedical Innovation, Health and Nutrition (Tsukuba, Japan). Total RNA was extracted and cDNA synthesized as previously described (53). RT-PCR was performed using Ex Taq polymerase (Takara Bio, Shiga, Japan) or Go Taq Colorless Master Mix DNA polymerase (Promega Corp., Madison, WI) and gene-specific primers with the following cycle conditions: 94°C for 4 min: then 30–38 cycles of 94°C for 30 s, 50–57°C for 30 s and 72°C for 30 s. The identity of amplified PCR products was verified by electrophoresis on 1.5% agarose gels containing 0.01% GelRed (Biotium, Fremont, CA). The primers used for RT-PCR are available on request.

Data Availability

Individual phenotype data, genotype data and summary statistics that support the findings of this study can be found at NDBC with the primary accession code hum0014 (http://humandbs.biosciencedbc.jp/; date last accessed July 3, 2017). Some access restrictions are applied to the individual data for ethical reasons.

URLs

BioBank Japan Project:https://biobankjp.org/english/index.html,

J-MICC study:http://www.jmicc.com/en/,

JPHC study:http://epi.ncc.go.jp/en/jphc/index.html,

TMM:http://www.megabank.tohoku.ac.jp/english/,

R:https://www.r-project.org/,

PLINK:https://www.cog-genomics.org/plink2,

MACH:http://csg.sph.umich.edu//abecasis/MaCH/,

HaploReg:http://www.broadinstitute.org/mammals/haploreg/haploreg.php,

MAGENTA:https://www.broadinstitute.org/mpg/magenta/,

MANTRA:http://datalib.edina.ac.uk/mantra/softwarepracticals.html/,

LDSC:https://github.com/bulik/ldsc/,

SHAPEIT:https://mathgen.stats.ox.ac.uk/genetics_software/shapeit/shapeit.html/,

Locuszoom:http://locuszoom.sph.umich.edu/locuszoom/,

1000 genome project:http://www.1000genomes.org/,

GWAS catalog:https://www.ebi.ac.uk/gwas/,

ROADMAP Epigenomics Project:http://www.roadmapepigenomics.org/.

All websites were last accessed date on July 3, 2017.

Supplementary Material

Supplementary Material is available at HMG online.

Acknowledgements

We would like to thank all the subjects who participated in this study. We are also grateful to the members of the BioBank Japan Project for supporting the study. We thank Professor Y. Okada (Department of Statistical Genetics, Osaka University Graduate School of Medicine) for his positive criticism during the preparation of this manuscript. We would also like to express our gratitude for the excellent technical support provided by M. Endo, T. Aoi, K. Ashikawa and other staffs in the Laboratory for Genotyping Development at RIKEN Center for the Integrative Medical Sciences.

Conflict of Interest statement. None declared.

Funding

This work was performed as part of the BioBank Japan Project, which is supported by the Japanese Ministry of Education, Culture, Sport, Science and Technology and by the Japanese Agency for Medical Research and Development. The work was also supported by the Japanese Society for the Promotion of Science (KAKENHI; grant numbers 24249083 and 26293372) and by an unrestricted grant from SENJU Pharmaceutical Co., Ltd. (Osaka, Japan) and TOPCON Co., Ltd. (Tokyo, Japan). The Tohoku Medical Megabank project is supported by the Japanese Ministry of Education, Culture, Sport, Science and Technology and the Japanese Agency for Medical Research and Development. The JPHC Study has been supported by the National Cancer Research and Development Fund since 2010 and was supported by a grant in aid of Cancer Research from the Ministry of Health, Labour and Welfare of Japan from 1989 to 2010. The J-MICC Study was supported by grants in aid of Scientific Research in Priority Areas of Cancer (17015018), Innovative Areas (221S0001) and a JSPS KAKENHI grant (16H06277) from the Japanese Ministry of Education, Culture, Sport, Science and Technology.

References

Author notes