Differential Expression Profile Prioritization of Positional Candidate Glaucoma Genes: The GLC1C Locus

Objectives  To develop and apply a model for prioritization of candidate glaucoma genes.

Methods  This Affymetrix GeneChip (Affymetrix, Santa Clara, Calif) study of gene expression in primary culture human trabecular meshwork cells uses a positional differential expression profile model for prioritization of candidate genes within the GLC1C genetic inclusion interval.

Results  Sixteen genes were expressed under all conditions within the GLC1C interval. TMEM22 was the only gene within the interval with differential expression in the same direction under both conditions tested. Two genes, ATP1B3 and COPB2, are of interest in the context of a protein-misfolding model for candidate selection. SLC25A36, PCCB, and FNDC6 are of lesser interest because of moderate expression and changes in expression. Transcription factor ZBTB38 emerges as an interesting candidate gene because of the overall expression level, differential expression, and function.

Conclusions  Only 1 gene in the GLC1C interval fits our model for differential expression under multiple glaucoma risk conditions. The use of multiple prioritization models resulted in filtering 7 candidate genes of higher interest out of the 41 known genes in the region.

Clinical Relevance  This study identified a small subset of genes that are most likely to harbor mutations that cause glaucoma linked to GLC1C.

High prevalence and potential for severe outcome combine to make adult-onset primary open-angle glaucoma (POAG) a significant public health problem.¹ Genetic components to glaucoma are suggested by high concordance for POAG among older monozygotic twins² as well as mapping of a large number of glaucoma risk factor loci^3,4 plus 12 POAG loci, called GLC1 loci.^5-14 No unifying cellular or biochemical themes have emerged from the identification of the first 3 GLC1 genes: myocilin (MYOC; Online Mendelian Inheritance in Man [OMIM] 601652) at the GLC1A locus,^15,16 optineurin (OPTN; OMIM 602432) at the GLC1E locus,¹⁷ and WD-repeat36 (WDR36) at GLC1G¹¹; these genes do not account for most of the familial cases of POAG.

The GLC1C locus is of special interest because the gene itself has not been found but the locus has been confirmed through identification of additional GLC1C families^7,18-20 who provide optimal samples for screening candidate genes for mutations.^7,18,20 The existence of 2 distinct GLC1C haplotypes suggests that mutations will not be limited to rare descendants of a single founder, so GLC1C mutations may be found in multiple current glaucoma genetics study populations.

Here we describe expression data for 41 genes within the GLC1C genetic inclusion interval and 7 high-priority candidate genes selected based on overall or differential expression. These data should help to inform the work of many groups in the field who are looking for genes that cause glaucoma.

For this institutional review board–approved study, eyes were obtained from the Midwest Eye Banks (Ann Arbor, Mich), which obtained informed consent and confirmed that no donors had been diagnosed with glaucoma. Fifth-passage primary cultures of human trabecular meshwork (TM) cells used in age experiments came from 12- and 60-year-old females and 16- and 74-year-old males. Cells were grown to confluence and maintained for 1 week before isolating RNA as previously described.²¹ In dexamethasone studies, fifth-passage TM cells were derived from young donors aged 12 and 16 years (as described in this article earlier) and a 17-year-old girl. The TM cells were grown with and without a 21-day course of 100nM dexamethasone as previously described.²¹

Labeling, hybridization to Affymetrix U133A GeneChips (Affymetrix, Santa Clara, Calif), and data extraction were done as previously described.²¹ The RNA was prepared from 2 separate flasks for each culture, and conditions and data from biological replicates were compared as previously described.²¹

Image analysis was performed using Affymetrix Microarray Analysis Suite version 5.1 with samples scaled to 1500. Signals below 300 were considered absent and signals between 300 and 750 were considered marginal. Values above 750 were considered present.

Fold change for each probe was calculated by dividing the mean “treated” signal by the mean “untreated” signal in the dexamethasone experiments and the mean “old” signal by the mean “young” signal in the age experiments. The inverse function was applied to values below 1 to indicate decreased expression. The significance of fold change was evaluated using the t test, with P<.005 considered statistically significant. Scatterplots were drawn using Spotfire DecisionSite 8.2 software (Spotfire, Inc, Cambridge, Mass).

Of the genes in the GLC1C interval between D3S3637 and D3S3694 obtained from the University of California, Santa Cruz, Genome Browser (http://genome.ucsc.edu),¹⁸ 28 were assayed by the U133A GeneChip and 13 were assayed by quantitative polymerase chain reaction (qPCR) using primers (Table 1) in iQ SYBR Green Supermix reactions (BioRad Laboratories, Hercules, Calif) as previously described.²¹ We did not perform qPCR confirmation of all of the GeneChip values because results from previous experiments show GeneChip data to be adequate in the signal range being evaluated here.

Data from qPCR were assigned an estimated Affymetrix signal value (EASV) using linear regression. The log₂ values of Affymetrix signals from 9 control genes (Table 1) were plotted against their corresponding qPCR threshold values normalized to glyceraldehyde-3-phosphate dehydrogenase (Figure 1). These data are from the same RNA and GeneChips used for the other aspects of this study. The plot had an R² value of 0.7459, and imputed values were calculated using the equation EASV = 2^[(C_t ^{× −0.3866) + 21.517]}, where C_t is the threshold value for qPCR. The validity of the regression was examined by comparing the actual and imputed signal values using qPCR data from 3 genes that were not used to construct the curve (Table 2). When comparing data for the 3 control genes under 4 conditions, we found that 8 of 12 data points showed less than a 2.5-fold difference between the imputed and actual Affymetrix values. In a separate experiment in which data from a different RNA sample tested on the Affymetrix U133Av2 version 2 chips were examined for the 6 experimental genes for which we need to impute values, we found that 5 of 6 measured data points fell within 2.5-fold of the imputed values.

The EASV gives us an approximate value for comparing data derived from GeneChips with qPCR data for genes not on the GeneChip, but only approximate and qualitative conclusions can be drawn. We have not compared variability across protocols, and we recommend that linear regression be performed using data in which the known and unknown Affymetrix and qPCR assays are carried out by the same person on the same instruments. Use of a larger number of replicate values is expected to improve the accuracy of the EASV. We consider the EASV in this case to be approximately accurate and not valid for making subtle distinctions between similar signal levels. Future development of this method, including compensation for factors such as PCR product size and guanine-cytosine content, may allow for increased precision but will not allow for extension of this approach to genes expressed at levels above the saturation point of the Affymetrix GeneChip. Imputed values are used for display purposes, but primary conclusions are drawn from qPCR where available.

The 7 879 854–base pair GLC1C interval between D3S3637 and D3S3694 contains 41 genes (Table 3). We evaluated gene expression under 4 conditions: TM cells cultured from 2 young donors, TM cells cultured from 2 older donors, and 3 young TM cell cultures treated with and without dexamethasone. Although some genes can show substantial individual variation in responses to dexamethasone between individuals, the largest difference in fold change between individuals for any gene in this region was only 0.46.

Adenosine triphosphatase, sodium-potassium–transporting, β3 polypeptide (ATP1B3) has the highest mean signal level across all 4 conditions tested (range, 16 000-26 000), which is 6.6- to 10.8-fold higher than the average signal of 2400 for all probes and all conditions on 27 GeneChips. Coatomer protein complex, subunit β2 (β prime) (COPB2) has the second highest signal (range, 6000-11 000). Two other genes, transmembrane protein 22 (TMEM22) and solute carrier family 25, member 36 (SLC25A36), reach signal levels above 5000 for at least 1 of the 4 conditions, and 15 other genes show signals above 1000 for at least 1 condition (Table 4). Signals for 16 of the 41 genes in the GLC1C region are considered present under all 4 conditions (Figure 2).

Eleven genes with Affymetrix signals below 750 in all 4 conditions were considered marginally expressed. Signals for 3 of these genes were considered absent in 1 or more conditions (Table 4). Probes corresponding to 4 genes, forkhead box L2 (FOXL2), armadillo repeat-containing 8 (ARMC8), sex-determining region Y–box 14 (SOX14), and calsyntenin 2 (CLSTN2), showed absent signals in all 4 experimental conditions (Table 4 and Figure 2).

For 13 genes not on the U133A GeneChip, surrogate Affymetrix microarray signals imputed from the qPCR cycle threshold (see the “Quantitative PCR” subsection of the “Methods” section) allowed for a qualitative comparison of signal levels for all GLC1C-interval genes on the same plot (Table 4 and Figure 2). Three genes encoding acid phosphatase–like 2 (ACPL2), hypothetical protein MGC40579 (MGC40579), and 5′-to-3′ exoribonuclease 1 (XRN1) were considered present in all 4 conditions but did not show large fold change values associated with age or dexamethasone treatment (Table 5). Signals for the remaining genes missing from the U133A GeneChip were considered marginal or absent because of failure to produce either a qPCR or GeneChip signal for all of the conditions (Table 4). Data on these same genes from a very limited number of Affymetrix U133 Plus 2.0 chips confirm absent expression (data not shown).

The 1 significant change large enough to seem biologically interesting is the 5.3-fold decrease with age for TMEM22 (Figure 2A and Table 4). Four genes responded to dexamethasone exposure with statistically significant decreases (propionyl coenzyme A carboxylase, β polypeptide [PCCB] and TMEM22) and increases (SLC25A36 and zinc finger and BTB domain containing 38 [ZBTB38]) (Figure 2B and Table 4). A decrease in TMEM22 was the only statistically significant expression differential in both experiments (Figure 2 and Figure 3).

In silico data show that most GLC1C-region genes expressed in the TM are also expressed in many tissues outside the eye.²²TMEM22 (expressed in the brain, heart, pancreas, kidney, and lung) and ZBTB38 (negligible expression in several tissues that were assayed) show a more restricted range of expression than most of the GLC1C-region genes (http://www.genecards.org/cgi-bin/carddisp.pl?gene=TMEM22&search=tmem22 and http://www.genecards.org/cgi-bin/carddisp.pl?gene=ZBTB38). TMEM22, COPB2, and ATP1B3 are expressed in the retina at more than 5 times the TM level whereas mitochondrial ribosomal protein S22 (MRPS22), phosphoinositide-3-kinase, catalytic, β polypeptide (PIK3CB), and Fas apoptotic inhibitory molecule (FAIM) all show lower levels of expression in the TM than in the retina. The remaining genes showed comparable expression for the retina and the TM (data not shown).

Because elevated intraocular pressure is characteristic of GLC1C glaucoma, prioritization of genes based on functional information from primary cultured TM cells is potentially relevant. The original GLC1C family had a mean maximum known intraocular pressure of 23.3 mm Hg (range, 14-36 mm Hg), and elevated intraocular pressure has been observed in subsequent GLC1C families as well.^7,18,20 Thus, while a GLC1C mutation might also affect retinal ganglion cells or the optic disc, the role of the TM in regulation of intraocular pressure makes it relevant to the study of the GLC1C gene.

MYOC offers us an intriguing example of a type of information that can be used in prioritizing candidate genes: differential expression. MYOC can initiate familial POAG when mutated¹⁵ or can play a role through changes in expression in response to the environment (dexamethasone exposure) or development (aging). Our positional differential expression profile model calls for mutation screening in positional candidate genes showing altered gene expression that potentially mimics changes in the gene product or activity level that can result from mutations.

TMEM22 is the only gene in the GLC1C interval that fits the positional differential expression profile model (Table 6). TMEM22 has a statistically significant decrease in response to both aging and dexamethasone exposure (5.3- and 1.6-fold, respectively) (Table 4 and Figure 3) and the third-highest expression in young (signal level = 4 535) and untreated (signal level = 6 458) donors. Expression of TMEM22 outside of the eye is much more restricted than expression of most of the genes in this interval (GeneCards, http://www.genecards.org); further, our microarray data (F.W.R., unpublished data, June 2006) show absent expression in the retina. Limited information suggests that the TMEM22 gene product might be involved in transport functions (Table 6). With a pair of DUF6 domains within a region of homology to the bacterial l-rhamnose symport transporter protein encoded by the RhaT gene,²³TMEM22 may be a member of the drug metabolite transporter superfamily of permeases, some members of which are involved in the transport and metabolism of carbohydrates, amino acids, toxins, or drugs.^24-26 If TMEM22 were the GLC1C gene, functional studies of its transport functions could point to key metabolites and pathways for further study.

Two genes, ATP1B3 and COPB2, show signal levels in excess of 10 000, a signal level seen for only 4.7% of the probes on the U133A GeneChip. ATP1B3 shows the highest expression level (signal range, 16 217-26 673) (Table 4) but no significant differential expression. This signal level can be compared with the higher average MYOC signal level of 78 368 in response to dexamethasone exposure. There is evidence of multiple different β-chain isoforms and substantial alternative splicing of this isoform (http://www.genecards.org/cgi-bin/carddisp.pl?gene=ATP1B3). ATP1B3 encodes the β3 subunit of the integral membrane sodium-potassium adenosine triphosphatase responsible for transport of sodium and potassium ions across the plasma membrane (Table 6).²⁷ The ability of several transporters to alter the effects of ion transport on cell volume fluid flow in glaucomatous eyes has recently been ruled out.^28-30 It remains to be seen whether ATP1B3 is playing a role in the regulation of cellular or extracellular volume.

COPB2, the second most highly expressed gene in the GLC1C interval (signal range, 6414-10 693), shows no significant differential expression (Table 4). COPB2 is widely expressed throughout the body, and there appear to be at least 5 different splice variants (http://www.genecards.org/cgi-bin/carddisp.pl?gene=COPB2&search=COPB2copb2&suff=txt). COPB2 encodes coatomer complex subunit β2, 1 of the 7 proteins that make up the approximately 700-kd Golgi coatomer complex that coats non–clathrin-coated vesicles³¹ (Table 6). The COPB1 coatomer complex selectively recruits proteins from the Golgi apparatus into budding vesicles that are transported back to the endoplasmic reticulum, a process that might be expected to affect a variety of proteins important to glaucoma, perhaps even myocilin.

Based on sequence homology, SLC25A36 encodes member 36 of solute carrier family 25 (Table 6) and shows moderate expression (signal range, 3814-5312) and significant expression change in response to dexamethasone exposure (P = .003) but not aging (P = .30) (Table 4). It is expressed in many tissues, and 5 different splicing variants have been observed (http://www.genecards.org/cgi-bin/carddisp.pl?gene=SLC25A36&search=SLC25A36). Other members of this family of proteins are involved in the transport of small molecules such as citrate, adenine nucleotide, ornithine, or malate.^32-35

Although the fold change of PCCB in response to dexamethasone exposure is statistically significant, the magnitude of the change is modest (Table 4, Figure 2, and Figure 3). Mutations in this gene have caused propionic acidemia,³⁶ a trait not known to occur in members of the known GLC1C families. This does not rule out PCCB as a candidate, as mild cases might go undiagnosed and some genes can cause drastically different phenotypes involving completely different tissue types and clinical features.^37-39

We rank genes with substantial signal levels ahead of the genes with absent expression or qPCR product (Figure 2 and Table 6). Changes observed for most of these genes are too small to be significant indicators of biologically important changes and are more likely an indicator of experimental and interindividual variation. The GLC1C gene could potentially be expressed in vivo but not in this cultured cell system. We cannot conclude that expression levels in these cells from adolescents and adults indicate levels that occur during other stages of life, including prenatal development. Thus, the nonexpressed genes are assigned low priority but cannot be removed from the list based on our experimental system.

One of the functional themes in a number of developmental glaucomas is that of transcription factors such as paired-like homeodomain transcription factor 2 (PITX2), forkhead box C1 (FOXC1), and LIM homeobox transcription factor 1β (LMX1B),^40-43 whose phenotypic range includes POAG in some individuals.^43,44 It is interesting to note 3 transcription factors in the GLC1C region. Two of these, FOXL2 and SOX14, are not expressed at all in our cell culture system. Mutations in FOXL2 are already associated with blepharophimosis, ptosis, and epicanthus inversus,⁴⁵ which do not occur in the GLC1C families, although this does not completely rule out involvement of this gene in more than 1 phenotype.

ZBTB38^46,47 is expressed at moderate levels (signal range, 2151-4543) and shows a 1.65-fold increase in response to dexamethasone exposure (P = .004) but only a small, nonsignificant decrease with age (1.19-fold decrease; P = .05) (Table 4). Given that haploinsufficiency or duplication of transcription factor genes can cause disease,^43,45,48 a change of 1.65-fold is potentially biologically important. The ZBTB38 gene product is a Kaiso-like transcriptional repressor containing a BTB (POZ) domain and multiple C2H2 zinc finger domains^46,49 (http://ca.expasy.org/uniprot/Q8NAP3 and http://www.genecards.org/cgi-bin/carddisp.pl?gene=ZBTB38&search=ZBTB38). ZBTB38 is the human homolog of the rat Zenon gene that is hypothesized to play a role in survival or phenotypic maintenance of postmitotic neurons,⁴⁷ suggesting a possible role in the survival of retinal ganglion cells. ZENON has also been shown to regulate expression of tyrosine hydroxylase (TH), which participates in a pathway leading to levodopa production,⁵⁰ especially in neurons.⁵¹ Induction of ZBTB38 does not appear to lead to suppression of TH in our TM cultures, as TH is uniformly absent (signal range, 16-108). Feeding levodopa to a mouse model of glaucoma during development ameliorates some aspects of the glaucoma phenotype,⁵⁰ so we cannot rule out ZBTB38 effects on the formation of structures of the angle and outflow pathway, different expression patterns of ZBTB38 and TH in native tissue, or a role for ZBTB38 in retinal ganglion cells and the optic nerve of the adult. Thus, the connection between ZBTB38 and the animal model of Libby et al⁵⁰ is intriguing, but dexamethasone-induced changes in ZBTB38 expression do not appear to invoke a direct response in TH expression in cultured adult TM cells. Even if ZBTB38 regulates TH in neurons, its effects in the TM might operate through other downstream regulatory targets. If ZBTB38 were the GLC1C gene, there would be questions about how to relate the beneficial effects of feeding levodopa in mice⁵⁰ to the development of therapeutic strategies for GLC1C families with glaucoma.

TMEM22 fits our positional differential expression profile model whereas other models recommend ATP1B3, COPB2, SLC25A36, PCCB, and transcription factor ZBTB38. TMEM22 and other candidates are being screened. Some glaucoma mutations may have function-specific effects (Table 6) rather than acting through altered expression, as previously seen for night blindness or retinal degeneration resulting from certain rhodopsin mutations.⁵² We do not expect the positional differential expression profile model to pinpoint every glaucoma gene; information on the expression of GLC1C-region genes described here should be of interest to those engaged in the search for glaucoma genes, even if some investigators use a different basis for prioritization of GLC1C candidate genes.

Correspondence: Julia E. Richards, PhD, Department of Ophthalmology and Visual Sciences, University of Michigan, 1000 Wall St, Ann Arbor, MI 48105 (richj@umich.edu).

Submitted for Publication: June 30, 2006; final revision received September 18, 2006; accepted September 19, 2006.

Financial Disclosure: None reported.

Funding/Support: This work was supported by grants EY07003 (core grant, Kellogg Eye Center), EY010572 (core grant, Casey Eye Institute), EY011650 (Dr Wirtz), and EY09580 (Dr Richards) from the National Institutes of Health, Bethesda, Md, unrestricted grants to the Kellogg Eye Center and the Casey Eye Institute from Research to Prevent Blindness, Inc, New York, NY, and the Van Arnam Glaucoma Research Fund of the Department of Ophthalmology and Visual Sciences, University of Michigan, Ann Arbor.