TMEM132D and VIPR2 Polymorphisms as Genetic Risk Loci for Retinal Detachment: A Genome-Wide Association Study and Polygenic Risk Score Analysis

Abstract

Purpose: Retinal detachment (RD) is a sight-threatening ocular disease caused by separation of the neurosensory retina from the underlying retinal pigment epithelium layer. Its genetic basis is unclear because of a limited amount of data. In this study, we aimed to identify genetic risk loci associated with RD in participants without diabetes mellitus and to construct a polygenic risk score (PRS) to predict the risk of RD.

Methods: A genome-wide association study was conducted using data from the Taiwan Biobank to identify RD risk loci. A total of 1533 RD cases and 106,270 controls were recruited, all of whom were Han Chinese. Replication studies were performed using data from the UK Biobank and Biobank Japan. To construct the PRS, a traditional clumping and thresholding method was performed and validated by fivefold cross-validation.

Results: Two novel loci with significant associations were identified. These two genes were TMEM132D (lead single nucleotide polymorphism [SNP]: rs264498, adjusted-P = 7.18 × 10⁻⁹) and VIPR2 (lead SNP: rs3812305, adjusted-P = 8.38 × 10⁻⁹). The developed PRS was effective in discriminating individuals at high risk of RD with a dose-response relationship. The quartile with the highest risk had an odds ratio of 1244.748 compared to the lowest risk group (95% confidence interval, 175.174–8844.892).

Conclusions: TMEM132D and VIPR2 polymorphisms are genetic candidates linked to RD in Han Chinese populations. Our proposed PRS was effective at discriminating high-risk from low-risk individuals.

Retinal detachment (RD) is a sight-threatening ocular disease with an annual incidence of approximately 1 in 10,000.¹ It is characterized by separation of the neurosensory retina from the underlying retinal pigment epithelium (RPE). If untreated, RD may result in cataracts, glaucoma, proliferative vitreoretinopathy, and complete blindness.² Prompt surgery is indicated for rhegmatogenous retinal detachment (RRD) and tractional retinal detachment (TRD), with the success rate ranging from 75% to more than 95% according to severity.³ Individuals who are older, have a giant retinal tear, wide retinal detachment, and macula-off detachment are less likely to have restoration of good sight.⁴ For high-risk individuals with retinal break, preventive strategies include prophylactic laser treatment.⁵ 

RD is classified into RRD, TRD, and exudative retinal detachment. The annual incidence of RRD in Taiwan is approximately 16.4 per 100,000 population.⁶ It starts with a retinal break followed by ingress of the vitreous body into the subretinal space, thus acquiring the name “rhegmatogenous” (from Greek “rhegma” meaning discontinuity or a break). TRD is a much less common type of RD. It is caused by the tractional force of membranes in the vitreous or over the retina that separates the retina and is most commonly seen in patients with diabetic retinopathy. 

Because of small sample sizes in previous studies, the genetics of RD is poorly understood. Familial studies have found that siblings and offspring of patients with RRD have twice the risk of developing the disease independently of age, sex, and myopia status.^7,8 In 2013, the first genome-wide association study (GWAS) was performed with 867 Scottish patients with RRD, and the results identified only one significant single nucleotide polymorphism (SNP), namely rs267738 in the CERS2 gene.⁹ Another study suggested that intronic variants in COL2A1 conferred a higher risk of RRD along with Stickler syndrome.¹⁰ In 2019, a GWAS performed on 3977 cases recruited from the UK Biobank (UKB) identified six risk loci: FAT3, COL22A1, TYR, BMP3, ZC3H11B, and PLCE1.¹¹ However, a limitation of the study was that it did not elucidate differences in genetics among different ethnicities. Hence, further studies are needed to complete the genetic landscape of RD. 

In this study, we present the largest pure Han Chinese GWAS of RD to date and propose the first polygenic risk score (PRS) model to predict the risk of RD using data from the Taiwan Biobank (http://www.twbiobank.org.tw/, accessed on July 2022), a databank with more than 180,000 ethnically Han Chinese enrollees. Through the GWAS, we identified two genetic loci associated with RD, namely TMEM132D and VIPR2. In addition, our PRS successfully stratified the risk of RD according to an individual's genetic profile. The quartile with the highest risk had an odds ratio of 1244.748 compared to the lowest-risk group. 

Method

Participants, Phenotyping, and Genotyping

Data on all participants in our GWAS and subsequent PRS construction were obtained from the Taiwan Biobank (TWB), which includes extensive data on more than 180,000 volunteers from the Taiwanese population (https://www.biobank.org.tw/, accessed in July 2022). The TWB primarily uses a self-reported mechanism to collect phenotypic data. A well-trained researcher conducts a structured interview with each participant once they agree to be enrolled in the TWB, and then the researcher fills out a questionnaire with information about the participant's demographics, lifestyle behaviors, environmental exposures, dietary habits, family history, and health-related information. In addition, all participants provide blood samples that are used for genetic studies.¹² Genome-wide SNP data were generated using a custom Axiom Genome-Wide Array Plate system, with more than 653,000 SNPs in the TWBv1.0 array and 750,000 SNPs in the TWBv2.0 array. This study was conducted according to the guidelines of the Declaration of Helsinki and was approved by the Institutional Review Board of Taipei Veterans General Hospital (ID no. 2020-07-008A and 2020-12-009AC). Informed consent was obtained from all subjects involved. 

Because comorbid diabetes mellitus (DM) represents a major distinctive entity in the pathophysiology of TRD and may distort genetic results, we excluded all participants with self-reported DM status in this study. Cases were defined as individuals with self-reported RD, whereas controls were defined as participants without self-reported RD. We performed quality control to exclude imputed SNPs from the TWB array with genotyping rates <95% and to remove SNP variants with a missing call rate >5%, minor allele frequency <1%, or deviation of heterozygosity with Hardy-Weinberg equilibrium P < 1.0 × 10⁻⁶ using PLINK version 1.9 (www.cog-genomics.org/plink/1.9/).¹³ For sample filtering, PLINK v1.9 software was also used to identify samples with genetic relatedness, indicating that they were from the same individual or first-, second-, or third-degree relatives. These determinations were based on evidence of cryptic relatedness from identity-by-descent status (pi-hat cutoff of 0.125). After removing first-, second-, and third-degree relatives, 107,803 independent samples (1533 RD cases and 106,270 controls) and 5,687,850 SNPs that passed quality control were used for subsequent association analysis and to construct the PRS model. 

Association Analysis

Genome-wide association analysis was performed using a logistic regression model with PLINK v1.9 under the assumption of additive allelic effects. The outcome variables were RD case and control statuses. We adjusted for age and sex as covariates and used a conventional genome-wide significance threshold of P < 5 × 10⁻⁸ to identify significant SNPs. 

To prevent bias resulting from an imbalanced case-control ratio, we performed a Scalable and Accurate Implementation of Generalized mixed model (SAIGE) using the SAIGE package in R (version 1.0.0).¹⁴ Furthermore, we replicated our results in two independent cohorts: Biobank Japan (hum0197.v3.BBJ.RD.v1, https://pheweb.jp/pheno/Retinal_Detachment) and the UKB (phenocode-20002_1281, https://pheweb.org/UKB-Neale/pheno/20002_1281). 

The Manhattan plot of the adjusted GWAS results was plotted using the qqman package (https://github.com/stephenturner/qqman).¹⁵ Locus zoom plots for gene loci that reached the genome-wide significance level were illustrated using the LocusZoom tool (http://csg.sph.umich.edu/locuszoom).¹⁶ 

Construction of the PRS Model

To estimate an individual's genetic risk of developing RD, we constructed PRS models with the clumping and thresholding method using PLINK v1.9. Hyperparameters included clump-P, namely the P value threshold for an SNP to be included as an index. We tried five clump-P values, namely P < 0.00025, P < 0.0005, P < 0.00075, P < 0.001, and P < 0.05. For each clump-P, SNPs within 25 Mb of the index SNP were considered for clumping, and SNPs with r² > 0.1 (we tried r² = 0.01 and 0.001) with the index SNPs were removed. For SNPs with a negative variant logistic regression coefficient, we mathematically converted their minor alleles into major alleles, thereby resulting in positive weight values for all variants. The PRS was calculated as the summation of the coefficients of logistic regression for all variants.¹⁷ Fivefold cross-validation was used in the model construction for each clump-P.¹⁸ We split the TWB dataset into five equal and mutually exclusive subsets using R software. Each subset was used for testing once, with the other four subsets used to train the model. The process was then repeated five times. The model with the best predictive performance, as measured by the area under the receiver operating characteristic curve (AUC of the ROC curve) using the R package “pROC” was chosen as the final model. 

To demonstrate the genetic risk of RD predicted by the PRS model, the participants were split into four groups based on their PRS. The (min, Q1) group included individuals with PRS values less than the lowest 25%. The (Q1, Q2) group included individuals with PRS values between the lowest 25% and the median. The same rule applied to the other two groups. The odds ratios of acquiring RD were calculated in each group, with the (min, Q1) group defined as the reference group. 

Results

Participant Profiles

We investigated the genetic features of RD in 1533 RD cases and 106,270 controls recruited from the TWB. None of the participants had DM because it contributes to RD with different pathogenesis as mentioned above.¹⁹ The participants’ demographics and clinical characteristics are shown in Supplementary Table S1. 

Genetic Features of RD

After adjusting for age and sex, the genetic features of RD were further investigated (Supplementary Table S2). Two significant genomic loci were identified (Fig. 1). One locus, located on chromosome 12 (q24.33), was within the transmembrane protein 132D (TMEM132D) gene. The lead SNP of this locus was rs264498 (adjusted-P = 5.92 × 10⁻⁹), with a minor allele frequency of 0.1739 in the cases and 0.1349 in the controls. The odds ratio for each adverse allele was 1.354 (standard error = 0.052). Another locus was identified at chromosome 7 (7q36.3), spanning across the vasoactive intestinal peptide receptor 2 (VIPR2) gene, with the major SNP rs3812305 (adjusted-P = 8.38 × 10⁻⁹, odds ratio = 1.408). Locus features of the two SNPs are shown in Figure 2. 

Figure 1.

 

Manhattan plot showing the association of 5,687,850 single nucleotide polymorphisms (SNPs) with retinal detachment in the Taiwan Biobank genome-wide association study. The −log10 (adjusted P) values for the 5,687,850 SNPs are shown against their corresponding chromosomes and at their genomic position. The horizontal red and blue lines are adjusted P values of 5.0 × 10⁻⁸ and 1.0 × 10⁻⁵, respectively. The two gene regions with their leading SNPs that reached a genome-wide significance level are labeled.

Figure 2.

 

Locus plots of the TMEM132D gene (A) and VIPR2 gene (B). The leading SNPs for the two genes were rs264498 and rs3812305, respectively. Recombination rates are shown in blue on the y-axis. The genomic positions are plotted on the x-axis. R² indicates the correlation coefficient.

Inflation of type I error is commonly seen in an unbalanced case-control ratio in biobank studies. To confirm the relatedness of SNPs and RD, we performed SAIGE analysis (Supplementary Table S2, column P.saige). The aforementioned two significant SNPs, after SAIGE adjustment, still passed the genome-wide significance threshold of P < 5 × 10⁻⁸, with a P value of 5.92 × 10⁻⁹ for rs264498 and 7.66 × 10⁻⁹ for rs3812305. Furthermore, we replicated the association analysis using data from the UKB and Biobank Japan, and several SNPs in VIPR2 were successfully replicated with P < 0.05 in the UKB (Supplementary Table S2, in the two rightmost columns). 

The Ability of the PRS To Predict Individuals With RD

To assess the RD risk in individuals without DM, we constructed PRS models based on the GWAS results. Five PRS models with different tuning parameters are shown in Table 1. The mean PRSs were significantly greater in the RD cases compared to the controls across all models. Considering the clinical significance, P value, and AUC in the training and testing datasets, the PRS model with 661 SNPs was selected as our final model. Using this model, the mean PRS in the cases was 0.0086 compared to 0.0032 in the controls. 

Table 1.

 

Comparison of the Predictive Performance of the PRS With Different Tuning Parameters

 

This PRS model was excellent in discriminating high-risk RD patients with a dose-response relationship. The participants were then divided into four quartiles according to their PRS values (Table 2). In the training dataset, compared to individuals in the lowest PRS quartile (min to Q1), those in the second lowest PRS quartile (Q1 to Q2) had a 7.002-fold higher risk of developing RD (95% confidence interval [CI], 0.861–56.921). Furthermore, individuals in the second highest PRS quartile (Q2 to Q3) had a 35.065-fold higher risk (95% CI, 4.804–255.978), and those in the highest quartile (Q3 to max) had a 1244.748-fold higher risk (95% CI, 175.174–8844.892). In the testing subset one-fifth of the TWB dataset was used, and individuals in the highest and second highest PRS quartiles still had significantly higher RD risks compared to those in the lowest PRS quartile (odds ratio = 30.750 [95% CI, 15.322–72.264] and 5.921 [95% CI, 2.768–14.520], respectively). Illustrations of the distribution of PRS among the cases and controls are shown in Supplementary Figure S1. 

Table 2.

 

Distribution of RD Cases and Controls Regarding PRS Quartiles in 661 SNPs

 

We also calculated the odds ratio of developing RD in the high-risk group (the top 5%–25% of PRS) compared with the remaining individuals (Supplementary Table S3). Compared to the individuals in the bottom 75% of PRS, those in the top 25% had an 88.844 (95% CI, 65.275–120.923) higher risk of developing RD. However, the odds ratio did not increase when comparing the top 10% to the remaining 90% or comparing the top 5% to the remaining 95%, manifesting the very high risk of developing RD in the high PRS groups regardless of their absolute PRS value. 

The performance of the PRS model was also evaluated according to AUCs (Fig. 3). Predicted by age and sex alone, the AUC was only 0.5119. When genetic information was added via the PRS, the AUC increased to 0.8370. The effectiveness of the PRS was also replicated in the testing dataset, with the AUC of clinical variables alone being 0.5212, increasing to 0.6710 with the addition of the PRS. 

Figure 3.

 

Receiver operating characteristic (ROC) curves for the PRS models in predicting retinal detachment development. Green: using PRS alone; Red: using clinical covariates (COV) of age and sex alone; Blue: using the combination of PRS and COV.

Discussion

In this study, we identified two independent loci for RD using GWAS, namely TMEM132D and VIPR2. In addition, we developed a PRS model capable of predicting the risk of RD in individuals without DM. 

The retina can be regarded as an extension of the central nervous system; hence, molecular pathways in the nervous system may resemble those in the retina.²⁰ We identified a strong association between TMEM132D polymorphisms and RD. TMEM132 is a transmembrane protein family containing five proteins (TMEM132A, B, C, D, and E). The extracellular portions of these proteins contain cohesin and immunoglobulin domains, indicating their function as neural adhesion molecules that bring together the extracellular matrix and intracellular cytoskeleton.²¹ Although functions of the TMEM132 family of genes remain mostly unknown, they are differentially expressed in the developing forebrain; TMEM132A, C, and E are also expressed in early neural progenitors and developing cochlea.²² With regards to TMEM132D, its polymorphisms have been linked to various psychiatric diseases, including panic, anxiety, and major depressive disorders, and also frontotemporal dementia.^23–25 The mechanism through which risk variants influence anxiety-related behavior depends on alterations of the binding affinity of the RNA polymerase II complex to TMEM132D promoter and the downstream gene expression level.²⁵ In addition, the development and anomalies of human eyes have been shown to be associated with CpG methylation.^26,27 Given that methylation of CpG2 of TMEM132D has been shown to orchestrate the relationship between physical abuse and panic disorder, the epigenetic status of TMEM132D is worth exploring in RD.²⁸ Furthermore, TMEM132D has been shown to mediate neuronal morphogenesis through the WAVE-regulatory complex.²⁹ A component of the WAVE-regulatory complex, NAP1, has been shown to interact with the cytoplasmic domain of TMEM132D, and to control Arp2/3 complex-mediated actin polymerization.³⁰ In addition, the gain-of-function of TMEM132D has been shown to suppress the WAVE-regulatory complex, impair actin nucleation, and hamper cell motility in mammalian cells.²⁹ Moreover, after RD, actin depolymerization has been shown to prompt the axonal retraction of rod photoreceptors and subsequently break the first synapse in the visual pathway via the RhoA-LIMK-Cofilin signaling pathway.³¹ 

We also identified an association between VIPR2 polymorphisms and RD. VIPR2 (or VPAC2) is expressed in the human retina and RPE.¹⁶ It is a G-protein coupled receptor (GPCR) that responds to vasoactive intestinal peptide (VIP) and pituitary adenylate cyclase-activating polypeptide, and transduces signals via an increase in adenylyl cyclase activity and the subsequent activation of protein kinase A.^32,33 VIP has been shown to promote the proliferation, differentiation, and development retinal pigment epithelium.³⁴ In addition, VIP has been shown to stimulate melanogenesis in the RPE in a concentration-depend manner and promote fluid transportation.^34,35 Furthermore, VIPR2 is associated with pachychoroid diseases, circadian rhythm, and schizophrenia.^36,37 Several reports focusing specifically on the Han Chinese population, including a genome-wide meta-analysis, have also identified an association between VIPR2 with high myopia.^38,39 However, these studies did not investigate the association between RD and VIPR2, and the lead SNPs identified were also different. In animal models, dysfunction of the VIP-VIPR2 axis has been shown to impair bipolar cell function and induce myopia.⁴⁰ Pathologic myopia is a risk factor for RD, as the retina in myopic patients is thinner and thus more prone to a retinal break and RD development.^41,42 In addition, VIPR2 has been shown to promote pachychoroid thickness and the development of central serous chorioretinopathy, which are also diseases arising from the choroid-RPE-neurosensory retina complex.⁴³ In the present study, we identified an association between VIPR2 polymorphisms and RD. Considering the other aforementioned studies on VIPR2, our findings suggest either a direct or indirect link between VIPR2 polymorphisms and RD. However, the interplay between the gene, RD, and myopia remains to be further explored. 

In this study, we also attempted to replicate our findings using data from the UKB databank. Two SNP loci in VIPR2 returned significant P values, namely rs76583299 (P = 0.018) and rs79725155 (P = 0.021). Because the UKB dataset has a much higher ethnic diversity compared to the TWB database used in our study, the results suggested that confounding factors caused by a difference in genetic landscapes may play a significant role in the pathogenesis of RD. Furthermore, we constructed a PRS model using 661 SNPs associated with RD, and the results showed that the odds ratio of RD in the (Q3, max) quartile was 1244.748 compared to the (min, Q1) quartile, effectively differentiating the high-risk from the low-risk individuals. Previously, prophylactic measures for RD, such as laser coagulation, have been recommended exceptionally in individuals presenting with clinical risk factors such as retinal breaks or retinal degeneration.⁴⁴ A literature review also showed a lack of evidence on this issue.⁴⁵ Our proposed scoring system provides a new way of looking at the risks of RD from a genetic standpoint, suggesting the possibility of investigating prophylactic measures for high-risk RD groups in the future. In addition, it enables the possibility to adjust for genetic confounding effects in future studies on RD. Moreover, use of the PRS for risk assessment may enable selective health education for high-risk groups, further enhancing risk awareness in high-risk individuals. Ultimately, the PRS model may also open up the possibility for precision medicine for RD based on genetic risks. 

The strengths of this study lie in the homogeneity and size of the dataset from the TWB. Similar previous studies have either used smaller datasets for analysis or the dataset has included various ethnicities.^9,11 In contrast, the TWB only recruits individuals of Han Chinese ethnicity, with a moderate case size and a homogeneous ethnic data profile. Therefore our study significantly reduces ethnicity-related confounding. Furthermore, this study represents the first large-scale genome-wide RD study in Han Chinese. 

There are also several limitations to this study. Because this is a moderate-sized retrospective study with self-reported cases, problems such as selection bias, recall bias, as well as nongenetic confounding factors may be present. In addition, we could not really differentiate TRD or exudative retinal detachment from RRD in this study. We tried to use the presence of DM as a way to exclude most TRD cases instead of classifying them through the actual RRD diagnosis. We believe that the proportion of participants with TRD or exudative retinal detachment enrolled in this study should be low. In addition, we did not screen all of the controls with fundoscopy; therefore some of the controls may have had lattice degeneration or retinal breaks and may have had RD in the following years. Furthermore, relying on genetic risk alone is inadequate, and confirmation through clinical fundus features is indispensable. Despite these limitations, this study may serve as the foundation for future, larger, prospective investigations on this topic. 

In conclusion, we identified TMEM132D and VIPR2 polymorphisms as genetic candidates associated with RD in a Han Chinese population. Our PRS model was successful in differentiating individuals with a high risk of RD from those with a low risk. Our results may provide a framework for future, larger, prospective studies on RD genetics, as well as studies on RD in general. 

Acknowledgments

The authors thank the Big Data Center (BDC) of Taipei Veterans General Hospital (VGHTPE), and the Department of Statistics, Tamkang University for technological support. 

Supported by the National Science and Technology Council (NSTC 111-2314-B-075-036-MY3, and NSTC 112-2321-B-A49-007), Taipei Veterans General Hospital (111VACS-007 and V112C-187), the Ministry of Education, Higher Education SPROUT Project for Cancer Progression Research Center and Cancer and Immunology Research Center; the “Center for Intelligent Drug Systems and Smart Bio-devices (IDS2B)” from The Featured Areas Research Center Program within the framework of the Higher Education Sprout Project by the Ministry of Education in Taiwan. 

Disclosure: H.-K. Chuang, None; A.-R. Hsieh, None; T.-Y. Ang, None; S.-W. Chen, None; Y.-P. Yang, None; H.-J. Huang, None; S.-H. Chiou, None; T.-C. Lin, None; S.-J. Chen, None; C.-C. Hsu, None; D.-K. Hwang, None