Gene Expression and Pathways Underlying Form Deprivation Myopia in the Guinea Pig Sclera

Abstract

Purpose: Posterior scleral remodeling accompanies myopia. In guinea pigs developing myopia, the region around the optic nerve (peripapillary zone, PPZ) rapidly expands followed by inhibition in eye size in the periphery. We studied the differential gene expression in the sclera that accompanies these changes.

Methods: Guinea pigs were form-deprived (FD) for 2 weeks to induce myopia, while the fellow eye served as a control. After 2 weeks, the PPZ and the peripheral temporal sclera were isolated in representative animals to extract the RNA. RNA sequencing was undertaken using an Illumina HiSeq 2000, with differential expression analyzed using Voom and pathways analyzed using the Ingenuity Pathway Analysis tool. RNA from additional PPZ and peripheral temporal sclera in FD and fellow eyes was used for validation of gene expression using quantitative real-time PCR (qRT-PCR).

Results: In myopic sclera, 348 genes were differentially expressed between PPZ and the peripheral temporal region (corrected P < 0.05), of which 61 were differentially expressed in the PPZ between myopic and control eyes. Pathway analyses of these gene sets showed the involvement of Gαi signaling along with previously reported gamma-aminobutyric acid (GABA) and glutamate receptors among numerous novel pathways. The expression pattern of three novel genes and two myopia-related genes was validated using qRT-PCR.

Conclusions: Gene expression changes are associated with the rapid elongation that occurs around the optic nerve region during the development of myopia. A prominent change in Gαi signaling, which affects cAMP synthesis and thus collagen levels, may be critical in mediating the regional changes in myopic sclera.

The sclera is a fibrous connective tissue that constitutes the outer covering of the eye. Turnover of components in the scleral matrix is essential for maintaining scleral structure and thus ocular shape. Structural remodeling in the posterior sclera accompanies the abnormal increase in eye size that characterizes the development of myopia.^1,2 The sclera in myopic human eyes is thinner, and loss of matrix components has been observed in myopic cadaveric eyes,^3,4 Induction of myopia in animal models also leads to a thinner sclera^1,2 (McFadden S, et al. IOVS 2008;49:ARVO E-Abstract 3713) and is associated with a decrease in the synthesis of collagen, proteoglycan, and other scleral matrix components as seen in mammalian sclera from tree shrews,⁵ guinea pigs,⁶ marmosets,⁷ and in the fibrous scleral layer of chicks.⁸ This reduction in the scleral matrix has been proposed to result in a thinner, weaker, and more extensible tissue, which may facilitate excessive ocular elongation during myopia.^1,2 

Various candidate genes that might underlie these scleral changes in myopia have been investigated primarily using quantitative real-time polymerase chain reaction (qRT-PCR). These include collagen subtypes (e.g., collagen I, III, V) (Gentle A, et al. IOVS 2002;43:ARVO E-Abstract 2449), proteoglycans (e.g., aggrecan, decorin), growth factors (e.g., transforming growth factor-beta, fibroblast growth factor), and matrix metalloproteinase (e.g., subtypes 1, 2, 3).¹ More recently the expression of genes involved in matricellular signaling (e.g., CTGF, FBLN1), transcription regulators (e.g., RARB, RXRB), and cell surface receptors (e.g., FGFR, TGFBR)^9,10 has also been identified using qRT-PCR. 

These animal studies compared myopic eyes (or eyes recovering from myopia) to fellow or untreated eyes, using either whole sclera or scleral punches taken from the posterior pole. However, it is also likely that specific regional differences occur in different parts of the sclera.^8,11,12 For example, human myopic eyes tend to be myopic on axis, but show relative hyperopia in the periphery,¹³ suggesting that the changes underlying ocular shape may be complex. In the guinea pig model of myopia, a complex sequence of changes occurs in eye shape. Using an accurate method to characterize eye shape from sectioned eyes,¹⁴ it was found that a relatively early response to myopia is enlargement in the zone surrounding the optic nerve,¹⁵ termed the peripapillary zone (PPZ), with little change in the adjacent periphery (Zeng G, et al. IOVS 2011;52:ARVO E-Abstract 3923). However, after 2 to 3 weeks of FD, while the PPZ is still relatively enlarged, the periphery inhibits its growth, reducing the distance from the lens center to the retina by up to 100 μm (Zeng G, et al. IOVS 2011;52:ARVO E-Abstract 3923). These regional differences in ocular growth are accompanied by a rapid increase in the perimeter length of the sclera at the PPZ (+5.3 μm/°) with reduced scleral perimeter lengths in the temporal periphery (−1.3 μm/°) compared to the non-FD eyes (Zeng G, et al. IOVS 2013;54:ARVO E-Abstract 5180). Since it is well known that in myopic eyes, scleral changes are generally greater at the posterior pole compared to the far periphery, there may also be specific changes in scleral gene expression occurring in the PPZ. 

Considerable knowledge regarding gene expression changes associated with myopic sclera has been gleaned mostly from qRT-PCR and some microarray studies, but limited genome-wide studies have been undertaken (Frost MR, et al. IOVS 2012;53:ARVO E-Abstract 3452). Typically, analysis of gene expression patterns using qRT-PCR is limited to one or several genes. On the other hand, microarrays can detect the expression of thousands of genes in a single experiment, but they are limited by the number of genes represented on each microarray chip, and the complexity of normalization methods makes comparison between transcription levels difficult.^16,17 RNA sequencing (RNAseq)^18,19 has the advantage that it allows identification of gene expression profiles across the whole genome with possible discovery of novel genes and concomitant pathways without being limited by a predetermined choice of genes.^20,21 

Using RNAseq, over 2000 genes were identified to be differentially expressed in the sclera of myopic tree shrews, of which 38 candidate genes were validated using qRT-PCR (R² = 0.93) (Frost MR, et al. IOVS 2012;53:ARVO E-Abstract 3452). Agreement between these two techniques was also evident in 11 out of 14 bidirectionally regulated genes in the choroid of tree shrews that were identified during the development of and recovery from myopia (He L, et al. IOVS 2013;54:ARVO E-Abstract 3675). While these studies used the whole scleral or choroidal tissues, the aim of the current study was to identify differences in the expression of genes unique to the PPZ region (relative to the temporal periphery that excludes the posterior pole) between myopic and nonmyopic eyes using whole-transcriptome analysis (RNAseq), with the ultimate aim of identifying functional pathways in the scleral region that are most affected during myopia development. 

Methods

Animals

Young guinea pigs (Cavia porcellus, n = 29) were bred at the University of Newcastle and each litter was housed with their mother in plastic containers (65 × 45 × 20 cm) covered with tops made of stainless open grill. Individual boxes were illuminated overhead with white light (2 × LED; Luxeon Star, Philips Lumileds Lighting Company, USA) evenly diffused through an opaque barrier made of polymethyl methacrylate (Opal Perspex; Lucite International, Mitchell Plastics, Melbourne, Australia) set 200 mm above the boxes (luminance was 400 lux at the center of each box) and set on a 12-hour day/12-hour night cycle. All procedures were approved by the University of Newcastle and were in accordance with the ARVO Statement for the Use of Animals in Ophthalmic and Vision Research. 

Induction of Myopia and Ocular Measures

Animals were form-deprived (FD) in one eye for 2 weeks (from 5 to 19 days of age), while the other eye served as an untreated control. Form deprivation was induced with a translucent diffuser as previously described,²² and diffusers were cleaned and checked daily. At the end of 2 weeks, one drop of 1% cyclopentolate (Cyclogyl; Alcon, Macquarie Park, Australia) was topically applied to each eye and after 1.25 to 1.5 hours, the refractive error was measured using an infrared autorefractor (AR-20Nidek; Gamagori, Aichi, Japan). Axial length was measured using high-frequency ultrasound in anesthetized animals as previously described.²² Immediately after optical measurements were completed, the diffusers were replaced for at least 24 hours prior to tissue dissection. 

Scleral Dissection and Samples Used

Three of the animals were selected for gene profiling from a larger pool of 19 FD animals based on being from separate litters, on developing the mean amount of relative myopia compared to the fellow eye, and on the quality of their scleral dissections and RNA. Additional FD guinea pigs (n = 8) were used in the validation of candidate genes selected from the RNAseq results. The selected animals were deeply anesthetized with 1.5% isoflurane in oxygen and euthanized with a lethal dose of pentobarbitone sodium (130 mg/kg, Lethabarb; Virbac, Milperra, Australia) into the heart. Eyes were rapidly enucleated on ice (within 1 minute of death). The orbital fat, vitreous, and other ocular components were cleaned off the sclera and punches were made around the optic nerve (Fig. 1) using a 3-mm trephine punch (Medshop, Preston, Australia), and the optic nerve was separated with a 1-mm punch and discarded. The remaining tissue was cut vertically in half. The temporal half of this peripheral tissue was used for comparison as it clearly excluded the central visual axis and posterior pole (Fig. 1). Therefore, the three animals provided 12 scleral samples for analysis: 3 PPZ FD and 3 PPZ fellow untreated; 3 peripheral temporal FD and 3 peripheral temporal from untreated fellow eyes. All tissues were completely immersed in RNAlater (Qiagen, Chadstone, Australia) and transferred to −80°C until RNA extraction. 

Figure 1

 

Posterior view of the guinea pig eye showing the two regions sampled: PPZ, peripapillary zone with optic nerve at center cut out; T, peripheral zone on temporal side (blue). Red dashed circle and asterisk show the position of the posterior pole. Upward red arrows denote the location of the long ciliary artery used for horizontal orientation. NR, nasal retina; TR, temporal retina. Image adapted with permission from Zeng G, Bowrey HE, Fang J, Qi Y, McFadden SA. The development of eye shape and the origin of lower field myopia in the guinea pig eye. Vision Res. 2013;76:77–88. Crown Copyright © 2012 Published by Elsevier Ltd.

RNA Extraction

Total RNA was isolated using RNeasy fibrous tissue mini kit (Qiagen, Chadstone, Australia) according to the manufacturer's protocol. Briefly, the frozen sclera was allowed to thaw at room temperature while still immersed in RNAlater, after which it was homogenized with a rotor-stator homogenizer (Qiagen, Chadstone, Australia) in RLT buffer for 1 minute at full speed. The homogenized samples were treated with proteinase K and incubated at 55°C for 10 minutes. The samples were then centrifuged at full speed for 3 minutes and treated with 100% ethanol. The samples were then transferred to an RNeasy spin column with a 2-mL collection tube, washed with buffer RW1, and treated with DNase I to remove any contamination from genomic DNA. Following a final wash with buffer RPE, the RNA was eluted with 50 μL RNase-free water and quantified using a nanodrop 2000 spectrophotometer (Thermo Scientific, Scoresby, Australia). All buffers used above were obtained from Qiagen (Chadstone, Victoria, Australia). The extracted RNA was stored at −80°C until used. 

RNA Sequencing

Samples were sequenced at the Australian Genome Research Facility (AGRF, Melbourne, Victoria, Australia) using an Illumina HiSeq 2000 platform (Illumina, Inc., San Diego, CA, USA) that provided 150 to 180 million (100 base pair) single-end reads per lane. The samples were prepared for sequencing according to the Illumina TruSeq RNA v2 protocol. Briefly, 1 μg total RNA containing poly-A tail was taken from each of the 12 samples and purified using oligo (dT) beads, following which they were fragmented into smaller pieces using divalent cations with heat. The mRNA fragments were then transcribed into first-strand cDNA with reverse transcriptase and random primers, which was followed by second-strand cDNA synthesis using DNA polymerase and RNase H. Primer sequences that were not annealed to the template RNA during cDNA synthesis (overhangs) were converted into blunt ends using end repair mix. The 3′ end of the cDNA fragments was then adenylated, to which complementary adapter sequences, specific for the Illumina platform, were added for subsequent amplification using the polymerase chain reaction (PCR) to create a cDNA library. 

Quantitative Real-Time PCR (qRT-PCR)

A separate set of FD animals was used to study the expression profiles of a few candidate genes from RNAseq using qRT-PCR, for validation purposes. Scleral RNA was extracted from the PPZ and peripheral temporal regions from the FD eyes and the PPZ of the corresponding fellow eyes (n = 3). Guinea pig–specific primers were designed using the custom TaqMan assay design tool. Total RNA (100 ng) from each scleral sample was converted to cDNA using a QuantiTect Reverse Transcription kit (Qiagen, Chadstone, Australia). Briefly, any residual genomic DNA contamination was removed from the template RNA using the gDNA wipeout buffer (7×), followed by incubation with RT primer mix, RT buffer, and the reverse transcriptase enzyme at 42°C for 15 minutes and heating up to 95°C for 3 minutes to complete the reverse transcription process. The resulting cDNA was diluted to 1:10 for the subsequent qRT-PCR. 

A multiplex PCR was set up with the target gene (FAM-labeled) and housekeeping gene, GAPDH (Cp03755742_g1; VIC-labeled, Primer-limited; Applied Biosystems, Thermo Fisher Scientific, Waltham, MA, USA), in a single reaction. Each reaction was performed in duplicate (StepOnePlus Real Time PCR, Applied Biosystems) with the following settings: 1 cycle of 95°C for 10 minutes followed by 40 cycles at 95°C for 15 seconds and 60°C for 1 minutes. The amplification data were analyzed for differential gene expression using the comparative Ct method.²³ 

Data Analysis

Sequencing output (in fastq format) was analyzed using an IBM supercomputer (Barcoo: IBM iDataplex x86 system; Melbourne Bioinformatics, University of Melbourne, Australia). The guinea pig reference genome (cavPor3.0) along with the gene annotation file (in Gene Transfer Format [.gtf]) were downloaded from ENSEMBL (http://asia.ensembl.org/Cavia_porcellus/Info/Index; in the public domain). Additionally, RNAseq data were analyzed using the Galaxy pipeline through the research cloud of the National eResearch Collaboration Tools and Resources (NeCTAR). Raw sequencing files were mapped to the guinea pig genome using the TopHat tool,²⁴ and the gene read counts across replicates were estimated using the HTSeq package.²⁵ Genes with an expression value, measured in counts per million (CPM), of more than 1 in at least two replicates were normalized using the trimmed mean of M-values (TMM) method.²⁶ A blocking factor was used to account for matched eyes for each guinea pig when fitting the linear model. Additionally, a pairwise comparison was carried out between (1) PPZ and peripheral temporal sclera of myopic eyes and (2) PPZ of myopic sclera and PPZ of fellow control eyes. The library size ranged between 7.2 and 8.9 million mapped transcripts. The generated read count file and the annotated gene list for guinea pigs (.gtf file), were used to identify the differential gene expression profile using Voom.²⁷ An empirical Bayes moderated t-statistic was used to determine the P value of differentially expressed genes, which was corrected for multiple testing using the Benjamini-Hochberg method. Paired t-tests were used when comparing either the PPZ of myopic and fellow control eyes or the PPZ and peripheral temporal region within the myopic eyes. Moderated t-tests with empirical Bayes model were then used to generate a list of differentially expressed genes. The mean logarithmic expression of each gene was estimated across the three replicates in the FD and fellow control eyes, and the final log fold change in expression was calculated by estimating the difference in the logarithmic mean of these two values (i.e., difference between mean of treated and mean of control groups) as presented in Supplementary Tables S1 and S2. Analyses using the above tools (and packages) were carried out using default parameters described in the Galaxy pipeline. A volcano plot was generated to pictorially represent the degree of fold change in differentially expressed genes from the large RNAseq dataset in the form of scatter plots. A heatmap was also created to present the pattern of gene expression across replicates and between two treatment paradigms in the form of color gradients. Both the volcano plot and heatmap were designed using R software.²⁸ Significantly differentiated genes were used in subsequent pathway analyses using the Ingenuity Pathway Analysis tool (IPA; Qiagen, Redwood City, CA, USA). 

Results

Induction of Myopia and Sample Quality

At the end of the FD period, the mean intraocular difference in the larger set of animals in refractive error was −3.98 ± 0.53 diopter (D) and in ocular length was 88 ± 12 μm (P < 0.0001 in both cases, mean ± SEM, n = 19). In the three representative animals selected for RNAseq, these same mean differences were −3.76 ± 0.6 D and 95 ± 13 μm, respectively. The mean difference in refractive error and ocular length in the qRT-PCR larger set of FD animals (n = 8) was −5.52 ± 1.17 D (P < 0.01) and 80 ± 16 μm (P < 0.001), respectively, and for the selected animals actually used was −3.3 ± 0.9 D and 91 ± 27 μm, respectively (mean ± SEM). The average weight of scleral tissues dissected from each guinea pig eye was approximately 13 to 18 mg for PPZ samples and 25 to 35 mg for the peripheral temporal samples. The quantity of total RNA obtained from these two areas was approximately 20 and 30 ng/μL, respectively. While the RNA integrity number (RIN) was over 7 for most samples, which was considered to be a quality measure for nondegraded RNA, some of the samples presented with values less than 5. This may be due to the lower concentration of RNA extracted, which limits the ability to obtain an accurate RIN number in samples with yields less than 50 ng/μL.^29,30 Thus, the purity (A260/A280 ratio) of RNA was also measured using nanodrop with estimated RNA purities between 1.8 and 2.0 for all samples. 

Differential Gene Expression Between the PPZ and the Peripheral Temporal Sclera

Over 16,000 genes in the annotated gene list from the guinea pig genome were mapped to the RNAseq data using HTSeq. The Voom tool identified many thousands of differentially expressed genes between the PPZ and peripheral temporal sclera within the FD (myopic) guinea pig eyes. The number of differentially expressed genes with a fold change (FC) of at least 1.25, 1.5, 1.75, and 2.0 was 6341, 3846, 2810, and 2233. A pictorial representation of the overall spread of the differentially expressed genes across FC and P values is shown as a volcano plot (Fig. 2). A heatmap of the top 1000 genes (based on P value) was generated using the normalized read counts obtained for each differentially expressed gene (Fig. 3A). Additionally, 348 differentially expressed genes were statistically significant (see Supplementary Table S1 for the list of 348 genes) after multiple correction (P < 0.05) with a FC of at least 1.51. Most of these genes (337/348) also presented with at least a 2-fold change difference between the PPZ and the peripheral temporal sclera. Of these, the top 20 statistically significant genes following correction are shown in Table 1. 

Figure 2

 

The x-axis shows the log fold change (LogFC) and y-axis shows the significance of the differentially expressed genes (denoted as numerous dots). Dots to the right of “0” (in the x-axis) denote upregulation and those to the left denote downregulation of genes. For the y-axis, data points higher than a −log₁₀ (P value) of 1.25 (corresponding to an uncorrected P value of 0.05) denote significance (approximately marked with a horizontal broken line).

Figure 3

 

(A, B) Heatmaps of the top 1000 genes between the PPZ and peripheral temporal sclera within myopic and fellow control eyes, respectively. (C) Heatmap of 61 significant genes (P < 0.05) differentially expressed between the PPZ of form-deprived and that of fellow control eyes. The x-axis shows the gene expression in the form of a color gradient. Genes showing a greater change in differential expression between two given conditions are presented with darker shades, and genes with a z-score of 0 were not differentially expressed. The y-axis presents a dendrogram, which groups together the genes that express similar fold changes. PPZ FD, sclera from the peripapillary zone of form-deprived (i.e., myopic) eyes; T FD, temporal peripheral sclera of FD eyes; PPZ C, sclera from the peripapillary zone of fellow control (nonmyopic) eyes; T C, temporal peripheral sclera of fellow control (nonmyopic) eyes.

Table 1

 

Top 20 Differentially Expressed Genes Between the PPZ and Peripheral Temporal Sclera of FD (Myopic) Guinea Pig Eyes

 

In contrast, no statistically significant changes after multiple correction were apparent at the gene level when comparing the PPZ and peripheral temporal region within the fellow (control) eyes. An inconsistent grouping of the color gradient was observed in the heatmap generated for the top 1000 genes (based on uncorrected P value) in these nonmyopic eyes (Fig. 3B). 

Differential Gene Expression in the Sclera Related to Myopia: FD Versus Fellow Eyes

When the above subset of 348 genes was compared for differential expression between the PPZ of myopic and PPZ of fellow control eyes, 61 genes (see Supplementary Table S2) showed a significant change in expression (P < 0.05) after multiple correction. Of these, the top 20 statistically significant genes are shown in Table 2. A heatmap of the 61 differentially expressed genes is presented in Figure 3C, which shows a difference in the grouping of color gradient between myopic and control eyes. On the other hand, no significant changes were observed when comparing the expression of these 348 genes between the peripheral temporal sclera of myopic and fellow control eyes. 

Table 2

 

Top 20 Differentially Expressed Genes Between the PPZ of FD (Myopic) and PPZ of the Fellow (Control) Eyes in Guinea Pigs

 

Pathway Analysis: PPZ Versus Peripheral Temporal Sclera Within Myopic Eyes

Pathway analysis of the 348 differentially expressed significant genes between the PPZ and peripheral temporal sclera within the myopic eyes identified over 150 different pathways, with the top 10 significant ones presented in Table 3. The top pathway identified was the Gαi signaling pathway (P = 5.13 × 10⁻⁴), where eight associated genes (Kras, Aplnr, Rgs7, S1pr3, Adra2a, Prkacb, Chrm2, Adora1) were identified as differentially regulated, including upregulation of Chrm2 (FC: +15.65, P < 0.05). Secondly, a total of eight novel genes (Kcnq2, Tac1, Gria3, Gria1, Gria2, Prkacb, Kcnq3, and Slc1a4) identified in our study are part of the glutamate receptor signaling pathway (P = 1.12 × 10⁻²), a pathway that is also part of the larger neuropathic pain signaling (P = 8.51 × 10⁻⁴) pathway. Thirdly, the genes (Kcnq2, Gabrg2, Slc6a1, Abat, and Kcnq3) from the GABA receptor signaling pathway (P = 3.55 × 10⁻³) were also identified as differentially regulated between PPZ and peripheral temporal sclera in myopic eyes. Finally, Fzd3 (FC: +3.08, P < 0.05), a gene previously associated with axial length,³¹ was identified as one of 13 genes differentially regulated that also overlapped with the axonal guidance signaling pathway (P = 1.66 × 10⁻²). The expression pattern of genes from the above pathways are presented in Supplementary Figure S1. 

Table 3

 

Top 10 Pathways Underlying Differential Gene Expression of the 348 Significant Genes Between the PPZ and Peripheral Temporal Sclera of Myopic Guinea Pigs

 

Pathway Analysis: PPZ of Myopic and PPZ of Fellow Control Eyes

Pathway analysis of the 61 significant genes differentially expressed between the PPZ of myopic and PPZ of fellow control eyes presented over 50 different pathways. The top 10 signaling pathways and their P values are presented in Table 4, which shows 5 pathways as statistically significant (P < 0.05 after multiple correction). Among these five significant pathways, three overlapped with the significant pathways that emerged from the comparison of PPZ and peripheral temporal sclera within the myopic eyes. They were Gαi signaling and catecholamine and L-dopachrome biosynthesis pathways. Genes identified in these three pathways and their FC include Aplnr and Adora1 (FC: −4.4 and −8.4, respectively, Gαi signaling), Dbh (FC: +8.5, catecholamine biosynthesis), and Tyr (FC: −7.4, L-dopachrome synthesis). In addition, two other pathways also presented with statistical significance (P < 0.05): eumelanin biosynthesis and hepatic fibrosis/hepatic stellate cell activation pathways. They harbored genes such as Tyr (FC: −7.4), Myh7 (FC: −25.9), Col20a1 (FC: −20.7), and Ccl21 (FC: −21.7). Supplementary Figure S2 shows the gene expression patterns in these significant pathways. 

Table 4

 

Top 10 Pathways Underlying Differential Gene Expression of the 61 Significant Genes Differentially Expressed Between the PPZ of Myopic and PPZ of Control Eyes

 

Validation Using qRT-PCR

To validate our findings, we undertook qRT-PCR on an independent set of samples. The four genes we chose were novel genes: Kif5c and Ncan from Table 1 along with Shox2 and Osgin1 from Table 2, which had not previously been shown to be involved in scleral changes during myopia. Additionally, genes that had been previously known to be involved in scleral extracellular matrix remodeling during myopia such as Col12a1³² and Timp3 (Gao H, et al. IOVS 2008;49:ARVO E-Abstract 1738), which do not feature in Tables 1 and 2, were also studied for their differential expression patterns. Supplementary Table S3 shows the primer sequences for these four genes. Table 5 shows that three out of four genes showed a similar direction of gene regulation between the RNAseq and qRT-PCR data. Such validation of a few genes, if not all, adds more confidence in the differential gene expression data obtained from RNAseq. Additionally, we report a downregulation of Col12a1 in the PPZ sclera of FD eyes compared to the PPZ of the fellow controls (−1.37) and when compared to the peripheral temporal sclera of the FD eyes (−1.78), whereas Timp3 was upregulated in the PPZ of FD compared to the PPZ of fellow controls (+1.84) and when compared to the peripheral temporal sclera of FD eyes (+1.62). 

Table 5

 

Differential Gene Expression Patterns From RNAseq and qRT-PCR

 

Discussion

A total of 61 genes were differentially expressed between myopic and control eyes (P < 0.05 after multiple correction) in a critical area of the sclera located around the optic nerve (PPZ) that shows significant and rapid expansion during myopia induction. Regional differences in gene expression between this PPZ region and the relatively unaffected peripheral temporal sclera existed in myopic eyes, with 348 genes being differentially expressed (P < 0.05 after multiple correction). 

Pathway analysis of the significant genes differentially expressed between the PPZ and peripheral temporal sclera within the FD eyes identified various pathways, including Gαi signaling and catecholamine and L-dopachrome biosynthesis pathways. Genes within these three pathways were also found to be significantly differentially expressed when comparing the PPZ sclera between the myopic and fellow control eyes. 

Activity within the Gαi signaling pathway likely plays an important role in myopia development. The Gi alpha subunit (Gαi, or Gi/G0 or Gi protein) inhibits the synthesis of cyclic adenosine monophosphate (cAMP) from adenosine triphosphate (ATP).³³ cAMP has been shown to inhibit the synthesis of collagen³⁴ and myofibroblast differentiation in cardiac fibroblasts.³⁵ In the sclera of normal guinea pigs, activation of cAMP by subconjunctival treatment with forskolin resulted in a myopic refraction accompanied by downregulation of collagen mRNA.³⁶ Additionally, human scleral fibroblasts treated directly with forskolin resulted in reduced collagen synthesis.³⁶ 

Our finding that the Chrm2 gene from Gαi signaling was upregulated (FC: +15.7) in the PPZ of myopic sclera relative to the periphery is also in agreement with the upregulation found in the sclera of myopic mice and the finding that knockout mice for this gene showed inhibited ocular growth compared to the wild-type mice.³⁷ However, unique to the current study, many other genes from the Gαi signaling pathway were affected in the expanding scleral zone (Supplementary Figs. S1A, S2E). Many G protein–coupled receptors bind to the Gi subunit, one class of which are the adenosine receptors. Adora1 is one such receptor for adenosine with activity mediated by Gαi, and we found significant downregulation of this gene in myopic PPZ (FC: −8.4) relative to the same region in untreated eyes. Given that our RNAseq data showed significant downregulation of one of the isoforms of collagen (Col20a1) in the PPZ of myopic compared to the fellow control eyes (FC: −20.7), this suggests that downregulation of genes from the Gαi signaling pathway such as Aplnr (another a member of the G protein–coupled receptor gene family) and Adora1 may impact collagen synthesis via activation of cAMP. 

However, one of the major subtypes of collagen, Col1a1, that was consistently downregulated in the sclera of myopic guinea pigs⁶ and other animal models,² was not found to be altered in our study. This could be because the scleral tissues that are being compared were from two different regions within already myopic eyes (PPZ and peripheral temporal) as opposed to the sclera from myopic and nonmyopic eyes that have previously been shown to exhibit altered scleral Col1a1 expression. If Col1a1 is downregulated to a similar extent across different regions of myopic sclera of guinea pigs, not only would we not expect to detect it using the current approach, but the implication is that such changes do not underlie the differential expansion between the PPZ and peripheral regions. 

The current study focused on genes that are specifically different between these two zones, which presumably is the cause of not identifying a vast number of differentially expressed scleral matrix genes. Additionally, the sclera must possess contractile ability to mediate rapid expansion at the PPZ accompanied by inhibited growth at the periphery. This may be facilitated by cell types other than fibroblasts that are capable of rapid contractions, such as smooth muscle cells. Our data showed that genes involved in cell growth, proliferation and smooth muscle contraction were distributed across the 348 significantly expressed genes at the PPZ in response to form deprivation, supporting such an argument. In addition, our finding of differential expression of genes involved in neuronal activity, as evidenced by genes from the neuropathic pain signaling pathway, supports the earlier finding that showed a possible role of neuronal developmental process toward myopia progression in humans,³⁸ thus proving a common signaling mechanism between humans and animal models. 

In agreement with our findings, parts of the GABA receptor signaling pathway have been previously implicated in myopia development in the fibrous and cartilaginous sclera of chicks.³⁹ Earlier studies have shown antagonists to GABA receptors (GABA_A, GABA_A0r, and GABA_B) to inhibit myopia development whereas agonists increased myopia susceptibility in chicks.^40,41 In guinea pigs, the selective GABAc inhibitor, TPMPA, inhibits ocular growth and myopia (McFadden SA, et al. IOVS 2011;52:ARVO E-Abstract 6306). More recently it was shown that one of the major scleral matrix components, glycosaminoglycan (GAG), was reduced in cultured chick scleral fibroblasts treated with GABA agonists (baclofen and muscimol) and increased in the presence of GABA antagonists (CGP46381, bicuculline, and TPMPA).⁴² Additionally, GABA receptors such as GABA_A, GABA_B, and GABA_C receptors have been identified in the human sclera.⁴³ However, ours is the first study, to our knowledge, to identify other genes involved in the GABA receptor signaling pathway (Gabrg2, Abat, Kcnq3, Slc6a1, and Kcnq2) as being affected in myopic sclera. The reduced scleral matrix content in the myopic sclera and the effect of GABA receptors in altering matrix components suggest a possible role of these additional genes involved in the GABA pathway in mediating these extracellular matrix changes. 

The Gria1 gene that was identified in neuronal pain signaling and glutamate receptor signaling has been reported in the human sclera,⁴⁴ thus presenting a common scleral gene between animal models and humans. The Fzd3 gene from the axonal guidance signaling is part of the frizzled gene family that interacts with the Wnt signaling pathway, which are associated with ocular growth.^45,46 The above genes were significantly upregulated at the PPZ compared to the peripheral temporal sclera in FD guinea pig eyes (Supplementary Table S1). While the normal ocular development in guinea pigs follows a prolate pattern with increased axial and scleral growth along the central visual field and particularly about the PPZ region,¹⁴ the development of myopia tends to accelerate this process via molecular mechanisms as evidenced by numerous significant gene expression changes particular to this region. Additionally, a large-scale genome-wide associated study (GWAS) carried out in humans by the Consortium for Refractive Error and Myopia (CREAM) identified significant associations between genes involved in Wnt signaling and axial length elongation,⁴⁷ suggesting a possible involvement of a common pathway between human and animal models of myopia. 

Previous studies on the sclera of animal models such as tree shrews showed that 4 days of myopia induction resulted in downregulation of genes involved in extracellular matrix secretion (TGF-β1 and TGF-β2), major extracellular matrix genes such as collagen subtypes (Col1a1, Col12a1), and proteoglycans (ACAN, DCN, FMOD, KERA, NYX, OGN). This was accompanied by upregulation of the protease enzyme MMP14 and downregulation of its endogenous inhibitor TIMPs (TIMP1 and TIMP3).^9,10 While changes to major scleral extracellular matrix components have been reported in myopic guinea pig sclera,^6,48 we did not identify significant regulation of these genes at the PPZ sclera compared to the peripheral temporal sclera in FD guinea pig eyes. If similar magnitudes in gene expression in response to form deprivation occurred between these two regions, they were excluded from our analysis. Specifically, the current study restricted the gene analysis to the unique gene set associated with the PPZ, while earlier studies analyzed all scleral genes encountered at the posterior pole. It is noted that some, but not all, of the posterior pole is included in the PPZ region. This in turn suggests that different mechanisms/pathways may mediate the specific changes that readily occur in the zone around the optic nerve, which in human myopia is known to be ultimately vulnerable to the formation of staphyloma.⁴⁹ 

In the current study, three animals per experimental condition were used to perform RNAseq. Similar small sample sizes have been previously used in RNAseq analysis of three retinas from leucine zipper knockout mice and three retinas from wild-type mice, with differential gene expression from RNAseq successfully validated using qRT-PCR,⁵⁰ which shows confidence in the RNAseq data. Identification of the relation between sample size and the depth of sequencing showed that a sample size of 4 and 10 million reads per sample could detect over 80% of annotated genes in chicken lungs.⁵¹ In the current study, the depth of coverage provided by the Illumina HiSeq 2000 platform was between 12.5 and 15 million reads per sample, which provided enough coverage to detect gene expression profiles even in transcripts of relatively low abundance. Given that the advantage of RNAseq is to identify the expression of novel gene transcripts, we chose a few such genes from our study, which have not been previously implicated in scleral remodeling during myopia, for validation purposes using qRT-PCR. A similar direction of gene regulation between the two techniques, in three out of four genes (75%) investigated, adds to confidence in the RNAseq data. While the expression of numerous collagen subtypes has been investigated in the past during myopia, we chose to study the expression of Col12a1 and the endogenous inhibitor of matrix metalloproteinase (MMP), Timp3, which are consistently shown to be differentially regulated across various myopia induction paradigms, including bidirectional changes in the sclera in response to hyperopic and myopic defocus.³² Downregulation of Col12a1 reported in our qRT-PCR data in both PPZ FD versus PPZ fellow eye and PPZ FD versus peripheral temporal FD fellow eye experiments is consistent with the above studies on matrix degradation observed during myopia. In the case of Timp3, upregulation was observed in both PPZ FD versus PPZ fellow eye and PPZ FD versus peripheral temporal FD fellow eye experiments, which is contrary to findings previously reported. A possible explanation is that this particular subtype of TIMP may undergo isoform-specific changes at the rapidly expanding PPZ region of the guinea pig sclera. However, further investigation of other subtypes of TIMPs is required to support this claim. 

The current study used samples taken after 2 weeks of FD, as this corresponds to a time point when significant and rapid increases occur in the PPZ region and opposing changes in the periphery have commenced, but it is possible that earlier molecular changes might be different from those observed here. Given that this is the first study to undertake RNAseq in different areas of the guinea pig sclera in response to myopia induction, further targeted investigation across a range of specific time points would therefore be needed to shed light on whether the gene expression profiles specific to myopia development in the PPZ sclera might differ at different development times. 

In summary and in conclusion, RNAseq provided a useful technique to begin to delineate gene expression from different regions of the sclera. Specifically, differential gene expression was found in the sclera surrounding the optic nerve, an area that demonstrates a rapid elongation during early myopia development in the guinea pig. The findings of this study highlight three potential pathways that may be critical in contributing to the excessive ocular growth that accompanies the development of myopia. It is worth noting that further biological validation will be necessary to fully explore the full list of genes identified through this RNAseq study, but it provides a starting point to interrogate potential genes of interest involved in early stages of myopia identified from across the genome. 

Acknowledgments

Supported by National Health and Medical Research Council of Australia Senior Research Fellowship 1028444 and 1138585 (PNB) and grants from the Coal Port Authority and Hunter Medical Research Institute (SAM), also a Research Higher Degree Studentship and Melbourne International fee remission scholarship from the University of Melbourne Australia (NS). The Centre for Eye Research Australia (CERA) receives Operational Infrastructure Support from the Victorian Government. 

Disclosure: N. Srinivasalu, None; S.A. McFadden, None; C. Medcalf, None; L. Fuchs, None; J. Chung, None; G. Philip, None; A. Richardson, None; M. Riaz, None; P.N. Baird, None