Molecular Vision 2019; 25:222-236
<http://www.molvis.org/molvis/v25/222>
Received 11 July 2018 | Accepted 21 April 2019 | Published 23 April 2019
Molecular Vision 2019; 25:222-236
<http://www.molvis.org/molvis/v25/222>
Received 11 July 2018 |
Accepted 21 April 2019 |
Published 23 April 2019
Ye Lu,1 Diana Zhou,2 Hong Lu,3 Fuyi Xu,2 Junming Yue,4 Jianping Tong,1 Lu Lu2
The first two authors contributed equally to this study.
1Department of Ophthalmology, The First Affiliated Hospital, Zhejiang University College of Medicine, Hangzhou, China; 2Department of Genetics, Genomics and Informatics, University of Tennessee Health Science Center, Memphis, TN; 3Department of Ophthalmology, Nantong Eye Institute, Affiliated Hospital of Nantong University, Nantong, China; 4Department of Pathology, University of Tennessee Health Science Center, Memphis, TN
Correspondence to: Lu Lu, Department of Genetics, Genomics and Informatics, University of Tennessee Health Science Center, Memphis, TN 38163; Phone: (901) 448-7557; FAX: (901) 448-3500; email: lulu@uthsc.edu
Purpose: Glaucoma is characterized by optic nerve damage and retinal ganglion cell loss. The glycoprotein neuromedin B-associated (Gpnmb) gene is well-known to be involved in the glaucoma disease process. The purpose of this study is to identify a downstream gene through which Gpnmb affects the glaucoma phenotypes using a systems genetics approach.
Methods: Retinal gene expression data for the BXD recombinant inbred (RI) strains (n=75) have previously been generated in our laboratory for a glaucoma study, and these data were used for genetic and bioinformatics analysis. Expression quantitative trait locus (eQTL) mapping and genetic correlation methods were used to identify a gene downstream of Gpnmb. Gene-set enrichment analysis was used to evaluate gene function and to construct coexpression networks.
Results: The level of Gpnmb expression is associated with a highly statistically significant cis-eQTL. Stanniocalcin 1 (Stc1) has a significant trans-eQTL mapping to the Gpnmb locus. The expression of Gpnmb and Stc1 is highly correlated in the retina and other tissues, as well as with glaucoma-related phenotypes. Gene Ontology and pathway analysis showed that Stc1 and its covariates are highly associated with apoptosis, oxidative stress, and mitochondrial activity. A generated gene network indicated that Gpnmb and Stc1 are directly connected to and interact with other genes with similar biologic functions.
Conclusions: These results suggest that Stc1 may be a downstream candidate of Gpnmb, and that both genes interact with other genes in a network to develop glaucoma pathogenesis through mechanisms such as apoptosis and oxidative stress.
Glaucoma is currently the second leading cause of blindness in the United States, and affects millions worldwide. In 2010, approximately 2.22 million people were affected by open-angle glaucoma (OAG) in the United States, and more than 70 million people worldwide had glaucoma [1-4]. Glaucoma is a genetically heterogenous group of eye diseases that results in optic nerve damage and retinal ganglion cell (RGC) loss usually associated with elevated intraocular pressure (IOP) and vision loss. Unfortunately, most cases of glaucoma are asymptomatic until after sight is lost [5]. Glaucoma damage is irreversible, and therapeutic interventions are effective only at early stages of the disease.
Pigmentary glaucoma (PG), a form of secondary open-angle glaucoma, is the second most common form of glaucoma in young adults. In this disease, pigment sloughs off the posterior iris and damages the drainage structure of the eye, which, in turn, attenuates or blocks the flow of aqueous humor [6]. In most cases, the resultant IOP increase leads to glaucoma, yet in some cases, the increase does not, indicating that additional factors, including environment, age, and genetics, may influence disease progression [7]. The exact pathogenesis of glaucoma is not fully understood. However, RGC damage, oxidative stress, ischemia and hypoxia, abnormal immunity and inflammatory reactions, and apoptosis have all been found to play a part in the progression of glaucoma [8-13].
The glycoprotein neuromedin B-associated (Gpnmb; Gene ID: 93695, OMIM: 604368) gene is a gene that is known to affect pigmentary glaucoma. Gpnmb is associated with iris pigment dispersion (IPD), which is the breakdown of the posterior iris pigment epithelium [6]. Mutations in Gpnmb and tyrosinase-related protein 1 (Tyrp1; Gene ID: 22178, OMIM: 115501) have been found to contribute to IPD and iris stromal atrophy (ISA) in DBA/2J (D2) mice [14]. D2 mice are homozygous for mutations in Gpnmb and Tyrp1, and have a higher chance of developing pigmentary glaucoma (PG) [15-17]. The eyes of D2 mice have no obvious abnormalities before 3 months of age. However, by 6 months, most D2 mice begin to develop iris depigmentation, increased IOP, and damage to the optic nerve head (ONH) [18]. Because the progression of glaucoma in D2 mice is similar to that of humans, this strain has been used as a mouse model to study pigment dispersion syndrome (PDS) and PG [18]. Using these mice, Gpnmb has been found to interact with a large network of genes to affect the presentation of glaucoma [19]. Despite these findings, the direct downstream pathway through which Gpnmb contributes to glaucoma remains unknown.
The purpose of this investigation is to identify genes and pathways through which Gpnmb works to affect the glaucoma phenotype, using a systems genetics method and BXD mouse genetic reference population (GRP), one parental strain of which is D2 mouse. After expression quantitative trait locus (eQTL) mapping and genetic correlation analysis, we identified stanniocalcin 1 (Stc1) as a downstream candidate gene of Gpnmb with high variable expression among the BXD mice. The two parental strains of BXD GRP, D2 and B6 mice, express Gpnmb and Stc1 statistically significantly differently, which provides us with an opportunity to study expression variance of Gpnmb and Stc1 among the BXD population. Using retinal gene expression data of BXD mice, we identified a genetic network involving Gpnmb and Stc1, and explored pathways through which the two genes interact to affect the glaucoma-related phenotypes.
GeneNetwork website is a public web source that houses phenotype and transcriptome data from various tissues of BXD RI mice, including the eye and the retina. The gene expression data used in this study can be accessed using the data set “Full HEI Retina Illumina V6.2 (Apr10) RankInv” that we generated through collaborative effort for a glaucoma study [20]. This data set uses data from 75 BXD strains, their parental strains B6 and D2, and both reciprocal F1 hybrids between B6 and D2. Almost all mice were young adults of 60 to 90 days of age. All animals were housed with free access to standard laboratory chow and water, maintained on a 12 h:12 h light-dark cycle, and were killed via rapid cervical dislocation as described previously [21]. This study adheres to the ARVO Statement for Use of Animals in Research. All experimental protocols were approved by the Institutional Animal Care and Use Committee (IACUC) at the University of Tennessee Health Science Center (Memphis, TN). More detailed information on this data set, including individual mouse, tissue harvest, RNA extraction, microarray hybridization, data normalization, etc., can be found at GeneNetwork267.
The heritability of Gpnmb and Stc1 expression was calculated using broad sense heritability that compares the genetic variation between strains to the environmental variance within strains [22]. The formula used for heritability calculation is 0.5Vg / (0.5Vg + Ve), where Vg is genetic variance (variances of strain means), and Ve is environmental variance. The 0.5 factor in this ratio was applied to adjust for the twofold increase of additive genetic variance among inbred strains relative to outbred populations [23].
Expression QTL (eQTL) mapping was performed using the WebQTL module on our GeneNetwork website with our published methods [21,24,25]. Regression analysis was used to determine the relation between differences in an expressed trait and differences in alleles at markers across the genome. Simple interval mapping identified potential eQTLs that regulate Gpnmb and Stc1 expression levels and estimated the statistical significance at each location using known genotype data for those sites. Composite interval mapping was also performed to control for genetic variance associated with major eQTLs, and therefore, identify any secondary eQTLs that may have been otherwise masked. Each analysis produced a likelihood ratio statistic (LRS) score, providing a quantitative measure of confidence of linkage between the observed phenotype—in this case, variation in the expression levels of Gpnmb and Stc1—and known genetic markers. The genome-wide significance for each eQTL was established using a permutation test that compared the LRS of our novel site with the LRS values for 1,000–10,000 genetic permutations [26].
Downstream genes, or genes whose expressions are affected by differential expression of Gpnmb, were identified using multiple criteria. First, the expression of the candidate downstream gene should be significantly correlated with the expression of Gpnmb by evaluation in both retina and other tissues. The candidate gene should also have an expression level greater than 8 units in retina tissue. The candidate gene should have a significant trans-eQTL at the Gpnmb locus, located in the 10 Mb interval between 44 Mb and 54 Mb of Chromosome (Chr) 6 where the Gpnmb gene is located. The gene expression should be represented by probes that hybridize in exons or 3’-UTR region. The expression of the candidate gene should be highly correlated with eye related phenotypes. Finally, published literature must support the gene’s biological function in the retina or eye.
Genetic correlation analysis was performed by computing Pearson product-moment correlations of the strain means between the expression of Gpnmb and the expression of all other transcripts across the mouse genome to produce a set of genetically correlated genes. Genes with an expression level greater than eight and a statistically significant correlation with Gpnmb (p<0.05) were selected for further analysis. Tissue correlation is an estimate of the similarity of expression of two genes across different tissues or organs from a single individual, which eliminates any genetic or environmental effect on the gene expression difference, and helps to find a true biologic connection between any pair of genes. Tissue correlation was performed by computing Spearman rank-order correlations between the expression of Gpnmb and its downstream candidates across 20 more multiple different tissues. Literature correlation was performed by examining the correlation coefficient (r) value of genes already described using similar terminology in published papers. The Semantic Gene Organizer was used to perform literature correlation, thus further filtering for a biologic correlation between Gpnmb and other genes based on literature reports [27]. Genes with r >0.5 were considered to have a high literature correlation. Riken clones, intergenic sequences, predicted genes, and probes not associated with functional mouse genes were eliminated to create a list of genes statistically significantly correlated with Gpnmb. Later, this process was used with the identified Gpnmb downstream gene Stc1 to calculate its highly correlated genes. These computations were all performed using tools on GeneNetwork.
The GeneNetwork website contains extensive phenotypic data sets ranging from behavioral to morphological to pharmacological. We queried the BXD published phenotypes database in GeneNetwork for eye-related traits, and focused the analysis on the traits that were statistically significantly correlated with Gpnmb and Stc1 expression in the retina (p<0.05). Then, we conducted 20,000 permutation tests as a more random and accurate measure of significance, and all traits that were no longer statistically significantly correlated after the permutation tests were eliminated.
To investigate Stc1 function and gene pathways through which Stc1 might play a role in retinal-related diseases, the top 500 genes with statistically significant genetic correlations (p<0.01) and literature correlations (r >0.58) with Stc1 were uploaded to the Webgestalt website for Gene Ontology (GO) and pathway analysis [28]. The p values generated from the hypergeometric test were automatically adjusted to account for multiple comparisons using the Benjamini and Hochberg correction [29]. The categories with an adjusted p value (adj P) of less than 0.05 indicated that the set of submitted genes were statistically significantly over-represented in those categories.
The gene network was constructed and visualized using the Cytoscape utility through the Gene-set Cohesion Analysis Tool (GCAT). Nodes in the network represent genes, and edges between two nodes represent the cosine score of latent semantic indexing (LSI) that determines whether the functional coherence of the gene sets is larger than 0.6. The significance of the functional cohesion is evaluated by the observed number of gene relationships above a cosine threshold of 0.6 in the LSI model. The literature p value (LP) is calculated using Fisher’s exact test by comparing the cohesion of the given gene set to a random one [27].
The expression levels of Gpnmb and Stc1 gene were verified with quantitative reverse transcription PCR (qRT-PCR). Retina tissue from D2 mice were harvested at 2–3 months (before the onset of glaucoma) and 6–7 months (after the onset of glaucoma) of age (≥3 samples per group). Total RNA was extracted with the RNeasy Mini Kit (Qiagen, Germantown, MD) following the manufacturer’s instructions. RNA concentration and purity were evaluated with NanoDrop (Thermo Scientific, Waltham, MA). The qRT-PCR reactions (20 μl total volume) were performed using the iTaq Universal SYBR Green One-Step Kit (Bio-Rad, Hercules, CA), including 4 μl RNA (20 ng/μl), 1.125 μl primers (5 μM), 7.5 μl 2X SYBR Green RT–PCR Reaction Mix, 0.3 μl iScript Reverse Transcriptase, and 2.075 μl Nuclease-free H2O. The RT step involved incubation at 50 °C for 10 min. The PCR cycling conditions included an initial denaturation of 95 °C for 5 min followed by 40 cycles of 95 °C for 10 s and 60 °C for 30 s. The qRT-PCR reaction was performed on the CFX Connect Real-Time System (Bio-Rad). The sequences of the PCR primers were as follows: Gpnmb, forward 5′-GCT GGT CTT CGG ATG AAA ATG A-3′, reverse 5′-CCA CAA AGG TGA TAT TGG AAC CC-3′; Stc1, forward 5′-ACG AGG CGG AAC AAA ATG ATT-3′, reverse 5′-TGC ACT TTA AGC TCT CTT TGA CA-3′; Gapdh forward 5′-GGA GCC AAA AGG GTC ATC AT-3′, reverse 5′-GTG ATG GCA TGG ACT GTG GT-3′. The relative expression of each gene was normalized to Gapdh by using the ΔΔ2Ct method, and data were presented as mean ± standard deviation (SD).
Gpnmb has highly variable expression among the BXD strains (Figure 1). The average expression across all BXD strains was 10.871±0.240 (mean ± standard error [SE]), with the parental B6 strain having the highest expression at 15.074±0.381 and the parental D2 strain having one of the lowest expressions at 7.560±0.473, indicating a statistically significant difference between the parental strains (p<0.01). The expression levels of most BXD strains are between B6 and D2 mice. The heritability of Gpnmb expression is 0.589, which indicates that genetic variability affects Gpnmb expression, and allows us to identify a locus that controls Gpnmb expression among BXD mice.
Simple interval mapping for Gpnmb revealed a statistically significant eQTL with an LRS of 43.0 on chromosome 6 at 49.66 megabases (Mb), which is close to where the Gpnmb gene itself is located chromosome 6 at 49.06 Mb (Figure 2). Composite interval mapping revealed no secondary locus that could modulate Gpnmb expression. This indicates that Gpnmb is cis-regulated, meaning that most of the variation in Gpnmb expression is caused by sequence variants in or near Gpnmb. Using our open access sequence data resources at GeneNetwork, we identified 37 single nucleotide polymorphisms (SNPs) in Gpnmb between the BXD parental strains (Table 1). Four SNPs are located in the coding region, including one nonsynonymous SNP, and the rest are located in the intron area. One or several of these SNPs are responsible for Gpnmb expression differences in the BXD mice.
Genetic mapping showed 11 statistically significant trans-eQTLs that were located within 5 Mb of Gpnmb (chromosome 6 from 44 to 54 Mb). After filtering with genetic correlation (p<0.05) and expression level (greater than 8), four genes were found to be statistically significantly correlated with Gpnmb and highly expressed in the retina. After further filtering with tissue correlation analysis, only Stc1 was found to be statistically significantly correlated with Gpnmb (p=0.0009, Figure 3). To validate expression of Gpnmb and Stc1 in the retina, we performed qRT-PCR for D2 mice before and after the onset of glaucoma. The results showed that expression of Gpnmb in D2 mice was decreased by approximately 60% in 6- to 7-month-old D2 mice when compared with 2- to 3-month-old D2 mice (p=0.039; Figure 4). For Stc1, the expression level decreased 20% in the 6- to 7-month-old D2 mice (p=0.059). In addition, correlation analysis revealed that the expression level of Stc1 was statistically significantly correlated with Gpnmb (p=0.044). Previous literature was consulted, and Stc1 was found to be the only gene with a biologic function connected to the retina [30]. Stc1 was chosen as the best downstream candidate gene of Gpnmb for further analysis.
Stc1 has highly variable expression among the BXD strains with the highest expression at 10.074±0.1480 in BXD 83 and the lowest expression at 8.715±0.089 in BXD51, indicating a 2.7-fold change in expression range (Figure 5). The parental B6 strain had expression of 10.074±0.1480, and the parental D2 strain had expression of 9.250±0.153, meaning that there is a statistically significant difference in expression between the two parental strains (p<0.05). Heritability of Stc1 was calculated to be 0.31, which indicates that Stc1 is affected by genetic variability, and allows us to identify a locus for Stc1.
Simple interval mapping for Stc1 found a statistically significant eQTL with an LRS of 19.6 (genome-wide p<0.05) on chromosome 6 at 48.92 Mb where the Gpnmb gene is located (Figure 6). The Stc1 gene itself is located on chromosome 14 at 69.04 Mb, which means that Stc1 has a trans-eQTL. Composite interval mapping showed no secondary locus that modulates Stc1 expression, suggesting that part of the expression variance in Stc1 is regulated by the DNA variants around the Gpnmb locus.
The top 500 probe sets were submitted to Webgestalt for gene function enrichment analysis. Of these probe sets, 26 could not be mapped to any Entrez gene ID or mapped to multiple gene IDs, resulting in 474 probe sets that could unambiguously map to 368 unique gene IDs. These genes were used for GO and pathway analysis.
The most statistically significant enrichment category for biologic processes was “cellular process” (n=321, adjP=2.26e-29) that includes “cell death” (n=66, adjP=1.84e-10), “programmed cell death” (n=64, adjP=1.16e-10), and “apoptotic process” (n=64, adjP=7.23e-11). The other most important and statistically significant enrichment categories included “mitochondrion” (n=67, adjP=2.01e-10) and “peroxisome” (n=8, adjP=1.86e-02) for cellular component; and “oxidoreductase activity” (n=30, adjP=0.0001), “ATP binding” (n=58, adjP=9.86e-08), “kinase activity” (n=44, adjP=1.68e-09), and “MAP kinase phosphatase activity” (n=5, adjP=3.15e-05) for molecular function. The GO results of all statistically significant enriched categories are listed in Appendix 1. Gene pathway analysis resulted in 21 statistically significant pathways (false discovery rate [FDR] p<0.05, Table 2), most of which are related to the electron transport chain, the mitochondrion, apoptosis, oxidative stress, the metabolism, the immune response, and inflammation.
To generate a Gpnmb and Stc1 coexpression network, and investigate a biologic relationship through which these two genes interact in the retina, we focused on the gene group that contained Gpnmb and Stc1 within the statistically significant enrichment. We found that Gpnmb and Stc1, as well as 27 other genes, are regulated by activator protein 1 (AP-1). We then uploaded them to GCAT for functional coherence analysis and gene network construction. These genes showed statistically significant functional cohesion with a literature p value of 3.267436e-16 (Figure 7). The network graph indicated that the expression levels of Gpnmb and Stc1 are directly linked to each other. Multiple resources, such as Chilibot, GeneCard, and PubMed, were used to study the function of each member in the Gpnmb and Stc1 coexpression network. In addition to Gpnmb and Stc1, some of genes in this network (Pex5 [Gene ID: 19305, OMIM: 600414], Rgs2 [Gene ID: 19735, OMIM: 600861], and Cd68 [Gene ID: 12514, OMIM: 153634]) are already known to be related to glaucoma. Most of the other genes are highly connected to apoptosis, oxidative stress, and mitochondria.
The results showed that Gpnmb expression is statistically significantly correlated with the optic nerve cross-sectional area, photoreceptor density, and iris pigmentation. Stc1 expression is statistically significantly correlated with optic nerve density and iris pigmentation. Among these correlated phenotypes, the iris pigmentation phenotype (graded from 0 (lowest) to 4 (highest) of iris transillumination) had a highly statistically significant negative correlation with the expression of Gpnmb and Stc1 (p<0.0005, Figure 8), indicating that Gpnmb and Stc1 play the same role in maintaining iris pigmentation.
Systems genetics is an approach to understanding the flow of biologic information that underlies complex traits [31]. This approach involves analysis of sets of causal interactions among DNA variants (such as SNPs), classic traits (such as IOP), and intermediate phenotypes (such as transcripts) in populations that vary for traits of interest. The BXD recombinant inbred (RI) strains are the best mouse genetic reference population for a systems genetics study of glaucoma, because the parental strains, the D2 mouse and the B6 mouse, differ in the presentation of glaucoma. Their progeny have different phenotypes and gene expression, and provide a unique and valuable opportunity to define novel genes and molecular networks involved in glaucoma. In this study, we used the systems genetics method and BXD RI strains to analyze genes that are highly correlated with and lie downstream of Gpnmb, a gene involved with glaucoma. We found Stc1 to be highly correlated with Gpnmb, involved in the glaucoma disease process, and the best downstream candidate gene of Gpnmb. We present a genetic network for Stc1 that involves genes related to optic neuropathy, gene enrichment categories Stc1 plays a role in, and pathways through which Stc1 interacts with other genes to affect glaucoma phenotypes.
Stanniocalcin-1, encoded by Stc1, is a glycoprotein found to have antioxidant and antiapoptotic properties. Stc1 inhibits the inflammatory cascade, induces antioxidant and antiapoptotic mechanisms, and reduces reactive oxygen species [32-35]. In glaucoma, these factors have been shown to reduce IOP and provide neuroprotection in human eyes [36,37], as well as prevent the ensuing apoptosis in retinal photoreceptors [38,39]. Stc1 delays RGC apoptosis in the rat model of intraorbital optic nerve transection and decreases oxidative stress by reducing the reactive oxygen species (ROS) in the eye [40].
After gene function enrichment analysis, one of the most statistically significant biologic processes found to be associated with Stc1 was “apoptotic process,” as well as its two related pathways “chemokine signaling pathway” and “apoptosis.” Apoptosis, the programmed death of a cell, is an important process in the development of the nervous system. Apoptosis also has been implicated in various neurodegenerative diseases, such as glaucoma [41]. RGCs have been shown to die via apoptosis in experimental glaucoma, and after the axon of the nerve is severed, as defined by light microscopy [42]. Further investigations have shown that elevated IOP in the extracellular matrix of the retina is strongly correlated with RGC apoptosis in glaucoma [43]. The expression of apoptosis and inflammation-related genes, such as Il18 (Gene ID: 16173, OMIM: 600953), NF-κB (Gene ID: 18033, OMIM: 164011), Mapk1 (Gene ID: 26413, OMIM: 176948), Mmp-2 (Gene ID: 17390, OMIM: 120360), Timp-1 (Gene ID: 21857, OMIM: 305370), and apoptotic signaling components, has been found to be elevated in the DBA/2J glaucoma mouse model [44]. The present GO and pathway results support the role of Stc1 in apoptosis, and most likely in glaucoma.
One of the most statistically significant cellular components associated with Stc1 was “mitochondrion.” Neurons, which have a high-energy requirement, are dependent on mitochondria for not only their source of energy but also for calcium signaling and apoptosis. IOP elevation has also been linked to mitochondrial damage in the optic nerve head through mitochondrial fission and cristae depletion, although currently it is unknown whether IOP elevation leads to mitochondrial alterations [45]. When there is a mitochondrial malfunction, free radicals are produced in excess, which can lead to oxidative stress, common in glaucomatous tissues. This oxidative damage can damage cellular macromolecules, resulting in neuronal degeneration and RGC death [46,47]. Oxidative DNA damage has also been found to be statistically significantly increased in patients with glaucoma [48]. Further support for mitochondrial and oxidative stress involvement in Stc1 was shown in the present pathway analysis, where we found the “electron transport chain,” “mitochondrial gene expression,” and “oxidative stress” are statistically significant pathways.
Transcription factor analysis found AP-1 regulates Gpnmb and Stc1 gene expression. To study the interaction of Gpnmb and Stc1, we constructed a Gpnmb and Stc1 coexpression network based on the group of genes targeted by AP-1. The AP-1 family of transcription factors consists of homo- or heterodimers encoded by the fos and jun gene families. AP-1 regulates cell growth and differentiation, and is affected by the oxidative state of the cell [49]. AP-1 has been suggested to play a role in apoptosis, although the exact mechanism is not clear [50-52]. Considering the relationship AP-1 has with Gpnmb and Stc1 and with oxidative stress and apoptosis, AP-1 may affect the pathogenesis of glaucoma through the regulation of Gpnmb and other genes that interact with Gpnmb in this gene network.
In addition to Gpnmb, we found three genes (Pex5, Rgs2, and Cd68) in the gene network are highly correlated with glaucoma based on previous literature reports. Peroxisome biogenesis factor 5, encoded by Pex5, plays an essential role in peroxisomal protein import. The Pex5-knockout mouse model has peroxisomal metabolism defects, leading to pleomorphic mitochondria and a marked increase in ROS [53]. Zellweger syndrome spectrum is a syndrome in which there are no functional peroxisomes due to deletions or mutations in PEX (Gene ID: 5830, OMIM: 600414) genes [54]. Patients with Zellweger syndrome may exhibit characteristics such as craniofacial dysmorphia, neonatal seizures, and ocular anomalies, like cataracts, pigmentary retinopathy, and glaucoma [55,56].
Regulator of G protein signaling 2, encoded by Rgs2, acts as a GTPase activating protein. Rgs2 is widely expressed in the rat retina, particularly in the ganglion cell layer [57]. Mice with null Rgs2 expression have also been associated with lower IOP and enhanced RGC survival [58]. This is possibly because Rgs2 is a negative regulator of Gαq, which affects the constriction of smooth muscle [59-61]. Mice with null Rgs2 also had increased width in the Schlemm’s canal, which is close to a possible outflow pathway for aqueous humor in patients with open-angle glaucoma [62]. However, the exact mechanism through which Rgs2 affects IOP is not clear.
Cd68 encodes cluster of differentiation 68, a protein highly expressed in macrophages, and has been used as a marker for phagocytic amoeboid microglia in axonal injury [63,64]. Normally, Cd68 is not present in the optic nerve head. However, in glaucoma, Cd68 has been found to be present in microglia, indicating that there is phagocytic activity at the ONH and inflammatory activity in the stroma of the iris and ciliary body [65,66].
Several of the genes in the present gene network are also highly correlated with apoptosis, oxidative stress, and mitochondria based on previous literature. Il24 (Gene ID: 93672, OMIM: 604136) and Mt3 (Gene ID: 17751, OMIM: 139255) are associated with oxidative stress, Mpv17 (Gene ID: 17527, OMIM: 137960) is a mitochondrial protein, and Eef1a2 (Gene ID: 13628, OMIM: 602959) has antiapoptotic effects.
Interleukin 24 (IL-24), encoded by Il24, is a proinflammatory cytokine in the IL-20 subfamily of interleukins [67]. Il24 has proapoptotic characteristics in cancer cells [68,69]. Il24 expression has been found to be expressed in the retina of DBA/2J mice [70].Il24 also induces the production of Il6, which is upregulated in the trabecular meshwork by oxidative stress, elevated IOP, and optic nerve head injury [71-73]. Mt3 encodes metallothionein 3, which is a protein that regulates zinc levels in the central nervous system in response to the oxidative status of the cell [74,75]. Under oxidative stress, zinc levels in the cytosol and lysosomes rise, resulting in increased cell death [76]. Mpv17 codes for a mitochondrial inner membrane protein, and its absence or malfunction has been found to cause oxidative phosphorylation depletion [77]. Mpv17 has also been implicated in mitochondrial DNA (mtDNA) depletion syndromes and oxidative phosphorylation activity in yeast [78,79]. Eukaryotic translation elongation factor 1 alpha 2 (Eef1a2) has been shown to have antiapoptotic effects [80-82].
Because these genes have been associated with either glaucoma or one of the many factors that affect the glaucoma phenotype, they may be associated with the glaucoma phenotype. However, further studies are needed to clarify the relationship between Gpnmb, Stc1, and these genes.
In summary, the present systems genetics analysis suggested that Stc1 could be a downstream candidate gene for Gpnmb. The gene function enrichment analysis of Gpnmb and Stc1 covariant genes suggests that Gpnmb may interact with Stc1 and other genes in the network to develop glaucoma pathogenesis through mechanisms of apoptosis and oxidative stress. However, further investigation is needed to elucidate the relationship among Gpnmb, Stc1, and the genes of the proposed genetic network.
Appendix 1. The gene ontology results of all significant enriched categories.
This work was supported by NIH Grant R01EY021200 (LL, PI) and Research to Prevent Blindness Stein Innovation Award (Robert Williams, PI). This study was supported in part by the by the Natural Science Foundation of Jiangsu China Grant BK2008186 and BK20161287 (HL), and Science Foundation of Nantong China Grant HS2011005 and HS2014005 (HL). We thank Dr. Robert Williams for his financial and bioinformatics support for this study. We thank Dr. David Ashbrook for his editing manuscript. Dr. Lu Lu (lulu@uthsc.edu) and Dr. Jianping Tong (tongjp2000@hotmail.com) are co-corresponding authors for this paper.