Network-based bioinformatics analysis of spatio-temporal RNA-Seq data reveals transcriptional programs underpinning normal and aberrant retinal development
- Devi Krishna Priya Karunakaran†1,
- Sahar Al Seesi†2,
- Abdul Rouf Banday1,
- Marybeth Baumgartner1,
- Anouk Olthof1, 3,
- Christopher Lemoine1,
- Ion I. Măndoiu2 and
- Rahul N. Kanadia1Email author
© The Author(s). 2016
Published: 31 August 2016
The retina as a model system with extensive information on genes involved in development/maintenance is of great value for investigations employing deep sequencing to capture transcriptome change over time. This in turn could enable us to find patterns in gene expression across time to reveal transition in biological processes.
We developed a bioinformatics pipeline to categorize genes based on their differential expression and their alternative splicing status across time by binning genes based on their transcriptional kinetics. Genes within same bins were then leveraged to query gene annotation databases to discover molecular programs employed by the developing retina.
Using our pipeline on RNA-Seq data obtained from fractionated (nucleus/cytoplasm) developing retina at embryonic day (E) 16 and postnatal day (P) 0, we captured high-resolution as in the difference between the cytoplasm and the nucleus at the same developmental time. We found de novo transcription of genes whose transcripts were exclusively found in the nuclear transcriptome at P0. Further analysis showed that these genes enriched for functions that are known to be executed during postnatal development, thus showing that the P0 nuclear transcriptome is temporally ahead of that of its cytoplasm. We extended our strategy to perform temporal analysis comparing P0 data to either P21-Nrl-wildtype (WT) or P21-Nrl-knockout (KO) retinae, which predicted that the KO retina would have compromised vasculature. Indeed, histological manifestation of vasodilation has been reported at a later time point (P60).
Thus, our approach was predictive of a phenotype before it presented histologically. Our strategy can be extended to investigating the development and/or disease progression of other tissue types.
The retina has been the most accessible part of the developing central nervous system with a wealth of information on detailed birth order of its cell types and on many genes involved in executing specific programs such as cell cycle regulation, cell fate determination, and neuronal differentiation. However, a comprehensive gene regulatory network is still not achieved as gene-centric approach can only go so far. To address this issue, transcriptome capture to identify co-transcriptionally regulated genes across retinal development has previously been attempted and was of great value . However, these efforts were hampered by the lack of depth of the captured transcriptome and lack of fractionation to gain higher resolution. Another concern was that at any given time the retina consists of different cell types with varied transcriptomes, which renders finding meaning from co-transcriptionally regulated genes difficult. We wanted to investigate whether higher depth of the captured transcriptome through RNA-Seq with minimal cross-compartment (nucleus-cytoplasm) normalization could resolve this issue.
Here we report analysis of RNA-Seq data from cytoplasmic and nuclear transcriptome of the developing retina. We show that combinatorial use of RNA-Seq with our custom bioinformatics strategy reveals the precise order of gene activation and transitions in processes during retinal development. Transition in gene expression was validated and resolved at the isoform level through our custom microarray. Importantly, we show proof of principle by extending our methodology to analyze RNA-Seq data from P21-Nrl-WT and KO retinae. Our approach which focuses on understanding the temporal progression in gene expression during normal/aberrant development can be extended to development and disease progression of other tissues.
All experiments used CD1 mice from Charles River Laboratory, MA. All mice procedures were compliant with the protocols approved by the University of Connecticut’s Institutional Animal Care and Use Committee (IACUC).
Retinae were dissected from E16 embryos and P0 pups followed by fractionation protocol as described previously . Once the fractions were obtained, Trizol (Invitrogen, CA, cat # 15596-026) was used as per the manufacturer’s instructions.
Library preparation for deep sequencing
After the total RNA was prepared from the two fractions, ribosomal RNA (rRNA) was removed using Ribozero Ribosomal RNA removal kit (Epicenter, WI, cat # RZH1046) by following the manufacturer’s instructions. The removal of rRNA was confirmed by gel electrophoresis and was used for RNA-Seq library preparation. RNA-Seq library was prepared using Script-seq mRNA seq library preparation kit (Cambio, UK, cat # SS10906). The library was deep sequenced in multiple runs using Illumina Hi-Seq 2000 platform at the University of Connecticut Health Center Deep sequencing core facility. P21 Nrl- WT and KO RNA-Seq data was shared with us by Dr. Anand Swaroop; National Eye Institute .
CD1 reference creation
Read mapping statistics and rRNA levels in the E16CE, P0CE, and P0NE samples
Percentage of transcriptome mapped read pairs
Percentage of rRNA reads
# mapped bases in Gb
Gene expression analysis
E16CE, P0CE and P0NE reads were mapped against the CD1 Ensembl 68 transcriptome reference. The P21 WT and KO single end reads were mapped against the C57BL6 reference transcriptome based on Ensembl version 68. Mapping was done using bowtie and allowed for one mismatch in an alignment seed of 30 bases. Gene expression levels were estimated using IsoEM , an expectation-maximization algorithm that estimates isoform frequency from single and paired RNA-Seq reads. IsoEM exploits read disambiguation information provided by the distribution of insert sizes generated during sequencing library preparation, and takes advantage of base quality scores, strand, and read pairing information. Isoform expression is reported as Fragment per Kilobase per million mapped reads (FPKM) units and gene expression is the sum of FPKM of its constituent isoforms. For gene differential expression, two methods were run, GFOLD  and Fisher’s exact test with house-keeping gene normalization as in . Gapdh was used as the housekeeping gene for this analysis. Genes were called differentially expressed if they showed ≥2 fold expression in one sample by both methods. GFOLD was run on the CD1 transcriptome aligned reads, with default parameters and a p-value of 0.01. Fisher’s exact test was run on estimated number of reads mapped per kilobase of gene length (calculated from IsoEM estimated FPKM values). Similar to GFOLD, a p-value of 0.01 was used for Fisher’s exact test.
Functional annotation analysis
We designed a custom Affymetrix microarray to en masse interrogate the presence/absence of unique exon/exon junctions in isoforms expressed in the RNA-Seq data. After mapping the RNA-Seq reads from the three samples (E16-CE, P0-CE and P0-NE) to the Ensembl 68 transcripts and running IsoEM to estimate the FPKM values for each of the three samples, genes expressed in any of the three samples were selected. Exon-exon junctions that are unique among expressed transcripts in genes that have more than one expressed transcript were selected (junctions with flanking sequences of length ≤12 bases were eliminated). As a result, we included probes for 28,574 unique junctions from 11,923 transcripts on the custom Affymetrix chip.
Custom affymetrix data analysis
Cytoplasmic RNA (1 μg) was prepared from retinae harvested from E12, E16, E18 embryos and P0, P4, P10 and P25 and processed at Yale Center for Genome Analysis. Expression levels of probe targets were computed from the raw intensity values using the Robust Multichip Average (RMA) algorithm [15, 16], which was performed with the affy R package . Subsequently, the data were processed through Gene Expression Similarity Investigation Suite (Genesis) (v1.7.6) for k-means clustering . Here, we ran 500 iterations to generate a total of 10 clusters that fell into three categories based on expression kinetics. These trends were defined as embryonic, postnatal, and embryonic + postnatal. The gene lists belonging to each cluster were separately analyzed for functional enrichment using DAVID.
RNA-Seq of fractionated retina
Skeletal muscle-specific genes expression. The table shows FPKM values of skeletal muscle-specific genes in E16CE, P0CE and P0NE samples
Validation of RNA-Seq data
We also used the same paradigm of genes with established expression kinetics to determine the level of cross-contamination between P0CE and P0NE RNA. To determine the level of nuclear RNA contaminating the cytoplasmic extract, we checked the expression of genes whose transcripts are predominantly nuclear, such as Xist, Malat1, Tsix and Neat1 [2, 35–37]. Indeed, for all four genes, the FPKM values were significantly higher in the NE compared to the CE (Fig. 3c). Determining the level of cytoplasmic RNA contamination in the nuclear extract presented a unique challenge as the majority of the RNA in the CE would be expected to be in the NE. For this, we examined expression of replication-dependent histone genes, as they are intronless and are known to be efficiently exported to the cytoplasm . Indeed, histone genes have higher FPKM in the CE than the NE (Fig. 3d). Furthermore, replication-dependent histone genes also serve to account for genomic DNA contamination. These genes lack introns and do not require splicing, so histone genomic DNA would be read as histone mRNA and inflate FPKM values in the NE, which was not the case (Fig. 3d), thus confirming minimal genomic contamination. Finally, genomic DNA contamination would result in FPKM value >0 for all genes; however, we observed 0 FPKM in the NE for a large number (804) of genes. In all, these controls suggest that there was minimal genomic DNA contamination in our fractionated NE RNA-Seq data.
RNA-seq revealed high-resolution transcription kinetics
Distribution of 2007 transcripts in P0NE_Only sample
Type of RNA
Protein coding genes
Transcriptionally coupled genes revealed molecular programs in the developing retina
One caveat to our binning strategy was that while it grouped genes based on transcription kinetics, it may have separated genes participating in a common biological process (identified by DAVID) into different bins. This in turn would prevent extraction of the expression kinetics of genes known to participate in executing a biological process of interest. To address this issue, we devised an iterative approach using GeneMANIA (Fig. 2c-II) to identify genes participating in a common biological process from the different bins. For example, 17 genes in OR_P0CE genes (E16CE-P0CE) that enriched for visual perception by DAVID were used as bait in GeneMANIA, followed by the aforementioned iterations (Fig. 2c; 1x-3x) to generate a final list of 36 genes (Fig. 2c, Right). Subsequently, each gene was assigned to its respective bin in both E16CE-P0CE and P0CE-P0NE comparisons (Fig. 2d-III). Redistribution of the genes into their respective bins is shown in Fig. 2e. For visual perception, the transcription of Rdh5, which converts all-trans retinal to 11-cis retinal  was initiated at/before E16, shut down just before birth and was initiated again at birth as we find it in P0NE_Only bin in the P0CE-P0NE comparison (Fig. 2e). In contrast, Crb1, Cngb3, Pcdh15 and Rgs9, which play a role in photo transduction [40, 41] and structural support/maintenance of photoreceptors [42, 43], were transcribed at/before E16 and upregulated at P0 as they were over-represented in P0NE in P0CE-P0NE (Fig. 2e). Rp1, which is a photoreceptor-specific microtubule-associated protein , was the only gene that was transcribed at P0 (P0CE_Only in E16CE-P0CE) that continued to be upregulated as it was over-represented in P0NE in P0CE-P0NE (Fig. 2e). Transcription of Guca1a was initiated between E16 and P0 (P0CE_Only in E16CE-P0CE), except it was turned off before P0 (P0CE_Only in P0CE-P0NE) (Fig. 2e). Through this method we were able to deconstruct the precise activation of genes involved in many aspects of vision acquisition/phototransduction during embryonic development.
The same analysis was performed on all the bins in P0CE-P0NE comparison. Specifically, genes in the OR_P0NE bin showed enrichment for 120 GO terms of which one of them was “synapse” (Additional file 3: Table S2.2). This functional enrichment agrees with studies showing that synaptogenesis occurs postnatally in the rodent retina . Further analysis of genes underlying the GOterm “synapse” showed transcription initiation of the AMPA receptor subunit genes including, Gria1, Gria2 and Gria4 before/at E16 (Non_DR in E16CE-P0CE) (Additional file 1: Figure S3). Similarly Grik2, which encodes for a subunit of the ionotropic kainate receptor, was also initiated before/at E16 (Additional file 1: Figure S3). Gad2, which is necessary for the production of the inhibitory neurotransmitter GABA, was transcribed before E16, while Gad1 transcription was initiated after E16 prior to birth (Additional file 1: Figure S3). Overall, genes involved in formation of the presynaptic activity were activated mostly during embryonic development (Additional file 1: Figure S3). In contrast, genes involved in postsynaptic activity had overlapping transcriptional activation with a subset of genes (Grid2, Grid1, Grik5, Grin3a, Ryr2 and Shank2) that were specifically activated in P0NE (Additional file 1: Figure S3). Finally, employing the same analysis for genes in P0NE_Only, which reflected de novo transcription, enriched for 14 GO terms, of which voltage-gated calcium ion channel activity was one of the top hits (Additional file 3: Table S2.2). Again, this enrichment did agree with previous studies where it has been shown that calcium channel activity is crucial for the construction of functional synapses that occurs postnatally .
Extending the analysis to other time points through the custom microarray
Next, we applied DAVID analysis to each cluster and found that embryonic clusters (Clusters 1 and 2) enriched for functions such as cell cycle regulation and cell projection organization (Fig. 5d). For the embryonic + postnatal clusters (Clusters 3–9), the functional GO terms that were enriched were those required for cell cycle regulation and terminal differentiation of neurons, such as vesicle-mediated transport, synapse formation, negative regulation of apoptosis, and axon guidance (Fig. 5e). Finally, the sole postnatal cluster (Cluster 10) was the only one that enriched for functions such as visual perception, photoreceptor cell differentiation and sensory perception of light (Fig. 5f). In all, the isoform-specific microarray confirmed RNA-Seq findings and further revealed the complexity of alternative splicing employed by the developing retina.
Comparison with other analysis methods
Temporal analysis combined with static analysis of Nrl WT and KO RNA-Seq is more informative
The loss of Nrl results in cell-fate switch from rod to cone photoreceptors . This made Nrl-KO an ideal system to test our hypothesis that temporal comparison would yield more information than the static analysis. First we performed static comparison between P21-Nrl-WT and P21-Nrl-KO data (Additional file 1: Figure S4, Additional file 5: Table S4.1), similar to the one previously reported . The objective of transcriptomics analysis of wild-type and knockout tissue is to find genes undergoing change to reveal the resulting biological change in the absence of that gene. Surprisingly, DAVID analysis of genes undergoing dynamic changes in gene expression in static comparison of P21-Nrl-WT and KO enriched for a couple of generic functions that did not give any meaning in terms of the knockout phenotype (Additional file 3: Table S2.5).
Binning vs. DE based analysis of Nrl WT and KO data
Here we have introduced a new strategy for RNA-Seq data analysis. To understand the advantages afforded by this approach we compared our binning method to the more commonly used DE method. First we analysed the static comparison between Nrl-WT vs. Nrl-KO data by the binning method and the DE-based method. Here we found that the DE-based method showed a large number (1116) of genes with differential expression (Additional file 6: Table S5.1). Similar data was obtained by the binning method except that the DE genes were now in bins labelled OR and Only, which reduced the large list of DE genes into manageable quanta. Moreover, the inherent value of the bins is inferred transcriptional kinetics. However, the binning method has a cost, which is reflected in the few statistically significant biological functions enriched by DAVID analysis. In contrast, the DE-based method generates a large list, which requires the investigator to decide the fold-change that might be relevant for his/her study. The one advantage of the DE method is that it produces a wide range of functional enrichments (12 GO terms) for DE gene-list for DAVID analysis (Additional file 7: Table S6.3).
Next, we analysed temporal comparison between P0 vs. P21 Nrl-WT and P0 vs.P21 Nrl-KO by both the binning method and the DE-based method (Additional file 6: Table S5.2, S5.3). Here the binning method revealed gene transcription kinetics across time, which is inherently valuable for deciphering developmental changes in normal and aberrant situations. The DE-based method produced a large DE gene-list, but did not yield any change in gene transcription over time (Additional file 7: Table S6.1, S6.2). Thus, it would require further deconstruction of the list, which is not necessary by the binning strategy. Also, DAVID analysis of the binned data provided enrichment of biological functions tethered to the specific gene expression profile. This is of great value as it is one of the central goals of performing RNA-Seq to capture transcriptome change over time (Additional file 3: Table S2.3, S2.4) . In all, the two approaches of RNA-Seq analysis have merits that can be leveraged to effectively analyse RNA-Seq data.
Amongst the co-transcriptionally regulated genes identified by our binning strategy, genes that remained transcriptionally unaltered employed a higher degree of alternative splicing than those undergoing dynamic regulation. This suggests that during development the major transcription initiations and terminations might lay down the foundation, while the proteome diversity generated through alternative splicing might be engaged in resolving the finer details, such as neuronal subtype specification and terminal differentiation. Another advantage of the binning strategy was that the large transcriptome data was quantized, allowing us to interrogate the dataset for genes with established expression kinetics. The purpose of this was to challenge the binned data generated based on the 1 FPKM threshold for known gene expression patterns (Fig. 3a-b).
An example of a “P0NE_Only” gene, Ces5a, whose FPKM unit in P0NE is comparable to those of other genes with known established expression kinetics (Nrl, Nr2e3, Gngt2)
The intrinsic value of identifying co-transcriptionally regulated genes is the expectation that they might reveal the biological processes being executed by the developing retina. Our bioinformatics pipeline can deconstruct the order of activation of specific genes engaged in executing a specific biological process so that one can begin to generate gene regulatory networks underlying retinal development. A key feature of our pipeline is the use of GeneMANIA to find potential partners of the core set of genes from a specific bin that enrich for a function in our DAVID analysis (Fig. 2c). A priori, one would predict a progressive increase in the number of genes with sequential application of the GeneMANIA part of the pipeline (Fig. 2c). However, we observed that there was quick convergence in the number of partner genes (Fig. 2c, Right). This suggests that leveraging RNA-Seq data to remove genes that were not expressed in the retina enriched for those genes relevant to retinal development and function at the time point under investigation.
Next we applied our analysis pipeline to find co-transcriptionally regulated genes in the P0 and P21-Nrl-WT comparison and the P0 and P21-Nrl-KO comparison (Additional file 1: Figure S4, Additional file 5: Table S4). One of the salient features of this analysis was that temporal analysis was more informative than static comparison. One explanation is that temporal analysis created bins that were developmentally regulated, which through DAVID analysis revealed changes in biological processes. For example, there is no cell cycle occurring at P21 so the majority of the cell cycle genes should be inactivated. Indeed, we observe cell cycle in the P0_Only bin in both P0 vs P21-Nrl-WT and P0 vs. P21-Nrl-KO analysis (Additional file 3: Table S2.3, S2.4). These genes in static analysis would show up as not expressed. Similarly, genes in P21_Only bin enriched for functions such as ion channel activity, ion transport, visual perception, synapse, voltage-gated ion channel activity, neurotransmission and others (Additional file 3: Table S2.3, S2.4). This was as expected as the retina is fully functional at P21 compared to P0. The advantage of our strategy is that it allowed us to understand the progression in gene expression kinetics in normal development and leverage that to understand how this progression deviates in the knockout retina. When P21_Only bin (either P21WT_Only or P21KO_Only) was analyzed through DAVID, we found many functions that were common to both sets of comparison, except examination of the number of genes underlying these functions revealed that there were subtle differences between the two bins (P21WT_only and P21KO_Only) (Additional file 3: Table S2.3, S2.4). This suggested that while many of the functions remain unaltered in the KO, there are subtle changes in the manner in which they might be executed. For example, “visual perception” showed up in the P0 vs. P21-Nrl-WT and P0 vs. P21-Nrl-KO comparisons in the P21_Only bins (Fig. 6). There were 24 genes underlying enrichment of this function in the WT comparison (Fig. 6), while there were 22 genes in KO comparison (Fig. 6). Upon comparing the gene identities from both sets, subtle differences emerged that allowed us to find the biological meaning from change in a single gene such as the rod photoreceptor-specific gene, Gnat1, that was absent in the Nrl-KO retina, which lack rod photoreceptors . Finding Gnat1 through temporal analysis raises the question whether it would have been found in static analysis. While Gnat1 was present in the P21WT_Only bin in static analysis, the rest of the genes that would normally be part of the GO term “visual perception” were in the Non_DR bin. Thus, without a priori knowledge one would not find this specific gene out of the entire list of genes in the P21WT_Only. Temporal analysis combined with our gene expression and binning strategy followed by our custom bioinformatics pipeline was able to find these subtle changes, which in case of static analysis was not possible (Additional file 1: Figure S4). While one could find these subtle changes in the static analysis by looking at specific genes, it requires a priori knowledge. The advantage of doing whole transcriptome analysis is that one could find patterns computationally, which can be leveraged to obtain new insights without the need for a priori knowledge. For example, GO terms such as potassium channel complex, sodium channel activity, synaptic vesicle, calcium ion transport and regulation of blood pressure regulation (Additional file 3: Table S2.4) were enriched by genes in the P21_Only bin in P0 vs. P21-Nrl-KO comparison, but were absent in the WT comparison. This finding suggests that there are specific anomalies in the Nrl-KO retina. Given that in the Nrl-KO retina, the majority of the rod photoreceptors have converted to cone photoreceptors, changes in ion transport and synaptogenesis are to be expected [53, 54]. However, regulation of blood pressure seemed out of place for the Nrl-KO retina. Indeed, closer examination of the genes underlying this function revealed the need to examine vasodilation in the Nrl-KO retina. Notably, previous reports showed that the Nrl-KO retina develops dilated retinal blood vessels and leakage at P60 . Thus, this confirmed the prediction made through shifts in the molecular signatures identified by our temporal analysis. Importantly, our analysis predicted an outcome based on gene expression pattern changes occurring between P0 to P21 that manifests at P60.
In summary, we were able to extract shifts in biological processes (Fig. 4a) governed by precise changes in gene expression through our unique RNA-Seq data acquisition/analysis platform. Importantly, it showed that the nuclear transcriptome was temporally shifted ahead of the cytoplasmic transcriptome at a developmental timepoint (Fig. 4c). Overlapping these discoveries with those made by extending our strategy to P21-Nrl-WT and KO analysis was most fruitful when ∆/time was extracted. This strategy identified perturbation in the molecular signature that enabled prediction of a phenotype that would manifest histologically at a later time (Fig. 4b). Indeed, this strategy would be effective toward deconstructing the progression of molecular changes during aberrant development or the progression of pathogenesis of the retinal diseases and can be extended to other tissues.
We would like to acknowledge the deep sequencing core facility at the University of Connecticut Health Center and the Yale Center for Genome Analysis. We would also like to thank Dr. Anand Swaroop from the National Eye Institute for sharing the P21-Nrl-WT and P21-Nrl-KO RNA-Seq data.
This article has been published as part of BMC Genomics Volume 17 Supplement 5, 2016. Selected articles from the 11th International Symposium on Bioinformatics Research and Applications (ISBRA '15): genomics. The full contents of the supplement are available online https://bmcgenomics.biomedcentral.com/articles/supplements/volume-17-supplement-5.'
This work was supported by P30 from the National Institute of Neurological Disorders and Stroke (5P30NS069266 to R.N.K) and K99-R00 from the National Eye Institute (4R00EY019547 to R.N.K) and Agriculture and Food Research Initiative Competitive Grant from the National Institute of Food and Agriculture (2011-67016-30331 to I.I.M) and award from National Science Foundation and a Collaborative Research Grant from Life Technologies (IIS-0916948 to I.I.M). The Physiology and Neurobiology Department at the University of Connecticut will defray the cost of publication through the corresponding author of this manuscript.
Availability of data and materials
The sequencing data can be obtained by contacting RNK at email@example.com.
RNK designed the experiment. DKPK, AB, MB, AO, CL, RNK performed the experiments. IIM and SA developed the bioinformatics pipeline used for the analysis. DKPK, SA, RNK analyzed the data. DKPK, RNK and SA prepared the manuscript. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
- Blackshaw S, Harpavat S, Trimarchi J, Cai L, Huang H, Kuo WP, Weber G, Lee K, Fraioli RE, Cho SH et al. Genomic analysis of mouse retinal development. PLoS Biol. 2004;2(9):E247.View ArticlePubMedPubMed CentralGoogle Scholar
- Banday AR, Baumgartner M, Al Seesi S, Karunakaran DK, Venkatesh A, Congdon S, Lemoine C, Kilcollins AM, Mandoiu I, Punzo C et al. Replication-dependent histone genes are actively transcribed in differentiating and aging retinal neurons. Cell Cycle. 2014;13(16):2526–41.View ArticlePubMedPubMed CentralGoogle Scholar
- Brooks MJ, Rajasimha HK, Roger JE, Swaroop A. Next-generation sequencing facilitates quantitative analysis of wild-type and Nrl(-/-) retinal transcriptomes. Mol Vis. 2011;17:3034–54.PubMedPubMed CentralGoogle Scholar
- Schmieder R, Lim YW, Edwards R. Identification and removal of ribosomal RNA sequences from metatranscriptomes. Bioinformatics. 2012;28(3):433–5.View ArticlePubMedGoogle Scholar
- Flicek P, Amode MR, Barrell D, Beal K, Brent S, Chen Y, Clapham P, Coates G, Fairley S, Fitzgerald S et al. Ensembl 2011. Nucleic Acids Res. 2011;39(Database issue):D800–806.View ArticlePubMedGoogle Scholar
- Fujita PA, Rhead B, Zweig AS, Hinrichs AS, Karolchik D, Cline MS, Goldman M, Barber GP, Clawson H, Coelho A et al. The UCSC Genome Browser database: update 2011. Nucleic Acids Res. 2011;39(Database issue):D876–882.View ArticlePubMedGoogle Scholar
- Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10(3):R25.View ArticlePubMedPubMed CentralGoogle Scholar
- Duitama J, Srivastava PK, Mandoiu II. Towards accurate detection and genotyping of expressed variants from whole transcriptome sequencing data. BMC Genomics. 2012;13 Suppl 2:S6.View ArticlePubMedPubMed CentralGoogle Scholar
- Nicolae M, Mangul S, Mandoiu II, Zelikovsky A. Estimation of alternative splicing isoform frequencies from RNA-Seq data. Algorithms for molecular biology. 2011;6(1):9.View ArticlePubMedPubMed CentralGoogle Scholar
- Feng J, Meyer CA, Wang Q, Liu JS, Shirley Liu X, Zhang Y. GFOLD: a generalized fold change for ranking differentially expressed genes from RNA-seq data. Bioinformatics. 2012;28(21):2782–8.View ArticlePubMedGoogle Scholar
- Bullard JH, Purdom E, Hansen KD, Dudoit S. Evaluation of statistical methods for normalization and differential expression in mRNA-Seq experiments. BMC bioinformatics. 2010;11:94.View ArticlePubMedPubMed CentralGoogle Scholar
- da Huang W, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc. 2009;4(1):44–57.View ArticleGoogle Scholar
- da Huang W, Sherman BT, Lempicki RA. Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists. Nucleic Acids Res. 2009;37(1):1–13.View ArticleGoogle Scholar
- Warde-Farley D, Donaldson SL, Comes O, Zuberi K, Badrawi R, Chao P, Franz M, Grouios C, Kazi F, Lopes CT et al. The GeneMANIA prediction server: biological network integration for gene prioritization and predicting gene function. Nucleic Acids Res. 2010;38(Web Server issue):W214–220.View ArticlePubMedPubMed CentralGoogle Scholar
- Irizarry RA, Hobbs B, Collin F, Beazer-Barclay YD, Antonellis KJ, Scherf U, Speed TP. Exploration, normalization, and summaries of high density oligonucleotide array probe level data. Biostatistics (Oxford, England). 2003;4(2):249–64.View ArticleGoogle Scholar
- Irizarry RA, Bolstad BM, Collin F, Cope LM, Hobbs B, Speed TP. Summaries of Affymetrix GeneChip probe level data. Nucleic Acids Res. 2003;31(4):e15.View ArticlePubMedPubMed CentralGoogle Scholar
- Gautier L, Cope L, Bolstad BM, Irizarry RA. affy--analysis of Affymetrix GeneChip data at the probe level. Bioinformatics. 2004;20(3):307–15.View ArticlePubMedGoogle Scholar
- Sturn A, Quackenbush J, Trajanoski Z. Genesis: cluster analysis of microarray data. Bioinformatics. 2002;18(1):207–8.View ArticlePubMedGoogle Scholar
- Hart T, Komori HK, LaMere S, Podshivalova K, Salomon DR. Finding the active genes in deep RNA-seq gene expression studies. BMC Genomics. 2013;14:778.View ArticlePubMedPubMed CentralGoogle Scholar
- Consortium SM-I, Consortium SM-I. A comprehensive assessment of RNA-seq accuracy, reproducibility and information content by the Sequencing Quality Control Consortium. Nat Biotechnol. 2014;32(9):903–14.View ArticleGoogle Scholar
- Kurose H, Bito T, Adachi T, Shimizu M, Noji S, Ohuchi H. Expression of Fibroblast growth factor 19 (Fgf19) during chicken embryogenesis and eye development, compared with Fgf15 expression in the mouse. Gene Expr Patterns. 2004;4(6):687–93.View ArticlePubMedGoogle Scholar
- Behesti H, Papaioannou VE, Sowden JC. Loss of Tbx2 delays optic vesicle invagination leading to small optic cups. Dev Biol. 2009;333(2):360–72.View ArticlePubMedPubMed CentralGoogle Scholar
- Hufnagel RB, Riesenberg AN, Quinn M, Brzezinski JA, Glaser T, Brown NL. Heterochronic misexpression of Ascl1 in the Atoh7 retinal cell lineage blocks cell cycle exit. Mol Cell Neurosci. 2013;54:108–20.View ArticlePubMedPubMed CentralGoogle Scholar
- Chiodini F, Matter-Sadzinski L, Rodrigues T, Skowronska-Krawczyk D, Brodier L, Schaad O, Bauer C, Ballivet M, Matter JM. A positive feedback loop between ATOH7 and a Notch effector regulates cell-cycle progression and neurogenesis in the retina. Cell Rep. 2013;3(3):796–807.View ArticlePubMedGoogle Scholar
- Jin Z, Zhang J, Klar A, Chedotal A, Rao Y, Cepko CL, Bao ZZ. Irx4-mediated regulation of Slit1 expression contributes to the definition of early axonal paths inside the retina. Development. 2003;130(6):1037–48.View ArticlePubMedGoogle Scholar
- Yoshida S, Mears AJ, Friedman JS, Carter T, He S, Oh E, Jing Y, Farjo R, Fleury G, Barlow C et al. Expression profiling of the developing and mature Nrl-/- mouse retina: identification of retinal disease candidates and transcriptional regulatory targets of Nrl. Hum Mol Genet. 2004;13(14):1487–503.View ArticlePubMedGoogle Scholar
- Corbo JC, Myers CA, Lawrence KA, Jadhav AP, Cepko CL. A typology of photoreceptor gene expression patterns in the mouse. Proc Natl Acad Sci U S A. 2007;104(29):12069–74.View ArticlePubMedPubMed CentralGoogle Scholar
- Cheng H, Khan NW, Roger JE, Swaroop A. Excess cones in the retinal degeneration rd7 mouse, caused by the loss of function of orphan nuclear receptor Nr2e3, originate from early-born photoreceptor precursors. Hum Mol Genet. 2011;20(21):4102–15.View ArticlePubMedPubMed CentralGoogle Scholar
- Omori Y, Katoh K, Sato S, Muranishi Y, Chaya T, Onishi A, Minami T, Fujikado T, Furukawa T. Analysis of transcriptional regulatory pathways of photoreceptor genes by expression profiling of the Otx2-deficient retina. PLoS One. 2011;6(5):e19685.View ArticlePubMedPubMed CentralGoogle Scholar
- Swaroop A, Xu JZ, Pawar H, Jackson A, Skolnick C, Agarwal N. A conserved retina-specific gene encodes a basic motif/leucine zipper domain. Proc Natl Acad Sci U S A. 1992;89(1):266–70.View ArticlePubMedPubMed CentralGoogle Scholar
- Swaroop A, Kim D, Forrest D. Transcriptional regulation of photoreceptor development and homeostasis in the mammalian retina. Nat Rev Neurosci. 2010;11(8):563–76.View ArticlePubMedGoogle Scholar
- Swain PK, Hicks D, Mears AJ, Apel IJ, Smith JE, John SK, Hendrickson A, Milam AH, Swaroop A. Multiple phosphorylated isoforms of NRL are expressed in rod photoreceptors. J Biol Chem. 2001;276(39):36824–30.View ArticlePubMedGoogle Scholar
- Walther C, Gruss P. Pax-6, a murine paired box gene, is expressed in the developing CNS. Development. 1991;113(4):1435–49.PubMedGoogle Scholar
- Hsieh YW, Yang XJ. Dynamic Pax6 expression during the neurogenic cell cycle influences proliferation and cell fate choices of retinal progenitors. Neural Dev. 2009;4:32.View ArticlePubMedPubMed CentralGoogle Scholar
- Donohoe ME, Zhang LF, Xu N, Shi Y, Lee JT. Identification of a Ctcf cofactor, Yy1, for the X chromosome binary switch. Mol Cell. 2007;25(1):43–56.View ArticlePubMedGoogle Scholar
- Kay GF, Penny GD, Patel D, Ashworth A, Brockdorff N, Rastan S. Expression of Xist during mouse development suggests a role in the initiation of X chromosome inactivation. Cell. 1993;72(2):171–82.View ArticlePubMedGoogle Scholar
- Sasaki YT, Ideue T, Sano M, Mituyama T, Hirose T. MENepsilon/beta noncoding RNAs are essential for structural integrity of nuclear paraspeckles. Proc Natl Acad Sci U S A. 2009;106(8):2525–30.View ArticlePubMedPubMed CentralGoogle Scholar
- Marzluff WF. Metazoan replication-dependent histone mRNAs: a distinct set of RNA polymerase II transcripts. Curr Opin Cell Biol. 2005;17(3):274–80.View ArticlePubMedGoogle Scholar
- Simon A, Lagercrantz J, Bajalica-Lagercrantz S, Eriksson U. Primary structure of human 11-cis retinol dehydrogenase and organization and chromosomal localization of the corresponding gene. Genomics. 1996;36(3):424–30.View ArticlePubMedGoogle Scholar
- Kohl S, Baumann B, Broghammer M, Jagle H, Sieving P, Kellner U, Spegal R, Anastasi M, Zrenner E, Sharpe LT et al. Mutations in the CNGB3 gene encoding the beta-subunit of the cone photoreceptor cGMP-gated channel are responsible for achromatopsia (ACHM3) linked to chromosome 8q21. Hum Mol Genet. 2000;9(14):2107–16.View ArticlePubMedGoogle Scholar
- He W, Cowan CW, Wensel TG. RGS9, a GTPase accelerator for phototransduction. Neuron. 1998;20(1):95–102.View ArticlePubMedGoogle Scholar
- Pellikka M, Tanentzapf G, Pinto M, Smith C, McGlade CJ, Ready DF, Tepass U. Crumbs, the Drosophila homologue of human CRB1/RP12, is essential for photoreceptor morphogenesis. Nature. 2002;416(6877):143–9.View ArticlePubMedGoogle Scholar
- Ahmed ZM, Riazuddin S, Ahmad J, Bernstein SL, Guo Y, Sabar MF, Sieving P, Riazuddin S, Griffith AJ, Friedman TB et al. PCDH15 is expressed in the neurosensory epithelium of the eye and ear and mutant alleles are responsible for both USH1F and DFNB23. Hum Mol Genet. 2003;12(24):3215–23.View ArticlePubMedGoogle Scholar
- Liu Q, Zuo J, Pierce EA. The retinitis pigmentosa 1 protein is a photoreceptor microtubule-associated protein. J Neurosci. 2004;24(29):6427–36.View ArticlePubMedPubMed CentralGoogle Scholar
- Wachtmeister L. Oscillatory potentials in the retina: what do they reveal. Prog Retin Eye Res. 1998;17(4):485–521.View ArticlePubMedGoogle Scholar
- Jhaveri DJ, O’Keeffe I, Robinson GJ, Zhao QY, Zhang ZH, Nink V, Narayanan RK, Osborne GW, Wray NR, Bartlett PF. Purification of neural precursor cells reveals the presence of distinct, stimulus-specific subpopulations of quiescent precursors in the adult mouse hippocampus. J Neurosci. 2015;35(21):8132–44.View ArticlePubMedGoogle Scholar
- Mears AJ, Kondo M, Swain PK, Takada Y, Bush RA, Saunders TL, Sieving PA, Swaroop A. Nrl is required for rod photoreceptor development. Nat Genet. 2001;29(4):447–52.View ArticlePubMedGoogle Scholar
- Calvert PD, Krasnoperova NV, Lyubarsky AL, Isayama T, Nicolo M, Kosaras B, Wong G, Gannon KS, Margolskee RF, Sidman RL et al. Phototransduction in transgenic mice after targeted deletion of the rod transducin alpha -subunit. Proc Natl Acad Sci U S A. 2000;97(25):13913–8.View ArticlePubMedPubMed CentralGoogle Scholar
- Roepman R, Bernoud-Hubac N, Schick DE, Maugeri A, Berger W, Ropers HH, Cremers FP, Ferreira PA. The retinitis pigmentosa GTPase regulator (RPGR) interacts with novel transport-like proteins in the outer segments of rod photoreceptors. Hum Mol Genet. 2000;9(14):2095–105.View ArticlePubMedGoogle Scholar
- Karan S, Frederick JM, Baehr W. Novel functions of photoreceptor guanylate cyclases revealed by targeted deletion. Mol Cell Biochem. 2010;334(1-2):141–55.View ArticlePubMedGoogle Scholar
- Ivanova E, Muller U, Wassle H. Characterization of the glycinergic input to bipolar cells of the mouse retina. Eur J Neurosci. 2006;23(2):350–64.View ArticlePubMedGoogle Scholar
- Roger JE, Ranganath K, Zhao L, Cojocaru RI, Brooks M, Gotoh N, Veleri S, Hiriyanna A, Rachel RA, Campos MM et al. Preservation of cone photoreceptors after a rapid yet transient degeneration and remodeling in cone-only Nrl-/- mouse retina. J Neurosci. 2012;32(2):528–41.View ArticlePubMedPubMed CentralGoogle Scholar
- Strettoi E, Mears AJ, Swaroop A. Recruitment of the rod pathway by cones in the absence of rods. J Neurosci. 2004;24(34):7576–82.View ArticlePubMedGoogle Scholar
- Cheng H, Aleman TS, Cideciyan AV, Khanna R, Jacobson SG, Swaroop A. In vivo function of the orphan nuclear receptor NR2E3 in establishing photoreceptor identity during mammalian retinal development. Hum Mol Genet. 2006;15(17):2588–602.View ArticlePubMedPubMed CentralGoogle Scholar