Molecular Vision 2019; 25:345-358
Received 20 February 2019 | Accepted 03 July 2019 | Published 05 July 2019
1Emory Eye Center, Department of Ophthalmology, Emory University, 1365B Clifton Road NE, Atlanta GA; 2Center for Neuropathology and Prion Research, Ludwig Maximilian University of Munich, Germany; 3Department for Translational Brain Research, German Center for Neurodegenerative Diseases (DZNE), Munich, Germany
Correspondence to: Felix L. Struebing, Center for Neuropathology and Prion Research, Ludwig Maximilian University of Munich, Germany; Feodor-Lynen-Str. 23; 81671 Munich, Germany; Phone: +49 89 2180 78028; FAX: +49 89 2180 78037; email: firstname.lastname@example.org
Purpose: The goal of the present study is to provide an independent assessment of the retinal transcriptome signatures of C57BL/6J (B6) and DBA/2J (D2) mice, and to enhance existing microarray data sets for accurately defining the allelic differences in the BXD recombinant inbred strains.
Methods: Retinas from B6 and D2 mice (three of each) were used for the RNA sequencing (RNA-seq) analysis. Transcriptome features were examined for both strains. Differentially expressed genes between the two strains were identified, and bioinformatic analysis was performed to analyze the transcriptome differences between the B6 and D2 strains, including Gene Ontology (GO) analysis, Phenotype and Reactome enrichment, and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis. The RNA-seq data were then directly compared with one of the microarray data sets (Department of Defense [DoD] Retina Normal Affy MoGene 2.0 ST RMA Gene Level Microarray Database) hosted on GeneNetwork.
Results: RNA-seq provided an in-depth analysis of the transcriptome of the B6 and D2 retinas with a total of more than 30,000,000 reads per sample. More than 70% of the reads were uniquely mapped, resulting in a total of 18,100 gene counts for all six samples. A total of 1,665 genes were differentially expressed, with 858 of these more highly expressed in the B6 retinas and 807 more highly expressed in the D2 retinas. Several molecular pathways were differentially active between the two strains, including the retinoic acid metabolic process, endoplasmic reticulum lumen, extracellular matrix (ECM) organization, and the PI3K-Akt signaling pathway. The most enriched KEGG pathways were the pentose and glucuronate interconversions pathway, the cytochrome P450 pathway, the protein digestion and absorption pathway, and the ECM-receptor interaction pathway. Each of these pathways had a more than fourfold enrichment. The DoD Normal Retina Microarray Database provided expression profiling for 26,191 annotated transcripts for B6 mouse, D2 mouse, and 53 BXD strains. A total of 13,793 genes in this microarray data set were comparable to the RNA-seq data set. For the B6 and D2 retinas, the RNA-seq data and the microarray data were highly correlated with each other (Pearson's r=0.780 for the B6 mice and 0.784 for D2 mice). These results suggest that the microarray data set can reliably detect differentially expressed genes between the B6 and D2 retinas, with an overall accuracy of 91.1%. Examples of true positive and false positive genes are provided.
Conclusions: Retinal transcriptome features of B6 and D2 mouse strains provide a useful reference for a better understanding of the mouse retina. Generally, the microarray database presented on GeneNetwork shows good agreement with the RNA-seq data, but we note that any allelic difference between B6 and D2 mice should be verified with the latter.
The development of high throughput RNA sequencing (RNA-seq) [1,2] is revolutionizing transcriptome analysis. It provides for a direct assessment of gene expression with relatively low background, a potentially unlimited dynamic range, and the ability to detect low abundance transcripts . RNA-seq also provides opportunities to revisit, verify, and improve existing microarray data sets by reinforcing valid findings, and eliminating false positive results that can arise due to hybridization artifacts . Large-scale transcriptome data sets using microarray technology were developed, and have been widely used for decades [4-6]. With the help of the retinal transcriptome microarray databases generated from BXD mouse strains hosted on GeneNetwork [7-11], our group has examined molecular networks contributing to various glaucoma-related phenotypes. We used the BXD strains and tools hosted on GeneNetwork to identify Aldh7a1 as a gene modulating retinal ganglion cell (RGC) death following elevated intraocular pressure (IOP) . We also identified Cdh11 as a gene modulating normal IOP . Finally, we identified Pou6f2 modulating central corneal thickness (CCT) in the mouse as a potential risk of glaucoma in humans . Subsequently, two independent genome-wide association studies in humans identified CDH11 (Gene ID 1009; OMIM 600023) and POU6F2 (Gene ID 11281; OMIM 609062) as genes modulating IOP in humans, and as factors for glaucoma risk [15,16], supporting our findings, and proving the power of this systems biology approach and the microarray data sets hosted on GeneNetwork. The BXD strains represent a recombinant inbred strain set derived from C57BL/6J (B6) and DBA/2J (D2) parents, whose genomic variance forms the foundation of the BXD mice.
In the present study, we generated RNA-seq data from the retinas of B6 and D2 mice. The transcriptome features and differences were investigated. We also compared the RNA-seq data to our microarray data set (Affymetrix MouseGene 2.0 ST Array Gene Level) hosted on GeneNetwork. With the parental RNA-seq data, allelic differences in the BXD recombinant inbred strains that may arise due to hybridization artifacts could be precisely defined. The results suggest that a combination of microarray and RNA-seq data can lead to more valid results, compared to either data set alone.
For the RNA-seq analysis, a total of six samples were analyzed, three samples from B6 mice (two males and one female) and three samples from D2 mice (two males and one female). All the animals were between 60 and 100 days of age. Mice were maintained on a 12 h:12 h light-dark cycle in a parasite-free facility with food and water ad libitum. All procedures involving animals were approved by the Animal Care and Use Committee of Emory University, and were in accordance with the ARVO Statement for the Use of Animals in Ophthalmic and Vision Research.
Mice were deeply anesthetized with 15 mg/kg of xylazine and 100 mg/kg of ketamine, and then euthanized with rapid cervical dislocation. Retinas from both eyes were quickly removed under a dissection microscope and directly placed into 160 U/ml Ribolock® (Thermo Scientific, Walton, MA) in Hank’s Balanced Salt Solution (Sigma, St. Louis, MO) on ice. Tissue was stored at −80 °C. RNA was isolated in batches using a Qiacube and the RNeasy Mini Kit (Qiagen, Hilden, Germany) according to the manufacturer’s instructions. The isolation included on-column DNase1 treatment to remove contaminating genomic DNA. All tissue was harvested between 10 AM and 12 PM to minimize circadian differences in gene expression. RNA integrity (RIN) was assessed on a Bioanalyzer 2100 (Agilent, Santa Clara, CA). The RIN score for all samples ranged from 9.1 to 9.5. RNA was then quantified with spectrophotometry, and 260/280 ratios for all samples were >2.1.
The retinas were dissected from the eye, and during the process, there are traces of RPE still attached to the retina. These traces were small, and to reveal the lack of significant RPE contamination, we examined the expression of RPE genes in the samples. The expression of genes in Figure 1 represents several different ocular cell types, including RPE. These data can be used to determine the makeup of the samples and the low level of RPE contamination.
Six samples of purified retinal RNA underwent library preparation with the Illumina (San Diego, CA) TruSeq Stranded Total RNA kit. Paired-end sequencing (100 bp) was performed by Hudson Alpha (Huntsville, AL). FastQC of the files was performed for quality control. Reads from each replicate were not trimmed before mapping, as this was recently shown to decrease the ability to map, and correlation with microarray data . To avoid mapping issues arising from mouse strain heterogeneity, separate genome indices were built for the B6 and the D2 mouse genome. This was achieved by manually substituting the mm10 B6 reference genome (Ensembl version 86) with D2-specific single nucleotide polymorphisms (SNPs) and indels (from dbSNP 142). The reads were then mapped to each genome using the two-pass mapping mode in RNA-STAR with default parameters. For each sample, strand-specific bam files were sorted and indexed with samtools. Counts per gene were tabulated with the “summarizeOverlaps” internal R function from uniquely mapped reads only. Differential expression analysis was performed using a generalized linear model in edgeR  after any genes with fewer than 0.5 counts per million (cpm) in more than three samples were removed. All p values were adjusted using the false discovery rate (FDR). Statistical significance was assumed for genes with FDR<0.001 and a twofold change between strains. For comparison with the microarray data, log2-transformed cpm counts were z-scaled, and multiplied by 2, before a constant of 8 was added to avoid negative values. This is essentially the same normalization procedure that was applied to the microarray data sets hosted on GeneNetwork [6,19,20]. All raw data are deposited in the Gene Expression Omnibus (GEO), accession GSE127942.
Correlation tests were performed on genes extracted from the RNA-seq and microarray data sets. Gene symbols were used for the comparisons. In cases where there were duplicate gene symbols represented by different probes in the microarray annotation file, an average value of all probes capturing different parts of the transcripts was calculated to represent the expression level of the gene. The absolute differential expression values for each gene between the two data sets were calculated as the expression level of B6 minus D2 . When the RNA-seq and microarray data sets were compared in detecting differentially expressed genes, the original probe value was used instead of the average value of all probes for a gene.
The differentially expressed genes between the B6 and D2 retinas were submitted for enrichment analysis to databases, including Gene Ontology (GO) [22,23], MGI mammalian phenotype ontology , Human Phenotype Ontology , and Reactome biologic pathways  with the R package EnrichR [27,28]. The top hits of the package output ‘combined scores’ and adjusted p values were used to determine significance. The differentially expressed genes were also submitted for Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis [29-31]. Redundancy control was run for each database. For GO analysis, “REVIGO”  was used for redundancy reduction, before the results were plotted. For all other enrichment analysis databases, “ReCiPa”  was used to reduce redundancy. The redundant pathways were combined into the pathway with a broader spectrum of genes and greater significance.
Quality evaluation (FastQC) of the RNA-seq files showed that >90% of reads had a Phred quality score of >30 (Appendix 1). Each replicate had >30 million raw reads, of which more than 70% could be uniquely mapped, resulting in a total of 18,100 gene counts for all three samples in the B6 group and three samples for the D2 group. Any genes with fewer than 0.5 cpm in the combined data were removed. Because the retina is a highly heterogeneous tissue composed of various cell types, this cutoff was chosen to include RNA from less common cells in the retina, such as specific RGC subtypes. As the first survey of the data quality, we examined the expression levels of signature genes for different retinal cell types (four of each), including RGCs, amacrine cells, bipolar cells, Müller glia, and photoreceptors [11,34-36] (Figure 1). These genes were compared with marker genes of the cornea and RPE, which should show negligible expression in clean retinal tissues . Commonly used housekeeping genes  were also examined as reference control (Figure 1). Most of the retinal cell type-specific signature genes were highly expressed, compared with the housekeeping genes, especially the gene for rhodopsin (Rho) and other photoreceptor markers. The retinas were dissected from the eye, and during the process, there were traces of RPE still attached to the retina. These traces were small, and to reveal the lack of significant RPE contamination, we examined the expression of RPE signature genes in the samples. The expression of genes in Figure 1 represents several different ocular cell types, including RPE. These data are used to determine the makeup of the samples and the low level of RPE contamination.
We also examined two genes with known mutations in the DBA/2J strain contributing to eye phenotypes, tyrosinase-related protein 1 (Tyrp1) and glycoprotein non-metastatic melanoma protein B (Gpnmb) [38,39]. Gpnmb showed a difference in expression between the two strains, while Tyrp1 and the other listed genes were similarly expressed between the B6 and D2 retinas (Figure 1).
Among the housekeeping genes examined, Actb, Gapdh, and Ppia were equally expressed across samples (low standard error of the mean [SEM] levels); however, Rn18s (Encoding 18S rRNA) displayed a surprising variability (Figure 1). This high variability of Rn18s is known to be associated with cell differentiation, forming organs and their maturation . These results suggest that the Actb, Gapdh, and Ppia are better choices as housekeeping genes for the retina.
Differentially expressed genes between the B6 and D2 retinas were extracted from RNA-seq data, and the p values were adjusted for multiple comparisons with the false discovery rate . We found 1,665 genes that were differentially expressed, with 858 of these more highly expressed in the B6 retinas and 807 more highly expressed in the D2 retinas (Appendix 2). The differentially expressed genes were submitted for enrichment analysis, including GO [22,23], Phenotypes in both Mouse  and Human , and Reactome  enrichment using EnrichR [27,28]. The top five pathways for each enrichment are shown in Figure 2. The pathways with statistical significance included the retinoic acid metabolic process and extracellular matrix organization in the GO biologic process (Figure 2A), the endoplasmic reticulum lumen in the GO cellular component (Figure 2C), and the extracellular matrix organization and collagen formation in the Reactome database (Figure 2F). The genes involved in each pathway are summarized in Appendix 3 and Appendix 4.
The differentially expressed genes were also submitted for KEGG analysis, in which the top 30 pathways are listed (Figure 3). The greatest number of differentially expressed genes fell into the metabolic pathway (n=107). The second most involved pathway was the PI3K-Akt signaling pathway (n=41, Appendix 5). There were 19 genes involved in the retinol metabolism pathway, with 16 genes expressed at a higher level in the B6 retinas and only two genes (Rdh9 and Cyp4a12b) showing higher expression in the D2 retinas (Appendix 6). The most enriched pathways were the pentose and glucuronate interconversions, the drug metabolism - cytochrome P450, the ECM-receptor interaction, alpha-linolenic acid metabolism, and protein digestion and absorption, with more than fourfold enrichment for each pathway (Figure 3, Appendix 7- Appendix 8-Appendix 9, Appendix 10, and Appendix 11).
RNA sequencing directly quantifies read counts, while microarray technology relies on the binding of a transcript to a predesigned probe sequence. To investigate these technical differences, and understand how they might be reflected in the mouse retinal transcriptome, we contrasted the RNA-seq data with the Department of Defense (DoD) Congressionally Directed Medical Research Programs (CDMRP) Normal Retina Database, which is publicly available on GeneNetwork. This previously established data set uses the Affymetrix MouseGene 2.0 ST Array, which provides expression profiling for 26,191 well-established annotated transcripts, 9,049 non-coding RNAs, and more than 600 microRNAs . Details about the data processing can be found in King et al. . This data set consists of 52 BXD strains, as well as B6, D2, and an F1 cross between B6 and D2 (B6D2F1). A total of 222 microarrays were used to describe the retinal transcriptomes of these combined 55 mouse strains, which were usually represented by four biologic replicates per strain.
A total of 13,793 genes (represented by 16,445 probe sets) in the microarray data set were also identified in the RNA-seq data set (Appendix 12). In the microarray database, some genes were represented by multiple annotated probe sets. When we encountered this situation, the mean value of all probe sets was used to represent the expression level of the gene. The RNA-seq data and the microarray data showed a high correlation; the Pearson’s r value for the B6 retina was 0.780 (Figure 4A) and for the D2 retina was 0.784 (Figure 4B). Thus, there was good agreement between the two data sets, suggesting that gene expression values from either RNA-seq or microarrays are reliable and comparable, except a few genes with either high or low expression.
The absolute difference in the mean expression of genes between the B6 and D2 retinas from the microarray data set was calculated (see Methods and Appendix 13), and compared to the RNA-seq data, as shown in Figure 5. Pearson’s r value was 0.408 for the differentially expressed values of all genes (p<0.001), indicating moderate agreement between the filtered data sets. When the expression level with differences above 0.5 on a log2 scale (|B6-D2| >0.5) was considered the threshold for the microarray data, 848 genes were identified as differentially expressed, of which 387 (45.6%) were confirmed by the RNA-seq data as true positive, while 461 (54.4%) were identified as false positive (Figure 5). Simultaneously, 15,597 genes did not show differential expression in the microarray data set, of which 14,598 (93.6%) genes were confirmed as true negative, and 999 (6.4%) genes were identified as false negative by the RNA-seq data (Figure 5). Thus, the sensitivity and specificity of detecting differentially expressed genes using the microarray data set were 27.9% and 96.9%, respectively. The precision, i.e., the positive predictive value of the microarray data in detecting differentially expressed genes, was 45.6%, and the negative predictive value was 93.6%. The overall accuracy of the microarray data sets relative to the RNA-seq data was 91.1%.
Although the microarray methodology might overestimate the number of differentially expressed genes, the results suggest that the retinal transcriptome is generally well represented with both technologies. Expression values from microarrays depend on probe hybridization, which can be affected by differences in the genome sequence observed between strains. Thus, genomic variants overlapping the probe targeting area could directly affect the affinity of the probe to the transcript [41,42]. RNA-seq is not dependent upon hybridization of a target to a probe, and thus, this technology can serve to validate the results obtained using a microarray system.
An example for validation of the microarray data by RNA-seq is Stab2. In the microarray data set, clear differential expression of the gene Stab2 was identified between the parental strains, as well as across the BXD strains (Probe_17243763, Max LRS=27). The likelihood ratio statistic (LRS) is used as a measurement of the association or linkage between differences in RNA expression , where values exceeding 17 usually approximate statistical significance at the p=0.05 level after correction for multiple testing. The Stab2 probe coverage area contained no SNP, and it was also identified as a differentially expressed gene in the RNA-seq data (Figure 6).
In another case, a clear differential expression of the gene Anapc1 was identified by the microarray data set (Probe_17391461, Max LRS=124.5). However, the RNA-seq data showed no statistically significant differential expression for this gene. Probe verification using the UCSC Genome Browser showed an SNP (rs27412956) in the probe coverage area that may lead to mis-binding toward the D2 sequence. This potentially explains why Anapc1 showed lower expression levels in the D2 mice compared to the B6 mice in the microarray data set (Figure 7), suggesting that an SNP in a probe hybridization locus alters probe binding. Moreover, this mis-binding affected not only the D2 parental strains but also all the BXD strains with D2-derived alleles for this genomic region, resulting in a false-high LRS. In such cases, the RNA-seq data set is helpful in identifying false positive results. Similar false positive cases caused by D2 SNPs were seen in a total of 211 genes (222 probes, Appendix 14).
In the present study, we examined the retinal transcriptome of two different mouse strains, C57BL/6J and DBA/2J. Although it is known that there are clear sex differences in gene expression, we chose mixed genders for this study, as this most closely mirrored the sex distribution of BXD mice from the GeneNetwork databases . The B6 mouse is the most widely used mouse strain in research. In many ways, this strain serves as the ultimate inbred mouse strain, not only being the background strain for most transgenic lines but also functioning as the reference genome for the mouse. In contrast, the D2 mouse holds a specific interest for the neuroscience and the vision research community, as this strain is used in behavioral addiction research , and develops pigment dispersion glaucoma with onset at approximately 6 months of age [39,45]. When we examined the retinal transcriptome between these two strains using RNA-seq, several differences became apparent. Two mutations are known to be responsible for the DBA/2J eye phenotype: GpnmbR150X, a nonsense mutation, and Tyrp1isa, an allele leading to two amino acid substitutions. However, at the mRNA level, this must not necessarily mean that these two genes are differentially expressed, because we found differential expression only for Gpnmb. Although more research on this topic is warranted, it can be speculated that the nonsense mutation in Gpnbm leads to increased nonsense-mediated decay (NMD), and therefore, early mRNA removal by the NMD machinery , which can be seen as different expression between D2 and B6 mice. In contrast, the amino acid substitutions in the Tyrp1 gene do not lead to early transcriptional termination, and therefore, are tolerated at the mRNA level.
In the Gene Ontology analysis for differentially expressed genes, the most significantly enriched pathway was the retinoic acid metabolic process (Figure 2A). The main contributing genes (CYP26A1, CYP26B1, and CYP26C1) come from the CYP26 family, also known as retinoic acid hydroxylase. All three genes are statistically significantly higher expressed in B6 mice than in D2 mice. These genes are also the main contributors to the significant enrichment of the cytochrome P450 pathway and the retinol metabolism pathway in KEGG analysis (Figure 3). Another significantly enriched pathway is the ECM organization, which is also the only significantly enriched pathway of the Reactome analysis (Figure 2A,F). The ECM is known to play multiple roles in the retina, including retinal development [47,48], regulation of retinoid transport, oxygen and nutrient transport , vascular homeostasis , glia-neuron interactions , and even regeneration of RGCs . This difference in the ECM between the B6 and D2 retinas may potentially contribute to the differences in RGC survival , and other retinal phenotypes between the two strains, which makes the biologic approach using BXD strains (offspring inbred strains of B6 and D2) more powerful in revealing genetics of retinal diseases, such as glaucoma and age-related macular degeneration.
Pathway analysis using KEGG demonstrated differences in the PI3K-Akt signaling pathway, with 54 genes that were differentially expressed (Figure 3). The PI3K-Akt signaling pathway is an important intracellular signaling pathway for cell cycle regulation. Overactivation of this pathway can reduce apoptosis, allow proliferation, and is associated with specific cancers . This difference in the PI3K-Akt signaling pathway could affect changes in the ability of axons to regenerate. In the retina, Pten silencing is known to increase retinal ganglion cell survival, and promote axonal regeneration after injury . There is a strong difference in the ability of optic nerve regeneration across the BXD strains . Part of this may be due to transcriptome variance between the parental strains (B6 and D2) in this PI3K-Akt signaling pathway, for it is known that AKT is involved in the effects of Pten knockdown . An interesting feature of the Pten-dependent regeneration response is that the parental strains (B6 and D2) are not at the extremes of the phenotype . Thus, axon regeneration displays genetic transgression, and is a complex trait modulated by multiple genomic loci in the BXD strain set.
Furthermore, pathways involved in the innate immune system were also identified in the KEGG pathway analysis (Appendix 15). Many contributing genes were related to the complement system, including C1qb, C1qc, C1ra, C3, and C4b. The retina-intrinsic innate immunity network was reported to be activated by various types of retinal injures, including optic nerve crush  and ocular blast injury . The related genes were mostly higher expressed in B6 mice rather than D2 mice (Appendix 7), indicating a greater innate immune response in the retina of B6 mice compared to that of D2 mice.
The BXD strains are derived from the B6 and D2 progenitor strains, and their genomic variance allows us to map quantitative trait loci (QTLs) with high resolution, and economically, using microarray technology. Although this approach provides large amounts of information, it is not without limitations: Because the core principle behind microarrays is hybridization between two nucleic acid strands, the signal depends on the probe binding to the target sequence. Genomic variants between the B6 and D2 retinas, such as SNPs and indels, may influence hybridization kinetics, and cause incorrect detection of the expression level of genes [3,58]. A sequence that contains a single SNP may show a decreased binding affinity toward the probe, resulting in a less intense fluorescence signal, and thus, a lower readout for the probe capturing that particular transcript .
RNA-seq has developed rapidly over the past 10 years, and has now replaced the use of microarrays in transcriptomics [1,2,59,60]. The inherently lower background noise and higher dynamic range make RNA-seq more accurate in detecting low abundance transcripts . In this study, the RNA-seq data from the B6 and D2 mice showed high consistency when compared with the DoD Normal Retina Microarray Data set (Figure 4). Beyond that, a moderate to high correlation of the differentially expressed genes between the two data sets suggested that the microarray data set still performs well at detecting differences between different strains of mice (Figure 5).
Ultimately, the combination of the technologies benefits analysis, especially for the study of recombinant inbred strains: Microarrays are more affordable for high throughput experiments, while RNA-seq analysis of the parental strains can limit the high false positive QTL rate. Note that not all the SNPs have to cause mis-binding between a microarray probe and the cDNA. Microarrays that feature arbitrarily designed mismatch probes still detect approximately 70% of perfect match signals, indicating that hybridization is not a binary measure, compared to RNA-seq read counts . Thus, cross-checking results obtained from microarray data in recombinant inbred populations with RNA-seq data from the respective parental strains refines the analysis. After all, the genome sequences of the BXD strains are inherited from either the B6 or D2 mice. If a differentially expressed gene cannot be confirmed with RNA-seq data, and meanwhile, all the BXD strains with a D2 haplotype at that locus show lower expression, then there is a high chance for this finding to be false positive. If differential expression is not confirmed with the parental RNA-seq data, and the differential expression in BXD strains is not haplotype-dependent, then there are good chances that the finding is valid. These QTLs could be the consequence of inherited modulators in cis or trans (enhancers, suppressors, or other regulatory elements) that affect gene expression. Ultimately, if differential expression is confirmed by the RNA-seq data, then we can be more confident about a real change in gene expression that is allele-based.
Differential expression analysis of B6 and D2 mice using RNA-seq provided a useful reference for better understanding of the biologic differences between these strains. Although the majority of expression differences is mirrored in the microarray data from BXD strains, RNA-seq complements these data sets. Any allelic difference between B6 and D2 mice should be verified with RNA-seq data to avoid false positive QTLs.
Appendix 1. Fastqc plot of retina samples for RNA-seq analysis.
Appendix 2. RNAseq_B6D2-zscaled_Retina
Appendix 3. Enrichlist_DB
Appendix 4. KEGG LIST OF DEGs
Appendix 5. DE genes involvement of the PI3K-Akt signaling pathway.
Appendix 6. DE genes involvement of the Retinol metabolism pathway.
Appendix 7. DE genes involvement of the Pentose and glucuronate interconversions pathway.
Appendix 8. DE genes involvement of the Drug metabolism - cytochrome P450 pathway.
Appendix 9. DE genes involvement of the Alpha-Linolenic acid metabolism pathway.
Appendix 10. DE genes involvement of the ECM-receptor interaction pathway.
Appendix 11. DE genes involvement of the Protein digestion and absorption pathway.
Appendix 12. RNAseq_microarray_comparison_B6D2_Retina
Appendix 13. RNAseq_microarray_one-to-many-comparison_B6D2_Retina
Appendix 14. False positives caused by D2 SNPs
Appendix 15. DE genes involvement of the staphylococcus aureus infection pathway
We thank Rebecca King for her technical assistance with the mice and RNA isolation. This work was supported in part by the Vision Core Grant P030EY006360, the Emory Integrated Genomics Core (EIGC) Shared Resource of Winship Cancer Institute of Emory University and NIH/NCI under award number P30CA138292, the Owens Family Glaucoma Fund, and an Unrestricted Grant from Research to Prevent Blindness. FS is supported by the European Union’s Horizon 2020 Research and Innovation Program (grant agreement MSCA-IF-EF-RI 792832).