The zebrafish retina grows for a lifetime. Whether embryonic and postembryonic retinogenesis conform to the same developmental program is an outstanding question that remains under debate. Using single-cell RNA sequencing of ∼20,000 cells of the developing zebrafish retina at four different stages, we identified seven distinct developmental states. Each state explicitly expresses a gene set. Disruption of individual state-specific marker genes results in various defects ranging from small eyes to the loss of distinct retinal cell types. Using a similar approach, we further characterized the developmental states of postembryonic retinal stem cells (RSCs) and their progeny in the ciliary marginal zone. Expression pattern analysis of state-specific marker genes showed that the developmental states of postembryonic RSCs largely recapitulated those of their embryonic counterparts, except for some differences in rod photoreceptor genesis. Thus, our findings reveal the unifying developmental program used by the embryonic and postembryonic retinogenesis in zebrafish.

During embryonic development, the central nervous system (CNS) derives from a pool of embryonic neural progenitor cells through the neurogenesis process (Kriegstein and Alvarez-Buylla, 2009; Merkle and Alvarez-Buylla, 2006; Paridaen and Huttner, 2014; Taverna et al., 2014). In the mammalian CNS, such neurogenesis is rapidly diminished beyond embryonic development (Merkle and Alvarez-Buylla, 2006; Paridaen and Huttner, 2014). In contrast, fish and amphibians exhibit the remarkable capacity of postembryonic neurogenesis, enabling them to grow CNS tissues throughout their life. Postembryonic neural stem cells are the cell source for such neurogenesis (Marcucci et al., 2016; Perron and Harris, 2000), raising the question as to whether postembryonic neural stem cells employ the same developmental program as their embryonic counterparts.

The vertebrate retina, an accessible CNS tissue, is one of the most popular systems to study the development of a CNS tissue at the molecular, cellular and circuit levels (Amini et al., 2017; D'Orazi et al., 2014; Davis and Dyer, 2010; Hoon et al., 2014; Malicki et al., 2016; Wallace, 2011). The retina has a well-characterized cellular type and organization. Five major neuronal types [retinal ganglion cells (RGCs), amacrine cells (ACs), horizontal cells (HCs), bipolar cells (BCs) and photoreceptor cells (PRs)] are interconnected through two synaptic layers, the inner plexiform layer (IPL) and the outer plexiform layer (OPL). Müller glial cells (MCs) are the only glial type in the retina (Fig. 1A).

All retinal types are the products of multipotent retinal progenitor cells (RPCs). At the clonal level, individual RPCs give rise to the clones that are highly variable in terms of cell number and cell-type composition (Holt et al., 1988; Turner et al., 1990; Wetts and Fraser, 1988). This clonogenic ability is conserved by neural progenitors of other CNS tissues across different vertebrate species (Agathocleous and Harris, 2009; Livesey and Cepko, 2001). RPCs produce distinct retinal types in a sequential but overlapping order (‘retinal histogenesis’; Bassett and Wallace, 2012; Cepko, 2014; Marquardt and Gruss, 2002; Fig. 1A). A previous study has proposed that RPCs progress through a series of developmental states (‘competence states’), in which RPCs become responsive to produce particular retinal types (Livesey and Cepko, 2001; Fig. 1A). In past decades, molecular studies have identified many genes that are involved in defining such developmental states (Agathocleous and Harris, 2009; Bassett and Wallace, 2012; Cepko, 2014). Only very recently has systematic identification of the developmental states of the retina in different vertebrate species begun to be carried out, with the advent of single-cell RNA sequencing (scRNA-seq) technology (Clark et al., 2019; Lu et al., 2020).

The zebrafish retina grows beyond embryonic development. Postembryonic retinal stem cells (RSCs), located in the most peripheral region of the ciliary marginal zone (CMZ), are the cell source for this growth (Centanin et al., 2011; Fischer et al., 2013; Marcucci et al., 2016; Marcus et al., 1999; Reinhardt et al., 2015). At the cellular level, a recent study demonstrated that cells of the second or third row at the CMZ tip exhibited stem cell characteristics (Tang et al., 2017). At the molecular level, RSCs specifically expressed rx2 that is controlled by a set of transcription factors including Sox2, Tlx, Gli3 and Her9 in medaka fish (Reinhardt et al., 2015).

Similar to their embryonic counterparts, postembryonic RSCs are multipotent, capable of giving rise to almost all retinal types (Hitchcock et al., 2004; Moshiri et al., 2004; Otteson and Hitchcock, 2003; Perron and Harris, 2000; Reh and Fischer, 2001). Notably, genes of earlier stages of embryonic retinal development are expressed in the more peripheral region of the CMZ, whereas genes of later stages of retinal development are expressed in the CMZ region closer to the neural retina (Dixit et al., 2014; Raymond et al., 2006). At the clonal level, a recent study further showed that RSCs in the CMZ produced the lineages with a similar cell number and cell-type composition to those derived from embryonic RPCs (Wan et al., 2016). All these results suggest that postembryonic RSCs in the CMZ conform to a developmental program that is similar to that used in embryonic RPCs (Harris and Perron, 1998). However, a large-scale forward genetic screen in zebrafish revealed that some genes were required for maintaining RSCs but not for embryonic retinal development, which suggests embryonic RPCs and postembryonic RSCs employ different genetic programs (Wehman et al., 2005). To resolve this inconsistency, we conducted a systematic investigation on the developmental programs that are used by embryonic RPCs and postembryonic RSCs using scRNA-seq (10x Genomics) in zebrafish. Our results showed that embryonic RPCs and postembryonic RSCs mostly conform to the same set of developmental states, with only one major difference in the production of rod photoreceptor cells.

scRNA-seq of embryonic RPCs at different developmental stages

In zebrafish, the retina originates from the optic vesicle, which emerges at ∼12 h postfertilization (hpf) (Wilson and Houart, 2004). At ∼24 hpf, the optic cup forms (Kwan et al., 2012; Picker et al., 2009; Rembold et al., 2006). The production of the first retinal type, RGCs, peaks at ∼36 hpf (Kay et al., 2005; Nawrocki, 1985). At ∼48 hpf, RPCs are mostly undergoing terminal cell divisions (He et al., 2012). All retinal types are fully differentiated by 72 hpf.

Fig. 1.

scRNA-seq of embryonic RPCs at different developmental stages. (A) Schematic showing that the zebrafish retina is composed of structurally and functionally radial retinal segments (dashed sector), which undergo the same developmental program in a sequential order (indicated by the arrows: from the center to the periphery). During the developmental program, RPCs give rise to different retinal types in the stereotyped order (indicated on the right). (B) Schematic workflow for scRNA-seq of staged embryonic RPCs using 10x Genomics technology. (C) t-SNE plots show the clustering of qualified retinal cells at different developmental stages (24, 36, 48 and 72 hpf). Different colors represent distinct cell types. The numbers of cells are indicated. (D) Plot showing proportions of distinct retinal cell types (indicated as different colors) over time, from the data in C. AC, amacrine cell; BC, bipolar cell; CMZ, ciliary marginal zone; HC, horizontal cell; MC, Müller glial cell; PR, photoreceptor cell; RGC, retinal ganglion cell.

Fig. 1.

scRNA-seq of embryonic RPCs at different developmental stages. (A) Schematic showing that the zebrafish retina is composed of structurally and functionally radial retinal segments (dashed sector), which undergo the same developmental program in a sequential order (indicated by the arrows: from the center to the periphery). During the developmental program, RPCs give rise to different retinal types in the stereotyped order (indicated on the right). (B) Schematic workflow for scRNA-seq of staged embryonic RPCs using 10x Genomics technology. (C) t-SNE plots show the clustering of qualified retinal cells at different developmental stages (24, 36, 48 and 72 hpf). Different colors represent distinct cell types. The numbers of cells are indicated. (D) Plot showing proportions of distinct retinal cell types (indicated as different colors) over time, from the data in C. AC, amacrine cell; BC, bipolar cell; CMZ, ciliary marginal zone; HC, horizontal cell; MC, Müller glial cell; PR, photoreceptor cell; RGC, retinal ganglion cell.

To reveal the program of embryonic retinal development, we dissociated the retinas at four developmental time points (24, 36, 48 and 72 hpf) and performed scRNA-seq using the 10x Genomics Chromium 3′ v2 platform (Fig. 1B, Materials and Methods). Libraries of different stages were individually sequenced to a mean depth of 428,141,933, 364,681,265, 185,525,516 and 482,751,657 reads per library, with a median of ∼13,661, 4512, 4793 and 4130 unique molecular identifiers (UMI) and ∼2358, 1214, 1337 and 1250 genes per cell, respectively. We performed raw data processing using Cell Ranger and then removed low-qualified cells and non-retinal cells using the Seurat R package (Fig. S1A; for details see Materials and Methods; http://satijalab.org/seurat/). Consequently, we obtained 4819, 4772, 3587 and 6278 cells for 24-, 36-, 48- and 72-hpf retinas, respectively (Fig. 1C, Fig. S1A). Using high-variance genes, we then clustered retinal cells of different time points individually using t-distributed stochastic neighbor embedding (t-SNE) analysis (Fig. 1C; for details see Materials and Methods). RPCs and retinal neuron types were annotated according to the putative marker genes for individual retinal types (Figs S1C and S2; https://zfin.org). CMZ cells constitute most, if not all, progenitors in the 72-hpf retina (Centanin et al., 2011; Johns, 1977; Straznicky and Gaze, 1971). Thus, proliferative cells of the 72-hpf retina were mostly CMZ progenitor cells (Figs S1C and S2). Biological replicates of 72-hpf samples were largely overlapped in the t-SNE plot (Fig. S1B), indicating data reproducibility. Also, the temporal emergence of distinct retinal types in the single-cell data (from 24 to 72 hpf) was consistent with retinal histogenesis described in earlier studies (Fig. 1D; He et al., 2012).

Defining cellular heterogeneity of embryonic RPCs

The 24-hpf RPCs represent a pool of equivalent progenitors for the subsequent neurogenesis (He et al., 2012). However, using t-SNE analysis, we identified four distinct clusters of 24-hpf RPCs (Cluster 1 to Cluster 4-24 hpf; Fig. 2A, Fig. S3A,D). Each cluster had specific marker genes that exhibited significantly higher expression than others. In our analysis, cluster-specific marker genes were characterized as the ones that have the values of P_val_adj<1e-14 and are not cell cycle-related genes (for details see Materials and Methods).

Fig. 2.

Defining cellular heterogeneity of embryonic RPCs. (A) Visualization of distinct clusters of different-staged embryonic RPCs using t-SNE. Cells are colored by their cluster assignments. Clusters are indicated by their corresponding labels. (B) Dot plots for the expressions of cluster-specific marker genes at different stages. Genes are colored by their cluster assignments. Circle size is proportional to the percentage of cells expressing the marker, and its intensity depicts the average transcript count within expressing cells. (C) Hierarchal clustering of distinct clusters of RPCs from different stages (Jaccard distance metric, average linkage). (D) t-SNE plots show distinct clusters from different embryonic stages after CCA alignment. The left big plot shows the cell distributions of aligned clusters; the right big plot shows the cell distributions of different embryonic stages. The expression patterns of cluster-specific marker genes are indicated in seven small plots. (E) The proportions of distinct clusters across different developmental stages.

Fig. 2.

Defining cellular heterogeneity of embryonic RPCs. (A) Visualization of distinct clusters of different-staged embryonic RPCs using t-SNE. Cells are colored by their cluster assignments. Clusters are indicated by their corresponding labels. (B) Dot plots for the expressions of cluster-specific marker genes at different stages. Genes are colored by their cluster assignments. Circle size is proportional to the percentage of cells expressing the marker, and its intensity depicts the average transcript count within expressing cells. (C) Hierarchal clustering of distinct clusters of RPCs from different stages (Jaccard distance metric, average linkage). (D) t-SNE plots show distinct clusters from different embryonic stages after CCA alignment. The left big plot shows the cell distributions of aligned clusters; the right big plot shows the cell distributions of different embryonic stages. The expression patterns of cluster-specific marker genes are indicated in seven small plots. (E) The proportions of distinct clusters across different developmental stages.

The first distinct cluster (Cluster 1-24 hpf) specifically expressed fabp11a, npm1a, her9 and col15a1b (Fig. 2B, Table S1). This cluster of cells was likely to represent the primitive RPC population. Specifically, her9, the analog of mammalian Hes1, is involved in the maintenance of the stemness of neural stem cells (Reinhardt et al., 2015). Also, col15a1b has been identified as the marker gene for RSCs in zebrafish (Gonzalez-Nunez et al., 2010; Raymond et al., 2006).

The second distinct cluster (Cluster 2-24 hpf) had the characteristic expression of her4.1, her4.2, her4.3, and her4.4, the analogs of mammalian Hes5 genes (Fig. 2B, Table S1). This family of genes is the putative target of Notch signaling that is crucial for RPC fast proliferation and neurogenesis entry (Kageyama and Ohtsuka, 1999; Kageyama et al., 2008, 2015; Kobayashi and Kageyama, 2014). This cluster of cells also expressed fabp7a (also known as blbp), a putative marker for radial glial cells (Fig. 2B; Diotel et al., 2016; Pinto and Gotz, 2007; Podgornyi and Aleksandrova, 2009). Thus, this cluster of cells was likely to represent RPCs which exhibited the properties of fast proliferation as well as radial glial characteristics at the earlier stage of retinal neurogenesis.

The third distinct cluster (Cluster 3-24 hpf) began to express pro-neural factor genes dla, neurod4, and atoh7 (Fig. 2B, Table S1). The dla gene, encoding the ligand for Notch receptors, is associated with the Notch activation that is crucial for retinal differentiation (Kageyama et al., 2009; Tallafuss et al., 2009). Meanwhile, Neurod4 and Atoh7, basic helix loop helix (bHLH) pro-neural transcription factors, have been shown to be widely involved in neuronal differentiation (Brown et al., 2001; Cherry et al., 2011; Chiodini et al., 2013; Murai et al., 2011; Seo et al., 2007). In particular, Atoh7 is one of the first pro-neuronal transcription factors and is essential for the production of the first neurons, RGC, in the zebrafish retina (Kay et al., 2001; Masai et al., 2000). Thus, it suggested that cells of this cluster were undergoing initiation of the neuronal differentiation program (Chiodini et al., 2013; Maurer et al., 2018). This cluster of cells constituted only a small fraction of 24-hpf RPCs (4.7%, Fig. 2E), consistent with the fact that 24-hpf RPCs are mostly undergoing proliferative cell divisions (He et al., 2012). It was notable that Cluster 4-24 hpf had the specific expression of p4hb and ythdf2 (Fig. S3A,D, Table S1), which does not have any known function for retinal development as yet, and so would benefit from future functional study.

Next, we analyzed 36-hpf RPCs. Interestingly, the distinct clusters that we could identify from 36-hpf RPCs were similar to those of 24-hpf RPCs (three clusters; Cluster 1 to Cluster 3-36 hpf; Fig. 2A,B, Fig. S3B,E, and Table S2). Furthermore, the correlation analysis showed that Clusters 1, 2 and 3 of 36-hpf RPCs were equivalent to Clusters 1, 2 and 3 of 24-hpf RPCs, respectively (Fig. 2C). Consistently, canonical correlation analysis (CCA), which identifies gene-gene correlation patterns that are conserved across the datasets (Butler et al., 2018), also showed that clusters of 24-hpf and 36-hpf RPCs could be well correlated and integrated (Fig. 2D). Although three distinct clusters were conserved at both 24 and 36 hpf, the 36-hpf retina had a much higher proportion of neurogenic RPCs (Cluster 3; 34% of total cells; Fig. 2E). This was consistent with the fact that 36-hpf RPCs have entered the neuronal differentiation program.

The 48-hpf retina comprises the most diversity of RPCs because 48-hpf RPCs at the central region mostly undergo terminal cell divisions, whereas the ones at the peripheral region are still at the earlier stages of embryonic retinal development (He et al., 2012; Rapaport et al., 2004; Venters et al., 2015). Indeed, we were able to characterize seven distinct clusters from 48 hpf RPCs (Cluster 1 to Cluster 7-48 hpf; Fig. 2A, Fig. S3C,F). Correlation analysis of clusters showed that the first three distinct clusters of 48-hpf RPCs resembled the corresponding clusters of 24-hpf RPCs (Cluster 1 to Cluster 3-24 hpf) and 36-hpf RPCs (Cluster 1 to Cluster 3-36 hpf) (Fig. 2C,D). In addition, 48-hpf RPCs had four new distinct clusters. Clusters 4 and 5 were related to the development of BCs and HCs because of their specific expression of vsx1 and otx2 (also known as otx2b), and rem1 and onecut1, respectively (Fig. 2A,B, Table S3; Beby and Lamonerie, 2013; Shi et al., 2011; Wang et al., 2014; Wu et al., 2013). Cluster 6 cells specifically expressed otx5 and nr2e3, which are essential for specifying PRs (Haider et al., 2009; Sauka-Spengler et al., 2001; Webber et al., 2008), whereas Cluster 7 cells were characterized by the expression of rlbp1a, which is the putative marker of MCs (Fig. 2A,B, Table S3; Collery et al., 2008). We also analyzed ∼4200 aggregated RPCs of 48-hpf retinas from two biological replicates, and the cells exhibited the same set of distinct clusters (Fig. S4).

Taken together, we characterized seven distinct RPC clusters in zebrafish retina (Fig. 2D). Clusters 1, 2 and 3 likely represent the primitive RPCs, fast-cycling RPCs and earlier neurogenic RPCs, respectively. These three clusters exist throughout the entire course of the retinal development. Meanwhile, Clusters 4, 5, 6 and 7 are likely involved in the generation of BCs, HCs, PRs and MCs, respectively.

Temporal emergence of distinct cell clusters

We next wondered whether distinct RPC clusters represent different developmental states for RPCs or independent sets of RPCs. If it was the former, we would expect that the generation of the seven distinct clusters should exhibit a coherent temporal order. In the developing retina, the differentiation program travels in a wave-like manner, from the ventral-nasal to the dorsal-temporal region (He et al., 2012). In a particular area of the retina, RPCs sequentially express a series of developmental genes over time. In the whole retina, genes of earlier stages were expressed in the more peripheral regions, whereas genes of later stages were expressed in the more central areas (Hu and Easter, 1999; Perron and Harris, 2000). Therefore, the examination of the spatial expression of genes that characterize each distinct cluster using in situ hybridization can indicate the temporal sequence of the generation of distinct clusters.

We then performed in situ hybridization of fabp11a, her4.2, atoh7 and vsx1, which represented Clusters 1 to 4, respectively. fabp11a, a specific gene for Cluster 1 cells, resided at the peripheral retina at 24 hpf. It became restricted into the more peripheral retina at 36 hpf. Eventually, it was expressed at the very tip region of the peripheral retina, where putative retinal stem cells are thought to be, at 48 hpf (Fig. 3A, Fig. S5A; Fischer et al., 2013; Perron and Harris, 2000; Raymond et al., 2006).

Fig. 3.

Temporal emergence of distinct cell clusters. (A) In situ hybridization images of cluster-specific marker genes (Clusters 1 to 4). Genes are colored by their corresponding clusters. The retina boundaries are indicated by the white dashed lines. Yellow arrowheads indicate the areas showing in situ signals (green). The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (B) The pseudo-time trajectories of distinct clusters of 48-hpf RPCs. Cells are colored by their cluster assignments. (C) The expressions of cluster-specific marker genes across pseudo-time. Genes are colored by their cluster assignments. (D) Schematic summary of the evolution of different distinct clusters in a representative coronal retina section. The spatial distribution of individual cluster-specific marker genes suggests the emergence order of distinct developmental states from Clusters 1 to 7.

Fig. 3.

Temporal emergence of distinct cell clusters. (A) In situ hybridization images of cluster-specific marker genes (Clusters 1 to 4). Genes are colored by their corresponding clusters. The retina boundaries are indicated by the white dashed lines. Yellow arrowheads indicate the areas showing in situ signals (green). The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (B) The pseudo-time trajectories of distinct clusters of 48-hpf RPCs. Cells are colored by their cluster assignments. (C) The expressions of cluster-specific marker genes across pseudo-time. Genes are colored by their cluster assignments. (D) Schematic summary of the evolution of different distinct clusters in a representative coronal retina section. The spatial distribution of individual cluster-specific marker genes suggests the emergence order of distinct developmental states from Clusters 1 to 7.

Cluster 2 cells were characterized by her4.2. In the 24-hpf and 36-hpf retinas, her4.2 was expressed in the more central region of the retina than Cluster 1-specific gene fabp11a (Fig. 3A, Fig. S5A). It indicated that Cluster 2 cells emerge later than Cluster 1 cells. At 48 hpf, her4.2 was mostly expressed in the more peripheral retina (Fig. 3A, Fig. S5A).

The expression of atoh7 characterized Cluster 3 cells. It was initiated at 24 hpf and elevated significantly at 36 hpf and onwards (Fig. 3A, Fig. S5A), which was consistent with the scRNA-seq result of a much higher proportion in Cluster 3 of 36-hpf RPCs than that of 24-hpf RPCs (Fig. 2E). At 36 hpf, atoh7 was expressed more centrally than fabp11a and her4.2 (Fig. 3A, Fig. S5A). At 48 hpf, atoh7 was restricted to the peripheral region (Fig. 3A, Fig. S5A). The results suggest that the generation of the Cluster 3 cells occurred later than cells of Clusters 1 and 2.

The gene vsx1 characterized Cluster 4 cells. It was highly expressed at 48 hpf and mainly located in the more central region than fabp11a, her4.2 and atoh7 (Fig. 3A, Fig. S5A). Together, our in situ results indicated that Clusters 1 to 4 occurred sequentially (Fig. 3A,D).

Next, we performed pseudo-time analysis of 48-hpf RPCs, which contained the entire repertoire of distinct clusters, to examine the temporal order of the generation of all distinct clusters. Considering the specific expression of genes (fabp11a and col15a1b) for cellular stemness as well as the most peripheral expression of fabp11a revealed by in situ hybridization (Fig. 3A), we assigned Cluster 1 as the root state (Fig. S5B). The main branch was composed of three clusters of cells in order, starting with Cluster 1, followed by Cluster 2 and finally Cluster 3. Subsequently, the trajectory segregated into two branches: Cluster 3 (atoh7 and neurod4), Cluster 4 (vsx1 and otx2) and Cluster 5 (rem1 and onecout1) on one branch, and Cluster 6 (otx5 and nr2e3) and Cluster 7 (rlbp1a and efhd1) on the other branch (Fig. 3B,C, Fig. S5C).

In summary, the analysis indicated that all distinct cell clusters conformed to a coherent temporal order, thereby representing a series of developmental states of embryonic RPCs (Fig. 3D). In addition, we also characterized developmental states using the combined cells from 24-, 36- and 48-hpf retinas. The results were consistent with the states that we described from staged RPCs above (Fig. S6).

Individual developmental states represent various aspects of retinal development

To further examine these developmental states, we disrupted cluster-specific genes in F0 founders using sets of four of the CRISPR/Cas9 ribonucleoprotein complexes (Cas9 RNP) (Wu et al., 2018). PCR results showed the efficient deletion of genome sequences of targeted genes in the founders (Fig. S7A,B). The sequencing results further showed that the majority of targeted genes (95.1% for bhlhe40, n=41; 97.6% for mab21l2, n=41; 100% for fabp11a, n=18) exhibited various types of mutations, including single-site mutation, multi-site disruption or truncation, and site-spanning deletion, indicating the high efficiency of gene disruption in F0 founders (Fig. S7C-E).

We found that the disruption of Cluster 1-specific genes, including mab21l2, nr2f5, bhlhe40 and fabp11a, led to a small eye phenotype, with little influence on neural differentiation and retinal layer formation at 72 hpf (Fig. S8A,C), consistent with earlier reports (Gath and Gross, 2019; Yamada et al., 2004). However, in some severe cases, the disruption of Cluster 1-specific genes resulted in no eye. For example, in 192 mab21l2-disrupted embryos, 127 were small eye and 8 were eyeless (Fig. S8B). Earlier studies have also shown that downregulation of Cluster 1-specific genes mab21l2 and fabp11a resulted in a severe eyeless phenotype (Qi et al., 2016; Sghari and Gunhaga, 2018). Furthermore, we found that the eyes in the embryos with disruption of Cluster 1-specific genes had reduced their size at 30 hpf (Fig. 4A,B). This was confirmed using the knock-out system of duplex guided RNPs (dgRNPs; Fig. S8D,E; Hoshijima et al., 2019). All these results suggest that Cluster 1-specific genes affect early retinal development. As all these analyses were done in the founders, it is worth verifying the eye defects caused by the disruption of Cluster 1 genes using genetic mutant lines in the future. To examine whether the disruption of Cluster 1-specific genes decreased the proliferation of early RPCs, we performed clonal analysis of single 12-hpf RPCs marked by the photo-conversion of nucleus-localized Kaede from green to magenta. Surprisingly, the number of cells per clone at 30 hpf showed no significant difference between the control groups (blank control and scrambled control) and the groups with the disruption of individual Cluster 1-specific genes (mab21l2, nr2f5, bhlhe40 and fabp11a) (Fig. S8F,G), suggesting that decreased proliferation of early RPCs might not contribute to the small eye phenotype. However, immunostaining of Caspase 3, the marker of apoptotic cells, showed a significant increase in Caspase 3-positive cells (∼35% of total cells) in the 12-hpf retinas with disruption of Cluster 1-specific genes, whereas few apoptotic cells were observed in wild-type retinas (∼5% of total cells; Fig. S8H,I). In the 30-hpf retinas, there were also few apoptotic cells after disruption of Cluster 1-specific genes (Fig. S8J). This indicated that disruption of Cluster 1-specific genes could induce cell death of early RPCs at 12 hpf but not RPCs at 30 hpf. Together, our analyses suggest the involvement of Cluster 1-specific genes in the survival of early RPCs, thereby resulting in the phenotype of small eyes, and even no eyes.

Fig. 4.

Individual developmental states represent various aspects of retinal development. (A) Representative images showing small eyes at 30 hpf as the result of disruption of Cluster 1-specific marker genes (mab21l2, nr2f5, bhlhe40 and fabp11a). Scale bars: 120 μm. (B) Plot showing the proportions of small eyes at 30 hpf as the result of disruption of individual Cluster 1-specific marker genes. ***P<0.001 (Fisher's exact test compared with the blank control). ns, no significant difference. (C) Confocal images of the retinas (coronal view) of SoFa fish lines at 5 dpf after disruption of Cluster 3- (atoh7, foxn4), Cluster 4- (otx2, prox1a) and Cluster 5- (onecut1) specific marker genes. The cell types that decreased significantly in cell number are labeled with stars. Scale bar: 30 μm. (D) Quantitative plots showing the proportions of cell types that decreased significantly in cell number from corresponding gene disruption groups in C (data are mean±s.e.m., data points are shown). ***P<0.001 (Wilcoxon test compared with the control groups). KO, knockout.

Fig. 4.

Individual developmental states represent various aspects of retinal development. (A) Representative images showing small eyes at 30 hpf as the result of disruption of Cluster 1-specific marker genes (mab21l2, nr2f5, bhlhe40 and fabp11a). Scale bars: 120 μm. (B) Plot showing the proportions of small eyes at 30 hpf as the result of disruption of individual Cluster 1-specific marker genes. ***P<0.001 (Fisher's exact test compared with the blank control). ns, no significant difference. (C) Confocal images of the retinas (coronal view) of SoFa fish lines at 5 dpf after disruption of Cluster 3- (atoh7, foxn4), Cluster 4- (otx2, prox1a) and Cluster 5- (onecut1) specific marker genes. The cell types that decreased significantly in cell number are labeled with stars. Scale bar: 30 μm. (D) Quantitative plots showing the proportions of cell types that decreased significantly in cell number from corresponding gene disruption groups in C (data are mean±s.e.m., data points are shown). ***P<0.001 (Wilcoxon test compared with the control groups). KO, knockout.

We further performed the similar gene disruption for atoh7 (Cluster 3), foxn4 (Cluster 3), otx2 (Cluster 4), prox1a (Cluster 4) and onecut1 (Cluster 5) in the SoFa transgenic line, which fluorescently marks different types of retinal neurons (Almeida et al., 2014). The disruption of atoh7 and foxn4 resulted in the significant loss of RGCs and ACs, respectively (Fig. 4C,D, Fig. S8K), consistent with earlier studies (Brown et al., 2001; Li et al., 2004) that suggested Cluster 3 cells were required for the production of early-born retinal neurons. Meanwhile, the disruption of foxn4 could also significantly decrease the HC number, showing its essential role in the fate specification of HCs (Fig. 4C,D, Fig. S8K). The disruption of otx2 significantly reduced the numbers of BCs and PRs (Fig. 4C,D, Fig. S8K), which was consistent with previous findings in otx2 knockout mutants (Koike et al., 2007). Interestingly, although an earlier study in mice showed the requirement of prox1a for HC generation (Dyer et al., 2003), our results instead showed that the disruption of prox1a led to the loss of BCs (Fig. 4C,D, Fig. S8K). Thus, Cluster 4 cells were likely to be involved in late-born retinal cell type production. Furthermore, the disruption of onecut1 resulted in a substantial loss of HCs (Fig. 4C,D, Fig. S8K; Wu et al., 2013), showing the role of Cluster 5 cells in HC production.

Together, our results indicated that individual developmental states reflected different aspects of retinal lineage progression, from the maintenance of retinal progenitors to the generation of different retinal types.

The transition between distinct states

Interestingly, in addition to distinct clusters, we also observed that some cell clusters did express the specific genes of two distinct clusters (Fig. 5A,B). For example, in 36-hpf retinas, one cluster of RPCs (Cluster 1→2) expressed the specific genes of the Clusters 1 and 2, whereas another cluster of RPCs (Cluster 2→3) expressed specific genes of Clusters 2 and 3 (Fig. 5A,B). Gene expression patterns of both cell clusters suggest that they may represent cell states that mediated the transitions of two sequential distinct developmental states.

Fig. 5.

The transition between distinct states. (A) The t-SNE plot of 36-hpf RPCs. Cells are colored by their cluster assignments. Black dashed circles indicate the clusters that exhibit transitional gene expression pattern (termed as ‘transitional cluster’). (B) Heat-map for gene expression of the top six differentially expressed genes for each cluster in 36-hpf RPCs. Clusters 1- to 3-specific marker genes are depicted by pink, green and orange dashed squares, respectively. The genes highly expressed in transitional clusters are highlighted with blue. (C) Biological process enrichment analysis for genes highlighted by blue in B. The x-axis corresponds to the gene ratio quantifying the fraction of genes belonging to a particular biological process. The color code and size correspond to the corrected P-value and number of genes involved in the biological process, respectively. (D) Violin plots showing the expression levels of known S/G2/M markers in distinct and transitional clusters. (E) The proportions of cells within different cell cycle phases in transitional clusters.

Fig. 5.

The transition between distinct states. (A) The t-SNE plot of 36-hpf RPCs. Cells are colored by their cluster assignments. Black dashed circles indicate the clusters that exhibit transitional gene expression pattern (termed as ‘transitional cluster’). (B) Heat-map for gene expression of the top six differentially expressed genes for each cluster in 36-hpf RPCs. Clusters 1- to 3-specific marker genes are depicted by pink, green and orange dashed squares, respectively. The genes highly expressed in transitional clusters are highlighted with blue. (C) Biological process enrichment analysis for genes highlighted by blue in B. The x-axis corresponds to the gene ratio quantifying the fraction of genes belonging to a particular biological process. The color code and size correspond to the corrected P-value and number of genes involved in the biological process, respectively. (D) Violin plots showing the expression levels of known S/G2/M markers in distinct and transitional clusters. (E) The proportions of cells within different cell cycle phases in transitional clusters.

More intriguingly, such clusters, which exhibited transitional gene expression patterns, were characterized by cell cycle S/G2/M genes, such as mcm5, cdc20, aspm and ube2c (Fig. 5C,D). According to the cell cycle scores, which were based on the expression levels of G2/M and S phase marker genes in each cell, we could assign a cell cycle phase to each cell (for details see Materials and Methods). The results also showed that such clusters of cells were mostly in cell cycle stages of S/G2/M (Fig. 5E). These results suggest that the transition of two sequential distinct developmental states may be S/G2/M dependent.

Conserved developmental states between embryonic and postembryonic retinogenesis

To further characterize developmental states of CMZ progenitors derived from postembryonic RSCs in zebrafish (Centanin et al., 2011; Marcus et al., 1999; Tang et al., 2017), we sorted out green cells from 14 day postfertilization (dpf) retinas of Tg(PCNA:GFP) fish (Fig. S9A). Subsequently, we employed scRNA-seq using the 10x Genomics platform. Sorted green cells were mostly proliferating, which mostly represented the CMZ progenitors and derived newborn neurons (Fig. S9A,B). The library was sequenced to a mean depth of 313,577,164 reads per library, with a median of 5427 UMI and 1193 genes per cell. The raw data were processed as described above (for details see Materials and Methods). After excluding lens cells and retinal neurons according to marker genes, we selected 4593 qualified CMZ progenitors for cell cluster identification using t-SNE analysis (Fig. S9C).

We were able to identify seven distinct clusters (Clusters 1- to 7-CMZ), each of which expressed a set of specific genes (Fig. 6A, Fig. S9C,D, and Table S4). The pseudo-time analysis showed that these clusters followed a developmental trajectory. As we described above, we assigned the cluster that expressed fabp11a and col15a1b (Cluster 1-CMZ) as the root state owing to the expression of col15a1b, a putative marker for RSCs in zebrafish (Gonzalez-Nunez et al., 2010; Fig. S9E). Following the Cluster 1-CMZ, cells of Cluster 2-CMZ (expressing genes such as her4.1 and her4.2) and Cluster 3-CMZ (atoh7 and neurod4) occurred sequentially. Then, the trajectory segregated into two branches, Clusters 4- and 5-CMZ (otx2 and vsx1, and rem1 and oneout1, respectively) on one branch and Clusters 6- and 7-CMZ (mafba and pcp4a, and six7 and rcvrn2, respectively) on the other branch (Fig. 6B,C, Fig. S9F). Thus, all distinct cell clusters of CMZ progenitors conformed to a coherent temporal order and represented a series of developmental states of postembryonic RPCs.

Fig. 6.

Conserved developmental states between embryonic and postembryonic retinogenesis. (A) Visualization of distinct clusters of 14-dpf CMZ progenitor cells using t-SNE. Cells are colored by their cluster assignments. Clusters are indicated by their corresponding labels. (B) The pseudo-time trajectories of distinct clusters of 14-dpf CMZ progenitor cells. Cells are colored by their cluster assignments. (C) The expressions of cluster-specific marker genes across pseudo-time. Genes are colored by their cluster assignments. (D) In situ hybridization images of Clusters 1- to 4-specific marker genes (fabp11a, her4.2, atoh7 and vsx1) in the CMZ of the retinas at 14 dpf. Genes are colored by their corresponding clusters. The CMZ boundaries are indicated by the white dashed lines. The areas showing in situ signals (green) are indicated by the region between yellow dashed lines. The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (E) Schematic showing the spatial distribution of cells of Clusters 1 to 4 in the CMZ of the retina at 14 dpf. (F) Hierarchal clustering of distinct clusters from embryonic (48-hpf RPCs) and postembryonic (14-dpf CMZ progenitor cells) stages (Jaccard distance metric, average linkage). Yellow arrow indicates the cluster with specific expression of mafba. (G) In situ hybridization images of mafba in 48-hpf RPCs and the CMZ progenitor cells of 14-dpf retina. The 48-hpf retina and the CMZ of 14-dpf retina are indicated by the white dashed lines. Yellow arrowheads indicate the areas showing in situ signals (green). Scale bars: 30 μm.

Fig. 6.

Conserved developmental states between embryonic and postembryonic retinogenesis. (A) Visualization of distinct clusters of 14-dpf CMZ progenitor cells using t-SNE. Cells are colored by their cluster assignments. Clusters are indicated by their corresponding labels. (B) The pseudo-time trajectories of distinct clusters of 14-dpf CMZ progenitor cells. Cells are colored by their cluster assignments. (C) The expressions of cluster-specific marker genes across pseudo-time. Genes are colored by their cluster assignments. (D) In situ hybridization images of Clusters 1- to 4-specific marker genes (fabp11a, her4.2, atoh7 and vsx1) in the CMZ of the retinas at 14 dpf. Genes are colored by their corresponding clusters. The CMZ boundaries are indicated by the white dashed lines. The areas showing in situ signals (green) are indicated by the region between yellow dashed lines. The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (E) Schematic showing the spatial distribution of cells of Clusters 1 to 4 in the CMZ of the retina at 14 dpf. (F) Hierarchal clustering of distinct clusters from embryonic (48-hpf RPCs) and postembryonic (14-dpf CMZ progenitor cells) stages (Jaccard distance metric, average linkage). Yellow arrow indicates the cluster with specific expression of mafba. (G) In situ hybridization images of mafba in 48-hpf RPCs and the CMZ progenitor cells of 14-dpf retina. The 48-hpf retina and the CMZ of 14-dpf retina are indicated by the white dashed lines. Yellow arrowheads indicate the areas showing in situ signals (green). Scale bars: 30 μm.

Furthermore, we confirmed some distinct clusters using in situ hybridization of fabp11a (Cluster 1-CMZ), her4.2 (Cluster 2-CMZ), atoh7 (Cluster 3-CMZ) and vsx1 (Cluster 4-CMZ) (Fig. 6D, Fig. S10A). They were expressed in the CMZ with a characteristic spatial order: fabp11a was expressed in the most peripheral region, followed by her4.2 (Cluster 2-CMZ), atoh7 (Cluster 3-CMZ) and vsx1 (Cluster 4-CMZ), which was in the most central area (Fig. 6D,E).

To examine whether embryonic RPCs and postembryonic RSCs shared a similar developmental program, we performed the correlation analysis of developmental states that present in CMZ progenitors and 48-hpf RPCs. The results showed that developmental states from CMZ progenitors were highly correlated with their counterparts of embryonic RPCs (Fig. 6F). However, we also observed some differences: CMZ progenitors, but not embryonic RPCs, had the cell cluster with high expression of mafba (Cluster 6-CMZ) (Fig. 6F). The expression of mafba was mainly in the CMZ region close to the PR layer at 14 dpf (Fig. 6G, Fig. S10B). Moreover, previous results have shown that chick Mafa drives rod photoreceptor cell formation (Kim et al., 2016). The evidence suggests that cells of Cluster 6-CMZ may play an essential role in rod photoreceptor cell genesis.

Embryonic RPCs and postembryonic RSCs share similar expression dynamics of hes family genes

The hairy-related hes/her genes, encoding a set of DNA-binding transcription repressors belonging to basic helix-loop-helix (bHLH) transcription factor family, are involved in many aspects of neural stem/progenitor cells including self-renewal, proliferation, and differentiation (Kageyama et al., 2008, 2015). Therefore, we wondered whether CMZ progenitors and embryonic RPCs exhibited similar expression dynamics of hes genes. We first examined their expression dynamics in embryonic RPCs of 48-hpf retinas. Gene-to-gene correlation analysis showed that all expressed hes/her genes were clustered into 3 groups (for details see Materials and Methods): hes1 group (including her6 and her9, analogs of mammalian Hes1), hes5 group (including her4.1, her4.2, her4.3 and her4.4, analogs of mammalian Hes5) and hes6 group (including hes6 and hes2.2) (Fig. 7A). Genes of these three groups exhibited distinct but overlapping expression patterns. The expression levels of hes1 group genes were the most abundant in Clusters 1 and 2. The expression of hes5 group genes peaked at Cluster 2 and diminished from Cluster 3. Moreover, hes6 group genes were mostly expressed in Clusters 3 and 4 (Fig. 7B).

Fig. 7.

Embryonic RPCs and postembryonic RSCs share similar expression dynamics ofhesfamily genes. (A) Gene-to-gene correlation of different hes/her genes in 48-hpf RPCs. Yellow squares depict different hes clusters. (B) Heat-map showing the distribution of hes/her genes across different distinct clusters of 48-hpf RPCs. (C) In situ hybridization images of her9 (hes1), her4.2 (hes5) and hes6 in the coronal view of the 48-hpf retinas. The retina boundaries are indicated by the white dashed lines. The areas showing in situ signals (green) are indicated by the region between yellow dashed lines. The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (D) Schematic showing the evolution of her9 (hes1), her4.2 (hes5) and hes6 across distinct clusters at a coronal view of the embryonic retina. The emergence order of hes/her family genes is her9 (hes1), her4.2 (hes5) and hes6. (E) Gene-to-gene correlation of different hes/her genes in 14-dpf CMZ progenitor cells. Yellow squares depict different hes clusters. (F) Heat-map showing the distribution of hes/her genes across different distinct clusters of 14-dpf CMZ progenitor cells. (G) In situ hybridization images of her6 (hes1), her4.2 (hes5) and hes6 in the coronal view of the 14-dpf CMZ. The retina boundaries are indicated by the white dashed lines. The areas showing in situ signals (green) are indicated by the region between yellow dashed lines. The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (H) Schematic showing the evolution of her6 (hes1), her4.2 (hes5) and hes6 across distinct clusters at a coronal view of the CMZ. The emergence order of hes/her family genes is her6 (hes1), her4.2 (hes5) and hes6.

Fig. 7.

Embryonic RPCs and postembryonic RSCs share similar expression dynamics ofhesfamily genes. (A) Gene-to-gene correlation of different hes/her genes in 48-hpf RPCs. Yellow squares depict different hes clusters. (B) Heat-map showing the distribution of hes/her genes across different distinct clusters of 48-hpf RPCs. (C) In situ hybridization images of her9 (hes1), her4.2 (hes5) and hes6 in the coronal view of the 48-hpf retinas. The retina boundaries are indicated by the white dashed lines. The areas showing in situ signals (green) are indicated by the region between yellow dashed lines. The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (D) Schematic showing the evolution of her9 (hes1), her4.2 (hes5) and hes6 across distinct clusters at a coronal view of the embryonic retina. The emergence order of hes/her family genes is her9 (hes1), her4.2 (hes5) and hes6. (E) Gene-to-gene correlation of different hes/her genes in 14-dpf CMZ progenitor cells. Yellow squares depict different hes clusters. (F) Heat-map showing the distribution of hes/her genes across different distinct clusters of 14-dpf CMZ progenitor cells. (G) In situ hybridization images of her6 (hes1), her4.2 (hes5) and hes6 in the coronal view of the 14-dpf CMZ. The retina boundaries are indicated by the white dashed lines. The areas showing in situ signals (green) are indicated by the region between yellow dashed lines. The nuclei are labeled by DAPI (magenta). Scale bar: 30 μm. (H) Schematic showing the evolution of her6 (hes1), her4.2 (hes5) and hes6 across distinct clusters at a coronal view of the CMZ. The emergence order of hes/her family genes is her6 (hes1), her4.2 (hes5) and hes6.

To confirm this, we examined the gene expressions of hes1 (her9), hes5 (her4.2) and hes6 in the coronal sections of retinas at 48 hpf using in situ hybridization (Fig. 7C). Consistent with the scRNA-seq results, in situ hybridization also showed the characteristic spatial order of these three hes gene groups. At 48 hpf, hes1 expression was always at the most peripheral region, where Cluster 1 cells resided. Next to that region, hes5 was expressed at the region where Cluster 2 cells were located. hes6 was found in the most central area where cells of Clusters 3 and 4 supposedly resided (Fig. 7C,D, Fig. S10C).

Furthermore, gene-to-gene correlation analysis showed that all expressed hes/her genes were also clustered into three groups in CMZ progenitors (Fig. 7E). Remarkably, the temporal expression patterns of these three groups of hes genes in CMZ progenitors were mostly similar to those in embryonic RPCs (Fig. 7F). Spatial expression patterns of hes1 (her6), hes5 (her4.2) and hes6 in the CMZ viewed by in situ hybridization were also similar to their expression patterns during embryonic retinogenesis. That is, hes1 was expressed in the leading edge of CMZ, hes6 was expressed in the most central region of the CMZ, whereas hes5 was expressed in the area between these two regions (Fig. 7G,H, Fig. S10D). Together, the temporal expressions of hes family genes were also largely conserved by embryonic RPCs and postembryonic RSCs.

Distinct fate outputs of embryonic RPCs and postembryonic RSCs

Finally, we wondered whether there was any difference in terms of cell fate output between embryonic RPCs and postembryonic RSCs. To this end, we first clustered the differentiated cells of the 72-hpf retina (1772 cells), which were derived from embryonic RPCs, as well as the cells of 14 dpf CMZ (817 cells), which were derived from postembryonic RSCs. The result showed that cell clusters derived from embryonic RPCs and postembryonic RSCs were mostly similar Fig. S11A-C. However, there was only one difference in a cluster of rod photoreceptor cells. Rods were significantly more abundant in differentiated cells derived from postembryonic RSCs than those from embryonic RPCs (Fig. 8A-C).

Fig. 8.

Distinct fate outputs of embryonic RPCs and postembryonic RSCs. (A) Cell clusters of PR subtypes that were obtained from the pool of single cells of 72-hpf retinas and GFP-positive cells of 14-dpf retinas of Tg(PCNA:GFP). Cells are colored by their cluster assignments. Rod cluster is enlarged on the right, and it is shown by stages of 72 hpf (gray) and 14 dpf (black). (B) Heat-map of gene expression showing the top 10 differentially expressed genes for each distinct cell cluster of PRs. Specific marker genes of cone and rod are depicted by green and orange, respectively. (C) Bar plot showing the ratios of rods in total PRs of 72 hpf and 14 dpf. (D) Image of EdU staining of the 20-dpf retina of the Tg(xops:nfsB-mCherry) fish after a pulse of EdU at 5 dpf. A schematic of plasmid design for Tg(xops:nfsB-mCherry) is provided above the image. Scale bar: 30 μm. (E) Bar plot showing the density of rods from embryonic (EM) and postembryonic (post EM) retinogenesis. ***P<0.001 (Wilcoxon test). Data are mean±s.e.m. n=86 sections from six retinas.

Fig. 8.

Distinct fate outputs of embryonic RPCs and postembryonic RSCs. (A) Cell clusters of PR subtypes that were obtained from the pool of single cells of 72-hpf retinas and GFP-positive cells of 14-dpf retinas of Tg(PCNA:GFP). Cells are colored by their cluster assignments. Rod cluster is enlarged on the right, and it is shown by stages of 72 hpf (gray) and 14 dpf (black). (B) Heat-map of gene expression showing the top 10 differentially expressed genes for each distinct cell cluster of PRs. Specific marker genes of cone and rod are depicted by green and orange, respectively. (C) Bar plot showing the ratios of rods in total PRs of 72 hpf and 14 dpf. (D) Image of EdU staining of the 20-dpf retina of the Tg(xops:nfsB-mCherry) fish after a pulse of EdU at 5 dpf. A schematic of plasmid design for Tg(xops:nfsB-mCherry) is provided above the image. Scale bar: 30 μm. (E) Bar plot showing the density of rods from embryonic (EM) and postembryonic (post EM) retinogenesis. ***P<0.001 (Wilcoxon test). Data are mean±s.e.m. n=86 sections from six retinas.

To experimentally test whether postembryonic RSCs could generate more rods than embryonic RPCs, we created a transgenic line Tg(xops:nfsB-mCherry). In it, nfsB-mCherry is driven by the promoter of xops that was specifically expressed in rod photoreceptors (Morris et al., 2005). Moreover, the expression pattern of rods labeled by this transgenic line was consistent with the expression of endogenous rho, a putative gene marker for rods (Fig. S11D). To define the peripheral region of the retina derived from postembryonic RSCs, a pulse of EdU to Tg(xops:nfsB-mCherry) fish at 5 dpf resulted in the clones of EdU-positive cells at 20 dpf. Thus, cells in the peripheral region of the retina between the EdU-positive clone and the CMZ were derived from postembryonic RSCs, whereas cells in the central region of the retina between the two EdU-positive clones were largely derived from embryonic RPCs (Fig. 8D). Quantitative analysis showed that the density of rod photoreceptors in the peripheral region (24 cells/100 μm) was significantly higher than that in the central region (12 cells/100 μm; Fig. 8E). It suggests the preferential production of rod photoreceptor cells by postembryonic RSCs.

Using scRNA-seq, we characterized a repertoire of seven distinct states for embryonic RPCs in zebrafish. Each state expresses a set of specific genes. Gene disruption assays suggest that distinct states represent the RPCs at different stages of retinal lineage progression. Specifically, we identified primitive RPCs (Cluster 1), fast-proliferating RPCs (Cluster 2), earlier neurogenic RPCs (related to the production of RGCs and ACs; Cluster 3), later neurogenic RPCs (related to the generation of BCs and PRs; Cluster 4) and fate-committed precursors for HCs (Cluster 5), PRs and MCs (Clusters 6 and 7).

The developmental programs of embryonic and postembryonic retinogenesis in zebrafish

The previous in situ studies of different developmental genes suggest that embryonic and postembryonic retinogenesis may use a similar developmental program (Harris and Perron, 1998; Perron et al., 1998). Using scRNA-seq, we systematically characterized the developmental states of embryonic RPCs and postembryonic RSCs. Our results also indicate that embryonic and postembryonic retinogenesis share a similar developmental program (Figs 3 and 6). However, an earlier study using a large-scale forward genetic screen identified a list of genes that are required for the regulation of postembryonic RSCs. Disruption of these genes affects the postembryonic retina growth from the CMZ, but retains healthy development of the embryonic retina. It suggests that embryonic RPCs and postembryonic RSCs may utilize different genetic programs (Wehman et al., 2005). However, this may be interpreted by the possibility that the maternal deposition of molecules could allow normal embryogenesis (Abrams and Mullins, 2009; Gross et al., 2005; Pelegri, 2003). We defined developmental states using cluster-specific marker genes. The disruption of these marker genes resulted in various defects in cell-fate acquisition. Thus, the clusters derived from our analysis provide a better representation of retinal development.

scRNA-seq analysis of RPCs across species

Most recently, similar single-cell analyses of the developing retinas has also been conducted in mice and humans (Clark et al., 2019; Lu et al., 2020). These scRNA-seq studies provide new insights into RPC heterogeneity, as well as raise technical concerns: (1) Many clusters of embryonic RPCs are conserved across species, such as fast-proliferating RPCs, later neurogenic RPCs, and fate-committed precursors. It suggests that the retinas of vertebrate species share major neurogenesis stages. (2) We characterize the developmental state representing the primitive RPCs. Furthermore, these genes are also shared mainly by postembryonic RSCs in the CMZ. Lack of these genes, such as fabp11a and mab21l2, can result in small eyes, and an even more dramatic eyeless phenotype (Fig. 4A, Fig. S8A,B). Surprisingly, these genes do not affect the proliferation significantly but maintain the number of progenitor cells in the optic vesicles via apotosis. To our knowledge, it is the first time that such an RPC cluster is identified. (3) In light of lineage relationships of different retinal fates revealed by earlier lineage tracing and clonal studies in zebrafish (He et al., 2012), the identification of lineage relationships among differentiated retinal types using pseudo-time analysis is sometimes problematic. For example, the pseudo-time analysis showed that ACs and BCs were on one branch, whereas PRs were on the other branch (Fig. 3B, Fig. S5C). However, earlier clonal studies have clearly shown that ACs and PRs are closer in terms of lineage relationships (He et al., 2012). This suggests that the similarity of the transcriptional profile of two cells is not always a robust predictor of their lineage relationship, particularly for differentiated cell types. In other words, because two types of neurons have similar transcriptional profiles it does not mean that they are closer in lineage origin.

The progression of developmental states

Identification of seven developmental states led to an immediate question as to the mechanisms driving state transitions. There is some evidence that sheds light on state transition. We observed the clusters that express specific genes of the leading and following distinct developmental states (Fig. 5A,B). This strongly suggests that such states mediate the transition. Also, we observed that such clusters expressed cell cycle S/G2/M genes (Fig. 5C,D), suggesting these transitions may be cell cycle dependent. Furthermore, we showed the temporal expression of hes/her family genes (Fig. 7). This suggests a possible role of Notch signaling in driving the state transition.

Our findings raise another interesting question as to how individual RPCs progress through these developmental states to generate clones that vary significantly in cell numbers and cell-type compositions (Holt et al., 1988; Turner and Cepko, 1987; Turner et al., 1990; Wetts and Fraser, 1988). First, individual RPCs do progress through all states, but their staying in each state is highly variable. If so, it is crucial to figure out the relative contributions of intrinsic and environmental control of state transition. For instance, individual RPCs could either count on an intrinsic clock, such as rounds of cell divisions, or environmental cues to decide the entry and departure of a given developmental state. Second, individual RPCs enter a given developmental state stochastically, thereby resulting in final clonal variability. Earlier clonal analysis in zebrafish has shown that clonal variability could be well explained by a model involving independent firing of pro-neuronal transcription factors (Boije et al., 2015). Most recently, single-cell transcriptome analysis of RPCs in mouse developing retina also suggests a system of stochastic fate choice (Clark et al., 2019).

Abundant rod production in the peripheral retina

Whether embryonic RPCs and postembryonic RSCs use a similar developmental program still lacks consensus (Harris and Perron, 1998; Perron et al., 1998; Wehman et al., 2005). Our analysis provides systematic evidence of single-cell transcriptomics to support a unifying program used by embryonic RPCs and postembryonic RSCs. However, postembryonic RSCs produce many more rod PRs than their embryonic counterparts. Currently, it is thought that rod photoreceptor cells derive from MCs instead of CMZ cells (Lenkowski and Raymond, 2014; Raymond et al., 2006); however, to our knowledge, this lacks strong experimental evidence. In the future, lineage analyses of single postembryonic RSCs and their derived MCs are needed to solve this issue.

Zebrafish

Zebrafish were raised and bred under standard conditions at 28°C. The embryos were staged using hpf, as previously described (Kimmel et al., 1995). The following zebrafish strains were used: wild-type AB strains, SoFa lines Tg(Atoh7:gapRFP::Ptf1a:GFP::Crx:CFPcaax) (ZFIN ID: ZDB-FISH-150901-13478), Tg(PCNA:GFP), Tg(xops:nfsB-mCherry). Tg(xops:nfsB-mCherry) transgenic lines were produced using established methods (Stuart et al., 1990). The other fish lines were gifted from William A. Harris' lab, University of Cambridge, UK. All animal procedures performed in this study were approved by the Animal Use Committee of Institute of Neuroscience, Chinese Academy of Sciences (NA-045-2019).

Sample preparation for scRNA-seq

The embryonic retinas of wild-type zebrafish at 24, 36, 48 and 72 hpf were dissected out according to the earlier literature (Lopez-Ramirez et al., 2016) with tungsten needles under the stereomicroscope. The pigment cells and lens cells from the retinas of 24-hpf and 36-hpf embryos were not removed as they were tightly connected with retinas. However, most of the pigment epithelium and lens were removed from the retinas of 48-hpf and 72-hpf embryos. The dissected retinas were then dissociated using 200 μl papain solution [4 μl papain (Worthington), 4 μl DNase I (1%; Sigma-Aldrich) and 8 μl L-cysteine (12 mg/ml, Sigma-Aldrich) into 184 μl of DMEM/F12 (Invitrogen)] at 37°C for 15 min. Within these 15 min, cells were further dissociated via gentle up-and-down pipetting every 2-3 min. The digestion was terminated using 800 μl washing buffer [100 ml of washing buffer was prepared by adding 650 μl glucose (45%; Invitrogen), 500 μl HEPES (1M; Sigma-Aldrich) and 5 ml fetal bovine serum (Gibco) into 93.85 ml DPBS (1×; Invitrogen)]. The single-cell suspension was filtered using a 40 μm cell strainer (BD Falcon) and then centrifuged at 500 g for 5 min at 4°C. The supernatant was discarded and pellets resuspended with 1× PBS with 0.04% bovine serum albumin (BSA).

The postembryonic retinas of Tg(PCNA:GFP) at 14 dpf were dissected as above, and the lens and pigment epithelium were removed as cleanly as possible. The single-cell suspension was prepared as described above. Subsequently, the GFP- positive cells were sorted into the solution of 1× PBS with 0.04% BSA on MoFlo XDP (Beckman Coulter).

Single-cell library preparation and sequencing

Approximately 5000-10,000 cells were loaded onto the Chromium Single Cell Chip (10x Genomics) according to the manufacturer's protocol. In general, the live-cell percentage was over 85%, confirmed by Trypan Blue staining. scRNA-seq libraries were generated using the GemCode Single-Cell Instrument and Single Cell 3′ Library & Gel Bead kit v2 Chip kit (10x Genomics) according to the manufacturer's protocol. Library quantification and quality assessment were performed by Qubit fluorometric assay (Invitrogen) with a dsDNA High Sensitivity Assay Kit (Thermo Fisher Scientific) and the fragment analysis with a High Sensitivity Large Fragment 50kb Analysis Kit (Advanced Analytical Technologies; AATI) separately. Indexed libraries were sequenced using the Illumina NovaSeq 6000 sequencer with the S2 flow cell using paired-end 150×150 bp as the sequencing mode. The sequencing depth was ∼50-200 reads per cell.

scRNA-seq data processing

After sequencing, data were analyzed using Cell Ranger Single Cell Software Suite (v2.1.0) provided on the 10x Genomics website (https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/what-is-cell-ranger) and several R pipelines [e.g. Seurat (Butler et al., 2018; Stuart et al., 2019) and Monocle (Qiu et al., 2017a,b; Trapnell et al., 2014)].

Cell Ranger count aligned single-cell FASTQ files (Novogene) to reference zebrafish genomes (Zv10) to generate single-cell raw count matrices. On average, 65,787, 59,076, 51,235, 147,511 and 59,031 raw reads per cell were obtained from 24-, 36-, 48-, 72-hpf and 14-dpf single-cell samples, respectively. And approximately 6508, 6173, 3621, 2826 and 5312 cells were recovered from 24-, 36-, 48-, 72-hpf and 14-dpf samples, respectively. For replicates of 48-hpf, 72-hpf and 14-dpf retinas, 37,507, 57,498 and 201,343 raw reads per cell were obtained and 5595, 4056 and 1700 cells were recovered. Two 48-hpf and 72-hpf replicates were aggregated individually using Cell Ranger aggr script with the default normalization parameters.

All of the raw count matrices were imported into R (v3.4.3/v3.5.1) and analyzed with the Seurat 2.3.1 package (https://satijalab.org/seurat/) except for the raw count matrices of aggregated 48-hpf replicates, which was analyzed with the Seurat 3.1 package (https://satijalab.org/seurat/).

Initially, low-quality cells or empty droplets [detected in fewer than 200 genes (features)] and low-abundance genes (detected in fewer than three cells) were filtered before creating the Seurat objects. Later, cell doublets and poor libraries [with very high genes (features) and high mitochondrial percentage] were filtered again after the Seurat objects created. After that, the filtered data were log normalized with a default size factor of 10,000 and scaled while regressing out the differences between the G2/M and S phase scores. Highly variable genes were identified based on their expression and dispersion using the LogVMR dispersion function in Seurat. Highly variable genes fit the following criteria: 0.0125<average expression<3 and s.d.>0.5. Dimensionality reduction was then performed using principal component analysis (PCA). Significant principal components (PCs; P-values<1e-4) were used for clustering based on a shared nearest-neighbor modularity optimization-based clustering algorithm with default parameters and resolution of 0.6. The clustering results were visualized in reduced dimensions via t-SNE embedding using the same PCs. For aggregated 48-hpf replicates, non-linear dimensionality reduction was performed and visualized using Uniform Manifold Approximation and Projection (UMAP) using the Seurat 3.1 package.

Following that, the Seurat function ‘FindAllMarkers’ was used to find genes with significant differential expression (cluster-specific marker genes; P_val_adj<1e-14) of each cluster with default parameters. Clusters were annotated based on the expression datasets of these cluster-specific marker genes on ZFIN (https://zfin.org).

As retinal cells could be contaminated with cells from neighboring tissues during single-cell suspension preparation, these contaminated cell clusters were identified and removed after the initial clustering. Following these step-wise processes, we obtained the purified data of each sample (Fig. S1A).

Embryonic developmental state analysis

Progenitor clusters of individual time points

The purified data were analyzed again following the same procedure described above using the Seurat package (2.3.1). Consequently, different clusters with cluster-specific marker genes (P_val_adj<1e-14) were obtained. According to these cluster-specific marker genes, different clusters were annotated into RPCs and different retinal neurons (Fig. S1C).

Following that, the RPCs of 24-, 36- and 48-hpf retinal samples were selected for further analysis via Seurat (2.3.1). The pipeline for analysis of developmental states is described in Fig. S12A. For 24-hpf RPCs, the clusters were largely biased with retinal patterning [e.g. ventral (aldh1a3), dorsal (aldh1a2); Fig. S12B-D; Duester, 2009; Peters, 2002; Picker et al., 2009]. Therefore, significant PCs associated with retinal patterning (Fig. S12E) were excluded and the remaining significant PCs were used to generate final Seurat clusters with t-SNE plots for visualization. But there was no patterning-related PCs in 36- and 48-hpf progenitors. All clusters were generated based on the default parameters. Cluster-specific marker genes (P_val_adj<1e-14) of each cluster were found via ‘FindAllMarkers’ function under different resolutions. Individual distinct clusters were identified based on the cluster-specific marker genes. Among these distinct clusters, some clusters with cluster-specific marker genes of unknown functions were classified as unknown clusters and excluded from downstream analysis. Meanwhile, some clusters with only cell cycle-related genes were also excluded from downstream analysis. Also, there were some clusters with fewer cluster-specific marker genes but with cluster-specific marker genes of two distinct clusters. Clustering analysis results were different under the resolution of 0.4, 0.6 and 0.8. At the lower resolution of 0.4, some distinct clusters did not appear but merged with other clusters. Whereas at the higher resolution of 0.6 and 0.8, these distinct clusters were able to be identified and the number remains constant (Fig. S13). For this reason, we decided to perform our analysis with a resolution of 0.6. The lists of distinct cluster-specific marker genes of 24-, 36- and 48-hpf RPCs are provided in Tables S1-S3, respectively.

Pseudo-time analysis of distinct clusters

The Monocle 2 R package was used to perform pseudo-time analysis. Distinct clusters were picked up via Seurat ‘SubsetData’ function (2.3.1). Following that, these subset Seurat objects were imported into Monocle for cell trajectory analysis via Monocle ‘importCDS’ function. Cell ordering was performed by using the top 800 differentially expressed genes that were determined by the ‘differentialGeneTest’ function. ‘DDRTree’ was applied to reduce dimensional spaces with default parameters, and the ‘plot_cell_trajectory’ function was used for visualization of the cell trajectory. The clusters with specific expression of fabp11a or npm1a were identified as the root states based on the expression of several stem cell markers [e.g. her9 (Reinhardt et al., 2015) and col15a1b (Gonzalez-Nunez et al., 2010)].

Comparison and integration of distinct clusters of different embryonic stages

The ClusterMap R package (Gao et al., 2019) was used to observe the correlation of distinct clusters identified from different embryonic stages (24, 36 and 48 hpf). Cluster-specific marker genes (P_val_adj<1e-14) of different distinct clusters from different stages were generated via ‘FindAllMarkers’ from the Seurat package (2.3.1). Using the hierarchical clustering of these marker genes of each distinct cluster, ClusterMap matched these clusters across different time points via default parameters.

Seurat objects of the identified distinct clusters from different embryonic stages were integrated using CCA, the alignment method described previously (Butler et al., 2018). Specifically, it learned a shared gene correlation structure that is conserved among different datasets. The top 1000 variable genes were identified in each object, and those that were present in all objects (2089 genes) were used. Then CCA was performed using these selected genes and 17 dimensions for downstream analysis. Next, the combined object was clustered using the similar procedure described above using the Seurat package (2.3.1).

Combine cells from different embryonic stages

In parallel with the pipeline described above to show the embryonic developmental states, another strategy was applied. First, Cell Ranger aggr was used to combine data from 24-, 36- and 48-hpf libraries to generate a combined dataset with the normalized sequencing depth. Second, the raw count matrix of the combined data was analyzed using the Seurat package (2.3.1). The procedure was mostly similar to that described above. Third, after the cell clustering, the clusters of RPCs were selected to perform further cell clustering (Fig. S6A).

To eliminate the effect of retinal patterning, significant PCs (P-values<1e-4) without retinal patterning genes were used for generating final Seurat clusters, which were visualized using the t-SNE plot. Next, using cluster-specific marker genes of these final clusters, the correlation of these clusters with previously defined states were performed via the ClusterMap R package. Finally, these clusters were classified into appropriate states. After these analyses, these identified clusters were imported into Monocle object to carry out pseudo-time trajectory analysis, as described above.

Progenitor clusters of 48 hpf

To enrich the progenitor numbers analyzed at 48 hpf, we used the aggregated 48-hpf dataset from two replicates. The purified dataset was analyzed following the same procedure described above using the Seurat package (3.1). According to cluster-specific marker genes (P_val_adj<1e-14), RPCs of 48-hpf retinal cells were selected for further analysis via Seurat (3.1). The analyzed pipeline was similar to the procedure described above. A subset of highly variable features that exhibit high cell-to-cell variation was calculated with default parameters. Scaling was performed to ‘regress out’ heterogeneity associated with the cell cycle stages based on the scores of difference between the G2/M and S phases before the PCA. PCA was performed on variable features to reduce the dimensionality of the dataset. All clusters were established based on the default parameters at a resolution of 0.4. Cluster-specific marker genes of each cluster were found via the ‘FindAllMarkers’ function. Individual distinct clusters were identified based on the cluster-specific marker genes.

Postembryonic developmental state analysis

The purified data of 14-dpf cells were analyzed again following the procedure described above. After progenitors were selected, significant PCs (P-values<1e-4) were used for clustering. Distinct clusters were identified according to the criteria used in embryonic RPCs above. The distinct cluster-specific marker genes of 14-dpf CMZ progenitors are listed in Table S4. Subsequently, distinct clusters were imported into Monocle object for generating the developmental trajectories. The pseudo-time analysis was performed by using the top 680 differentially expressed genes, which are determined by the ‘differentialGeneTest’ function. The ‘plot_cell_trajectory’ function achieves visualization of the cell trajectory. As for the embryonic retinogenesis, the cluster with specific expression of fabp11a or col15a1b was identified as the root state. To compare the developmental states of embryonic and postembryonic stages, ClusterMap was used to generate the hierarchical clustering of 48-hpf and 14-dpf developmental states.

Cell cycle phase analysis

The cell cycle phase of each cell was assigned based on the cell cycle scores using Seurat R package (https://satijalab.org/seurat/v2.4/cell_cycle_vignette.html). In principle, the score was calculated according to the expression of G2/M and S phase markers. Cells expressing neither G2/M nor S phase markers were likely not cycling and in G1 phase. Along with these scores, each cell was predicted into G2/M, S or G1 phase.

Fate output analysis between embryonic and postembryonic retinogenesis

The combined Seurat object containing data of 72-hpf and the replicated 14-dpf was created via the ‘merge’ function of the Seurat Package (2.3.1). The following procedure in Seurat was mostly similar to that described above. After clustering, the clusters of PRs were selected to perform further analysis. All clusters were visualized by t-SNE based on the default parameters at the resolution of 0.6. Annotations of different clusters were based on their cluster-specific marker genes.

Gene correlation analysis

Highly correlated gene clusters were obtained by measuring the pairwise Pearson correlation distances. These correlation distances were calculated using the ‘bioDist’ R package.

GO analysis

Gene Ontology (GO) analysis of highly differentially expressed genes of transition clusters at 36 hpf was performed by the ‘clusterProfiler’ R package (Yu et al., 2012). Biological process (BP) was the chosen ontology for analysis and presentation.

Gene disruption via the CRISPR/Cas9 system

We used the CRISPR/Cas9 system as described previously to perform the specific and rapid gene disruption (Wu et al., 2018). This CRISPR/Cas9 genetic screening could produce gene disruptions in zebrafish with optimized, redundant gene targeting [mostly four single guide (sg)RNAs per gene]. The targeted coding sequences were divided into five equal sections, and the four sgRNAs were designed to target the first four sections individually (Fig. S7B). The sgRNA sequences were achieved from the candidate list from an earlier study (Wu et al., 2018). PCR was used to generate templates for sgRNA mixture production from the pSQT1313 (Addgene #53370) plasmid with a forward primer containing the T7 promoter and the guide sequence (forward primes are listed in Table S5), and a reverse primer containing the standard chimeric sgRNA scaffold (reverse primer: AAAAAAAGCACCGACTCGGTGCCAC). Then we used these templates for sgRNAs synthesis via T7 MEGAscript (Invitrogen). The concentration of sgRNA mixtures was kept at ∼500 ng/μl. The injection mixture solution contained 1 μl cas9 protein (final concentration of 400 ng/μl) and 1.5 μl sgRNA mixture (final concentration of 300 ng/μl). Microinjection was performed by injecting 1 nl of the mixture into one-cell stage embryos. As a negative control, 1 nl of mixture of scrambled sgRNA and Cas9 protein was injected. sgRNA mixtures of Cluster 1-specific marker genes were injected into wild-type embryos, and sgRNA mixtures of Clusters 3-to 5-specific marker genes were injected into SoFa transgenic embryos. We designed primers outside of the two distal sgRNA sites to test the gene disruption efficiency (primers are listed in Table S5). As this gene disruption system could cause a large fragment deletion, the predicted PCR results after gene disruption would show shorter or more diffused bands compared with the bands in the control groups. To further confirm the efficiency of gene disruption, PCR products from two separate experiments were cloned into blunt zero vectors (TransGen Biotech) and individual clones sequenced. Retinal sizes were observed at 1 dpf and 3 dpf, and the effects of different cell types were observed at 5 dpf.

We also used dgRNPs to carry out Cluster 1-specific gene disruption. dgRNPs were designed to mimic the native form of the CRISPR/Cas9 system in Streptococcus pyogenes. They are composed of the Cas9 protein and a duplex RNA formed by target-specific crispr (cr)RNA and trans-activating crispr (tracr)RNA. Both RNAs were chemically synthesized with Alt-R modifications (IDT) to avoid degradation in host cells. Target site specificity is provided by crRNA. The preparations of crRNA:tracrRNA duplex and dgRNPs were according to methods published recently (Hoshijima et al., 2019). Potential gRNA target sites were designed using the website provided by IDT (https://sg.idtdna.com/site/order/designtool/index/CRISPR_PREDESIGN). The target sites were: mab21l2 (CATCGCACCCAACGAGTTTG; CATCGCCAAGACCATACGAG); nr2f5 (GGGACACCGGGAACTCCCTC; GGTACGGCTGACCATTATAC); scrambled (CAGGCAAAGAATCCCTGCC).

Cell type counting after gene disruption

The SoFa transgenic embryos after disruption of Clusters 3- to 5-specific genes were fixed at 5 dpf at 4°C overnight. The embryos were then dehydrated in 30% sucrose and 20-μm-thick slices sectioned coronally. After staining with 4′,6-diamidino-2-phenylindole (DAPI), these slices were imaged under the confocal microscope (Olympus FV1200) with 60× (s, NA=1.30) objective. Visualization and measurement were performed using ImageJ (version 2.35). The sections with the optic nerve were selected. A square area, in which one side was aligned with the optic nerve and one corner was in the center of the lens, was cropped. We drew the regions of interest (ROIs) of each retinal type based on the signals of SoFa fishlines. After that, we counted each type using the ITCN plugin from ImageJ.

In situ hybridization

The digoxigenin (DIG)-labeled fabp11a, her4.2, atoh7, vsx1, mafba, her6, her9 and hes6 anti-sense and sense probes were prepared using the T7 RNA Polymerase kit and DIG RNA Labeling kit (Roche). The primers used for synthesizing probes are listed in Table S6. Zebrafish embryos (24 hpf, 36 hpf, 48 hpf and 14 dpf) were fixed in 4% paraformaldehyde (PFA; Electron Microscopy Services) at 4°C overnight. These fixed embryos were dehydrated in 30% sucrose and then 20-μm-thick slices were sectioned coronally. In situ hybridization was performed as described previously (Tang et al., 2017). Slices were incubated in TNB buffer [0.5% blocking reagent (Roche) in TN buffer] with sheep anti-DIG-POD antibody (1:500; Roche, #11207733910) at 4°C overnight. Finally, the signals were detected using the Tyramide Signal Amplification kit (Invitrogen). TN buffer was prepared with the final concentration of 0.1 M Tris-HCl (pH 7.5) and 0.15 M NaCl. Following staining, slices were imaged under an inverted confocal microscope (Olympus FV1200) with 60× (w, NA=1.20) objective. All visualizations were performed with FV10-ASW 4.0 Viewer (Olympus) and ImageJ (version 2.35).

Immunostaining

For Rho immunostaining, Tg(xops:nfsB-mcherry) fish at 14 dpf were fixed in 4% PFA (Electron Microscopy Services) and cryoprotected in 30% sucrose (4°C, overnight). For Caspase 3 immunostaining, wild-type embryos with or without disruptions of Cluster 1-specific genes at 30 hpf were fixed in 4% PFA (Electron Microscopy Services) and cryoprotected in 30% sucrose (4°C, overnight). Then the samples were flash-frozen and cryosectioned at a thickness of 14 μm coronally. Fluorescent immunochemistry was performed on sections containing zebrafish retinas as described previously (Tang et al., 2017).

Sections were washed three times with 1× PBS for 15 min and permeabilized in 1× PBS with 0.5% Triton X-100 for 30 min. After blocking with 5% BSA solution (Sigma-Aldrich) at room temperature for 1 h, sections were incubated with the primary antibody Rho (1:1000; Abcam, #ab98887) or Caspase 3 (1:1000; Cell Signaling Technology, #9661T) at 4°C overnight. Then sections were washed with 1× PBS and incubated with Alexa Fluor 488-conjugated secondary antibody (1:1000; Jackson ImmunoResearch, #33206ES60) or Alexa Fluor 568-conjugated secondary antibody (1:1000; Invitrogen, #A10042) at room temperature for 2 h. DAPI staining was performed according to the standard protocol. Slides were finally mounted using the fluorescent mounting medium (Sigma-Aldrich).

For Caspase 3 staining at 12 hpf, whole embryos were fixed at 4°C overnight in 4% PFA in PBS with 0.25% Triton X-100 (PBT). After that, embryos were washed with PBT and then immersed in 150 mM Tris-HCl (pH 9.0) for 5 min at room temperature, then heated at 70°C for 15 min. They were then washed with PBT and treated with 0.05% Trypsin-EDTA in ice for 5 min. After a quick wash in PBT, embryos were blocked in blocking buffer (2% normal fetal bovine serum, 1% BSA, 1% DMSO in PBT) for 1 h at room temperature. Then embryos were incubated with Caspase 3 (1:1000; Cell Signaling Technology, #9661T) at 4°C overnight. After washing with PBT, embryos were incubated with Alexa Fluor 568-conjugated secondary antibody (1:1000; Invitrogen, #A10042) at room temperature for 2 h and with DAPI at room temperature for 15 min.

Plasmid construction

The Xops:nfsB-mCherry plasmid was constructed as follows: the xops promoter was amplified from the plasmid pFIN-XOPS-tdTOM (xops-tdTOM) plasmid (Addgene #44359) with primers xopsF 5′-TATAGGGCGAATTGGtatagggcgaattggGGCCGCAGATCTTTATACATTGC-3′ and xopsR 5′-CCGGTGGATCCCAAACCCTCGAGATCCCTAGAAGCCTGTGAT-3′ and then subcloned into pT2uas:nfsB-mCherry plasmid (from Dr Toshio Ohshima's lab, Waseda University, Japan) to substitute the UAS promoter using a homologous recombination kit (ClonExpressMultiS One Step Cloning Kit, Vazyme).

EdU staining

Tg(xops:nfsB-mCherry) zebrafish at 5 dpf were injected with 1 nl 10 mM EdU into the yolk. These fish were fixed at 20 dpf at 4°C overnight. The fixed fish were dehydrated in 30% sucrose and then 14-μm-thick slices sectioned transversely. The EdU staining was performed according to the protocol of the Click-iT™ EdU Alexa Fluor™ 647 Imaging Kit (Invitrogen). In total six retinas and 86 sections were used for the statistical analysis. The nfsB-mCherry-positive cells in the central retina (between the EdU-positive clones) and the peripheral retina (between the CMZ and EdU-positive clones) were counted. The length of the corresponding outer plexiform layer was measured using ImageJ (version 2.35).

Statistical analysis

All statistical analyses were performed in R (version 3.4.3/3.5.1/3.6.1). Statistical analyses were performed using the Fisher's exact test and Wilcoxon test. Data are presented as mean±s.e.m. *P<0.05; **P<0.01; ***P<0.001; ns, no significance.

We thank Prof. William A. Harris for deep discussions and thoughtful suggestions during the entire manuscript preparation. We also thank Dr Toshio Ohshima at Waseda University for kindly providing the pT2uas:nfsB-mCherry plasmid. We thank Ms Xiaoying Qiu and Ms Xinling Jia from the Institute of Neuroscience, Chinese Academy of Sciences, for the help in gene disruption. We thank Dr Min Zhang from the Molecular Technology Facility at the Institute of Neuroscience, Chinese Academy of Sciences, for her help in 10x single-cell library construction.

Author contributions

Conceptualization: B.X., X.T., J.H.; Methodology: B.X., X.T., M.J., H.Z., L.D., S.Y.; Validation: B.X., X.T., M.J., J.H.; Formal analysis: B.X., X.T.; Investigation: B.X., X.T., M.J.; Resources: B.X., X.T., H.Z., L.D., S.Y.; Data curation: B.X., M.J.; Writing - original draft: J.H.; Writing - review & editing: B.X., X.T., J.H.; Visualization: B.X., X.T.; Supervision: J.H.; Project administration: J.H.; Funding acquisition: J.H.

Funding

This work was supported by the National Natural Science Foundation of China (31871035), the Strategic Priority Research Program of Chinese Academy of Sciences (XDB32000000), Shanghai Basic Research Field Project (18JC1410100), Science and Technology Commission of Shanghai Municipality Major Project (2018SHZDZX05).

Data availability

The datasets generated and analyzed for this study are available in GEO under accession number GSE122680.

Abrams
,
E. W.
and
Mullins
,
M. C.
(
2009
).
Early zebrafish development: it's in the maternal genes
.
Curr. Opin. Genet. Dev.
19
,
396
-
403
.
Agathocleous
,
M.
and
Harris
,
W. A.
(
2009
).
From progenitors to differentiated cells in the vertebrate retina
.
Annu. Rev. Cell Dev. Biol.
25
,
45
-
69
.
Almeida
,
A. D.
,
Boije
,
H.
,
Chow
,
R. W.
,
He
,
J.
,
Tham
,
J.
,
Suzuki
,
S. C.
and
Harris
,
W. A.
(
2014
).
Spectrum of Fates: a new approach to the study of the developing zebrafish retina
.
Development
141
,
1971
-
1980
.
Amini
,
R.
,
Rocha-Martins
,
M.
and
Norden
,
C.
(
2017
).
Neuronal migration and lamination in the vertebrate retina
.
Front. Neurosci.
11
,
742
.
Bassett
,
E. A.
and
Wallace
,
V. A.
(
2012
).
Cell fate determination in the vertebrate retina
.
Trends Neurosci.
35
,
565
-
573
.
Beby
,
F.
and
Lamonerie
,
T.
(
2013
).
The homeobox gene Otx2 in development and disease
.
Exp. Eye Res.
111
,
9
-
16
.
Boije
,
H.
,
Rulands
,
S.
,
Dudczig
,
S.
,
Simons
,
B. D.
and
Harris
,
W. A.
(
2015
).
The independent probabilistic firing of transcription factors: a paradigm for clonal variability in the zebrafish retina
.
Dev. Cell
34
,
532
-
543
.
Brown
,
N. L.
,
Patel
,
S.
,
Brzezinski
,
J.
and
Glaser
,
T.
(
2001
).
Math5 is required for retinal ganglion cell and optic nerve formation
.
Development
128
,
2497
-
2508
.
Butler
,
A.
,
Hoffman
,
P.
,
Smibert
,
P.
,
Papalexi
,
E.
and
Satija
,
R.
(
2018
).
Integrating single-cell transcriptomic data across different conditions, technologies, and species
.
Nat. Biotechnol.
36
,
411
-
420
.
Centanin
,
L.
,
Hoeckendorf
,
B.
and
Wittbrodt
,
J.
(
2011
).
Fate restriction and multipotency in retinal stem cells
.
Cell Stem Cell
9
,
553
-
562
.
Cepko
,
C.
(
2014
).
Intrinsically different retinal progenitor cells produce specific types of progeny
.
Nat. Rev. Neurosci.
15
,
615
-
627
.
Cherry
,
T. J.
,
Wang
,
S.
,
Bormuth
,
I.
,
Schwab
,
M.
,
Olson
,
J.
and
Cepko
,
C. L.
(
2011
).
NeuroD factors regulate cell fate and neurite stratification in the developing retina
.
J. Neurosci.
31
,
7365
-
7379
.
Chiodini
,
F.
,
Matter-Sadzinski
,
L.
,
Rodrigues
,
T.
,
Skowronska-Krawczyk
,
D.
,
Brodier
,
L.
,
Schaad
,
O.
,
Bauer
,
C.
,
Ballivet
,
M.
and
Matter
,
J.-M.
(
2013
).
A positive feedback loop between ATOH7 and a Notch effector regulates cell-cycle progression and neurogenesis in the retina
.
Cell Rep.
3
,
796
-
807
.
Clark
,
B. S.
,
Stein-O'Brien
,
G. L.
,
Shiau
,
F.
,
Cannon
,
G. H.
,
Davis-Marcisak
,
E.
,
Sherman
,
T.
,
Santiago
,
C. P.
,
Hoang
,
T. V.
,
Rajaii
,
F.
,
James-Esposito
,
R. E.
, et al. 
(
2019
).
Single-cell RNA-Seq analysis of retinal development identifies NFI factors as regulating mitotic exit and late-born cell specification
.
Neuron
102
,
1111
-
1126
.
Collery
,
R.
,
McLoughlin
,
S.
,
Vendrell
,
V.
,
Finnegan
,
J.
,
Crabb
,
J. W.
,
Saari
,
J. C.
and
Kennedy
,
B. N.
(
2008
).
Duplication and divergence of zebrafish CRALBP genes uncovers novel role for RPE- and Müller-CRALBP in cone vision
.
Invest. Ophthalmol. Vis. Sci.
49
,
3812
-
3820
.
Davis
,
D. M.
and
Dyer
,
M. A.
(
2010
).
Retinal progenitor cells, differentiation, and barriers to cell cycle reentry
.
Curr. Top. Dev. Biol.
93
,
175
-
188
.
Diotel
,
N.
,
Vaillant
,
C.
,
Kah
,
O.
and
Pellegrini
,
E.
(
2016
).
Mapping of brain lipid binding protein (Blbp) in the brain of adult zebrafish, co-expression with aromatase B and links with proliferation
.
Gene Expr. Patterns
20
,
42
-
54
.
Dixit
,
R.
,
Tachibana
,
N.
,
Touahri
,
Y.
,
Zinyk
,
D.
,
Logan
,
C.
and
Schuurmans
,
C.
(
2014
).
Gene expression is dynamically regulated in retinal progenitor cells prior to and during overt cellular differentiation
.
Gene Expr. Patterns
14
,
42
-
54
.
D'Orazi
,
F. D.
,
Suzuki
,
S. C.
and
Wong
,
R. O.
(
2014
).
Neuronal remodeling in retinal circuit assembly, disassembly, and reassembly
.
Trends Neurosci.
37
,
594
-
603
.
Duester
,
G.
(
2009
).
Keeping an eye on retinoic acid signaling during eye development
.
Chem-Biol. Interact.
178
,
178
-
181
.
Dyer
,
M. A.
,
Livesey
,
F. J.
,
Cepko
,
C. L.
and
Oliver
,
G.
(
2003
).
Prox1 function controls progenitor cell proliferation and horizontal cell genesis in the mammalian retina
.
Nat. Genet.
34
,
53
-
58
.
Fischer
,
A. J.
,
Bosse
,
J. L.
and
El-Hodiri
,
H. M.
(
2013
).
The ciliary marginal zone (CMZ) in development and regeneration of the vertebrate eye
.
Exp. Eye Res.
116
,
199
-
204
.
Gao
,
X.
,
Hu
,
D.
,
Gogol
,
M.
and
Li
,
H.
(
2019
).
ClusterMap: compare multiple single cell RNA-Seq datasets across different experimental conditions
.
Bioinformatics
35
,
3038
-
3045
.
Gath
,
N.
and
Gross
,
J. M.
(
2019
).
Zebrafish mab21l2 mutants possess severe defects in optic cup morphogenesis, lens and cornea development
.
Dev. Dyn.
248
,
514
-
529
.
Gonzalez-Nunez
,
V.
,
Nocco
,
V.
and
Budd
,
A.
(
2010
).
Characterization of drCol 15a1b: a novel component of the stem cell niche in the zebrafish retina
.
Stem Cells
28
,
1399
-
1411
.
Gross
,
J. M.
,
Perkins
,
B. D.
,
Amsterdam
,
A.
,
Egaña
,
A.
,
Darland
,
T.
,
Matsui
,
J. I.
,
Sciascia
,
S.
,
Hopkins
,
N.
and
Dowling
,
J. E.
(
2005
).
Identification of zebrafish insertional mutants with defects in visual system development and function
.
Genetics
170
,
245
-
261
.
Haider
,
N. B.
,
Mollema
,
N.
,
Gaule
,
M.
,
Yuan
,
Y.
,
Sachs
,
A. J.
,
Nystuen
,
A. M.
,
Naggert
,
J. K.
and
Nishina
,
P. M.
(
2009
).
Nr2e3-directed transcriptional regulation of genes involved in photoreceptor development and cell-type specific phototransduction
.
Exp. Eye Res.
89
,
365
-
372
.
Harris
,
W. A.
and
Perron
,
M.
(
1998
).
Molecular recapitulation: the growth of the vertebrate retina
.
Int. J. Dev. Biol.
42
,
299
-
304
.
He
,
J.
,
Zhang
,
G.
,
Almeida
,
A. D.
,
Cayouette
,
M.
,
Simons
,
B. D.
and
Harris
,
W. A.
(
2012
).
How variable clones build an invariant retina
.
Neuron
75
,
786
-
798
.
Hitchcock
,
P.
,
Ochocinska
,
M.
,
Sieh
,
A.
and
Otteson
,
D.
(
2004
).
Persistent and injury-induced neurogenesis in the vertebrate retina
.
Prog. Retin. Eye Res.
23
,
183
-
194
.
Holt
,
C. E.
,
Bertsch
,
T. W.
,
Ellis
,
H. M.
and
Harris
,
W. A.
(
1988
).
Cellular determination in the Xenopus retina is independent of lineage and birth date
.
Neuron
1
,
15
-
26
.
Hoon
,
M.
,
Okawa
,
H.
,
Della Santina
,
L.
and
Wong
,
R. O. L.
(
2014
).
Functional architecture of the retina: development and disease
.
Prog. Retin. Eye Res.
42
,
44
-
84
.
Hoshijima
,
K.
,
Jurynec
,
M. J.
,
Klatt Shaw
,
D.
,
Jacobi
,
A. M.
,
Behlke
,
M. A.
and
Grunwald
,
D. J.
(
2019
).
Highly efficient CRISPR-Cas9-based methods for generating deletion mutations and F0 embryos that lack gene function in zebrafish
.
Dev. Cell
51
,
645
-
657.e644
.
Hu
,
M.
and
Easter
,
S. S.
(
1999
).
Retinal neurogenesis: the formation of the initial central patch of postmitotic cells
.
Dev. Biol.
207
,
309
-
321
.
Johns
,
P. R.
(
1977
).
Growth of the adult goldfish eye. III. Source of the new retinal cells
.
J. Comp. Neurol.
176
,
343
-
357
.
Kageyama
,
R.
and
Ohtsuka
,
T.
(
1999
).
The Notch-Hes pathway in mammalian neural development
.
Cell Res.
9
,
179
-
188
.
Kageyama
,
R.
,
Ohtsuka
,
T.
and
Kobayashi
,
T.
(
2008
).
Roles of Hes genes in neural development
.
Dev. Growth Differ.
50
Suppl. 1
,
S97
-
S103
.
Kageyama
,
R.
,
Ohtsuka
,
T.
,
Shimojo
,
H.
and
Imayoshi
,
I.
(
2009
).
Dynamic regulation of Notch signaling in neural progenitor cells
.
Curr. Opin. Cell Biol.
21
,
733
-
740
.
Kageyama
,
R.
,
Shimojo
,
H.
and
Imayoshi
,
I.
(
2015
).
Dynamic expression and roles of Hes factors in neural development
.
Cell Tissue Res.
359
,
125
-
133
.
Kay
,
J. N.
,
Finger-Baier
,
K. C.
,
Roeser
,
T.
,
Staub
,
W.
and
Baier
,
H.
(
2001
).
Retinal ganglion cell genesis requires lakritz, a Zebrafish atonal Homolog
.
Neuron
30
,
725
-
736
.
Kay
,
J. N.
,
Link
,
B. A.
and
Baier
,
H.
(
2005
).
Staggered cell-intrinsic timing of ath5 expression underlies the wave of ganglion cell neurogenesis in the zebrafish retina
.
Development
132
,
2573
-
2585
.
Kim
,
J.-W.
,
Yang
,
H.-J.
,
Oel
,
A. P.
,
Brooks
,
M. J.
,
Jia
,
L.
,
Plachetzki
,
D. C.
,
Li
,
W.
,
Allison
,
W. T.
and
Swaroop
,
A.
(
2016
).
Recruitment of rod photoreceptors from short-wavelength-sensitive cones during the evolution of nocturnal vision in mammals
.
Dev. Cell
37
,
520
-
532
.
Kimmel
,
C. B.
,
Ballard
,
W. W.
,
Kimmel
,
S. R.
,
Ullmann
,
B.
and
Schilling
,
T. F.
(
1995
).
Stages of embryonic development of the zebrafish
.
Dev. Dyn.
203
,
253
-
310
.
Kobayashi
,
T.
and
Kageyama
,
R.
(
2014
).
Expression dynamics and functions of Hes factors in development and diseases
.
Curr. Top. Dev. Biol.
110
,
263
-
283
.
Koike
,
C.
,
Nishida
,
A.
,
Ueno
,
S.
,
Saito
,
H.
,
Sanuki
,
R.
,
Sato
,
S.
,
Furukawa
,
A.
,
Aizawa
,
S.
,
Matsuo
,
I.
,
Suzuki
,
N.
, et al. 
(
2007
).
Functional roles of Otx2 transcription factor in postnatal mouse retinal development
.
Mol. Cell. Biol.
27
,
8318
-
8329
.
Kriegstein
,
A.
and
Alvarez-Buylla
,
A.
(
2009
).
The glial nature of embryonic and adult neural stem cells
.
Annu. Rev. Neurosci.
32
,
149
-
184
.
Kwan
,
K. M.
,
Otsuna
,
H.
,
Kidokoro
,
H.
,
Carney
,
K. R.
,
Saijoh
,
Y.
and
Chien
,
C.-B.
(
2012
).
A complex choreography of cell movements shapes the vertebrate eye
.
Development
139
,
359
-
372
.
Lenkowski
,
J. R.
and
Raymond
,
P. A.
(
2014
).
Müller glia: stem cells for generation and regeneration of retinal neurons in teleost fish
.
Prog. Retin. Eye Res.
40
,
94
-
123
.
Li
,
S.
,
Mo
,
Z.
,
Yang
,
X.
,
Price
,
S. M.
,
Shen
,
M. M.
and
Xiang
,
M.
(
2004
).
Foxn4 controls the genesis of amacrine and horizontal cells by retinal progenitors
.
Neuron
43
,
795
-
807
.
Livesey
,
F. J.
and
Cepko
,
C. L.
(
2001
).
Vertebrate neural cell-fate determination: lessons from the retina
.
Nat. Rev. Neurosci.
2
,
109
-
118
.
Lopez-Ramirez
,
M. A.
,
Calvo
,
C.-F.
,
Ristori
,
E.
,
Thomas
,
J.-L.
and
Nicoli
,
S.
(
2016
).
Isolation and culture of adult zebrafish brain-derived neurospheres
.
J. Vis. Exp.
108
,
53617
.
Lu
,
Y.
,
Shiau
,
F.
,
Yi
,
W.
,
Lu
,
S.
,
Wu
,
Q.
,
Pearson
,
J. D.
,
Kallman
,
A.
,
Zhong
,
S.
,
Hoang
,
T.
,
Zuo
,
Z.
, et al. .
(
2020
).
Single-cell analysis of human retina identifies evolutionarily conserved and species-specific mechanisms controlling development
.
Dev. Cell
53
,
473
-
491
.
Malicki
,
J.
,
Pooranachandran
,
N.
,
Nikolaev
,
A.
,
Fang
,
X.
and
Avanesov
,
A.
(
2016
).
Analysis of the retina in the zebrafish model
.
Methods Cell Biol.
134
,
257
-
334
.
Marcucci
,
F.
,
Murcia-Belmonte
,
V.
,
Wang
,
Q.
,
Coca
,
Y.
,
Ferreiro-Galve
,
S.
,
Kuwajima
,
T.
,
Khalid
,
S.
,
Ross
,
M. E.
,
Mason
,
C.
and
Herrera
,
E.
(
2016
).
The ciliary margin zone of the mammalian retina generates retinal ganglion cells
.
Cell Rep.
17
,
3153
-
3164
.
Marcus
,
R. C.
,
Delaney
,
C. L.
and
Easter
,
S. S.
 Jr.
(
1999
).
Neurogenesis in the visual system of embryonic and adult zebrafish (Danio rerio)
.
Off. Vis. Neurosci.
16
,
417
-
424
.
Marquardt
,
T.
and
Gruss
,
P.
(
2002
).
Generating neuronal diversity in the retina: one for nearly all
.
Trends Neurosci.
25
,
32
-
38
.
Masai
,
I.
,
Stemple
,
D. L.
,
Okamoto
,
H.
and
Wilson
,
S. W.
(
2000
).
Midline signals regulate retinal neurogenesis in zebrafish
.
Neuron
27
,
251
-
263
.
Maurer
,
K. A.
,
Kowalchuk
,
A.
,
Shoja-Taheri
,
F.
and
Brown
,
N. L.
(
2018
).
Integral bHLH factor regulation of cell cycle exit and RGC differentiation
.
Dev. Dyn.
247
,
965
-
975
.
Merkle
,
F. T.
and
Alvarez-Buylla
,
A.
(
2006
).
Neural stem cells in mammalian development
.
Curr. Opin. Cell Biol.
18
,
704
-
709
.
Morris
,
A. C.
,
Schroeter
,
E. H.
,
Bilotta
,
J.
,
Wong
,
R. O. L.
and
Fadool
,
J. M.
(
2005
).
Cone survival despite rod degeneration in XOPS-mCFP transgenic zebrafish
.
Invest. Ophthalmol. Vis. Sci.
46
,
4762
-
4771
.
Moshiri
,
A.
,
Close
,
J.
and
Reh
,
T. A.
(
2004
).
Retinal stem cells and regeneration
.
Int. J. Dev. Biol.
48
,
1003
-
1014
.
Murai
,
K.
,
Philpott
,
A.
and
Jones
,
P. H.
(
2011
).
Hes6 is required for the neurogenic activity of neurogenin and NeuroD
.
PLoS ONE
6
,
e27880
.
Nawrocki
,
L. W.
(
1985
).
Development of the Neural retina in the Zebrafish, Brachydanio rerio
.
Eugine, Oregon
:
University of Oregon
.
Otteson
,
D. C.
and
Hitchcock
,
P. F.
(
2003
).
Stem cells in the teleost retina: persistent neurogenesis and injury-induced regeneration
.
Vision Res.
43
,
927
-
936
.
Paridaen
,
J. T. M. L.
and
Huttner
,
W. B.
(
2014
).
Neurogenesis during development of the vertebrate central nervous system
.
EMBO Rep.
15
,
351
-
364
.
Pelegri
,
F.
(
2003
).
Maternal factors in zebrafish development
.
Dev. Dyn.
228
,
535
-
554
.
Perron
,
M.
and
Harris
,
W. A.
(
2000
).
Retinal stem cells in vertebrates
.
BioEssays
22
,
685
-
688
.
Perron
,
M.
,
Kanekar
,
S.
,
Vetter
,
M. L.
and
Harris
,
W. A.
(
1998
).
The genetic sequence of retinal development in the ciliary margin of the Xenopus eye
.
Dev. Biol.
199
,
185
-
200
.
Peters
,
M. A.
(
2002
).
Patterning the neural retina
.
Curr. Opin. Neurobiol.
12
,
43
-
48
.
Picker
,
A.
,
Cavodeassi
,
F.
,
Machate
,
A.
,
Bernauer
,
S.
,
Hans
,
S.
,
Abe
,
G.
,
Kawakami
,
K.
,
Wilson
,
S. W.
and
Brand
,
M.
(
2009
).
Dynamic coupling of pattern formation and morphogenesis in the developing vertebrate retina
.
PLoS Biol.
7
,
e1000214
.
Pinto
,
L.
and
Gotz
,
M.
(
2007
).
Radial glial cell heterogeneity--the source of diverse progeny in the CNS
.
Prog. Neurobiol.
83
,
2
-
23
.
Podgornyi
,
O. V.
and
Aleksandrova
,
M. A.
(
2009
).
BLBP-immunoreactive cells in the primary culture of neural precursors from embryonic mouse brain
.
Bull. Exp. Biol. Med.
147
,
125
-
131
.
Qi
,
J.
,
Dong
,
Z.
,
Shi
,
Y.
,
Wang
,
X.
,
Qin
,
Y.
,
Wang
,
Y.
and
Liu
,
D.
(
2016
).
NgAgo-based fabp11a gene knockdown causes eye developmental defects in zebrafish
.
Cell Res.
26
,
1349
-
1352
.
Qiu
,
X.
,
Hill
,
A.
,
Packer
,
J.
,
Lin
,
D.
,
Ma
,
Y.-A.
and
Trapnell
,
C.
(
2017a
).
Single-cell mRNA quantification and differential analysis with Census
.
Nat. Methods
14
,
309
-
315
.
Qiu
,
X.
,
Mao
,
Q.
,
Tang
,
Y.
,
Wang
,
L.
,
Chawla
,
R.
,
Pliner
,
H. A.
and
Trapnell
,
C.
(
2017b
).
Reversed graph embedding resolves complex single-cell trajectories
.
Nat. Methods
14
,
979
-
982
.
Rapaport
,
D. H.
,
Wong
,
L. L.
,
Wood
,
E. D.
,
Yasumura
,
D.
and
LaVail
,
M. M.
(
2004
).
Timing and topography of cell genesis in the rat retina
.
J. Comp. Neurol.
474
,
304
-
324
.
Raymond
,
P. A.
,
Barthel
,
L. K.
,
Bernardos
,
R. L.
and
Perkowski
,
J. J.
(
2006
).
Molecular characterization of retinal stem cells and their niches in adult zebrafish
.
BMC Dev. Biol.
6
,
36
.
Reh
,
T. A.
and
Fischer
,
A. J.
(
2001
).
Stem cells in the vertebrate retina
.
Brain Behav. Evol.
58
,
296
-
305
.
Reinhardt
,
R.
,
Centanin
,
L.
,
Tavhelidse
,
T.
,
Inoue
,
D.
,
Wittbrodt
,
B.
,
Concordet
,
J. P.
,
Martinez-Morales
,
J. R.
and
Wittbrodt
,
J.
(
2015
).
Sox2, Tlx, Gli3, and Her9 converge on Rx2 to define retinal stem cells in vivo
.
EMBO J.
34
,
1572
-
1588
.
Rembold
,
M.
,
Loosli
,
F.
,
Adams
,
R. J.
and
Wittbrodt
,
J.
(
2006
).
Individual cell migration serves as the driving force for optic vesicle evagination
.
Science
313
,
1130
-
1134
.
Sauka-Spengler
,
T.
,
Baratte
,
B.
,
Shi
,
L.
and
Mazan
,
S.
(
2001
).
Structure and expression of an Otx5-related gene in the dogfish Scyliorhinus canicula: evidence for a conserved role of Otx5 and Crxgenes in the specification of photoreceptors
.
Dev. Genes Evol.
211
,
533
-
544
.
Seo
,
S.
,
Lim
,
J.-W.
,
Yellajoshyula
,
D.
,
Chang
,
L.-W.
and
Kroll
,
K. L.
(
2007
).
Neurogenin and NeuroD direct transcriptional targets and their regulatory enhancers
.
EMBO J.
26
,
5093
-
5108
.
Sghari
,
S.
and
Gunhaga
,
L.
(
2018
).
Temporal requirement of Mab21l2 during eye development in chick reveals stage-dependent functions for retinogenesis
.
Invest. Ophthalmol. Vis. Sci.
59
,
3869
-
3878
.
Shi
,
Z.
,
Trenholm
,
S.
,
Zhu
,
M.
,
Buddingh
,
S.
,
Star
,
E. N.
,
Awatramani
,
G. B.
and
Chow
,
R. L.
(
2011
).
Vsx1 regulates terminal differentiation of type 7 ON bipolar cells
.
J. Neurosci.
31
,
13118
-
13127
.
Straznicky
,
K.
and
Gaze
,
R. M.
(
1971
).
The growth of the retina in Xenopus laevis: an autoradiographic study
.
J. Embryol. Exp. Morphol.
26
,
67
-
79
.
Stuart
,
G. W.
,
Vielkind
,
J. R.
,
McMurray
,
J. V.
and
Westerfield
,
M.
(
1990
).
Stable lines of transgenic zebrafish exhibit reproducible patterns of transgene expression
.
Development
109
,
577
-
584
.
Stuart
,
T.
,
Butler
,
A.
,
Hoffman
,
P.
,
Hafemeister
,
C.
,
Papalexi
,
E.
,
Mauck
,
W. M.
, III
,
Hao
,
Y.
,
Stoeckius
,
M.
,
Smibert
,
P.
and
Satija
,
R.
(
2019
).
Comprehensive integration of single-cell data
.
Cell
177
,
1888
-
1902.e1821
.
Tallafuss
,
A.
,
Trepman
,
A.
and
Eisen
,
J. S.
(
2009
).
DeltaA mRNA and protein distribution in the zebrafish nervous system
.
Dev. Dyn.
238
,
3226
-
3236
.
Tang
,
X.
,
Gao
,
J.
,
Jia
,
X.
,
Zhao
,
W.
,
Zhang
,
Y.
,
Pan
,
W.
and
He
,
J.
(
2017
).
Bipotent progenitors as embryonic origin of retinal stem cells
.
J. Cell Biol.
216
,
1833
-
1847
.
Taverna
,
E.
,
Götz
,
M.
and
Huttner
,
W. B.
(
2014
).
The cell biology of neurogenesis: toward an understanding of the development and evolution of the neocortex
.
Annu. Rev. Cell Dev. Biol.
30
,
465
-
502
.
Trapnell
,
C.
,
Cacchiarelli
,
D.
,
Grimsby
,
J.
,
Pokharel
,
P.
,
Li
,
S.
,
Morse
,
M.
,
Lennon
,
N. J.
,
Livak
,
K. J.
,
Mikkelsen
,
T. S.
and
Rinn
,
J. L.
(
2014
).
The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells
.
Nat. Biotechnol.
32
,
381
-
386
.
Turner
,
D. L.
and
Cepko
,
C. L.
(
1987
).
A common progenitor for neurons and glia persists in rat retina late in development
.
Nature
328
,
131
-
136
.
Turner
,
D. L.
,
Snyder
,
E. Y.
and
Cepko
,
C. L.
(
1990
).
Lineage-independent determination of cell type in the embryonic mouse retina
.
Neuron
4
,
833
-
845
.
Venters
,
S. J.
,
Mikawa
,
T.
and
Hyer
,
J.
(
2015
).
Early divergence of central and peripheral neural retina precursors during vertebrate eye development
.
Dev. Dyn.
244
,
266
-
276
.
Wallace
,
V. A.
(
2011
).
Concise review: making a retina--from the building blocks to clinical applications
.
Stem Cells
29
,
412
-
417
.
Wan
,
Y.
,
Almeida
,
A. D.
,
Rulands
,
S.
,
Chalour
,
N.
,
Muresan
,
L.
,
Wu
,
Y.
,
Simons
,
B. D.
,
He
,
J.
and
Harris
,
W. A.
(
2016
).
The ciliary marginal zone of the zebrafish retina: clonal and time-lapse analysis of a continuously growing tissue
.
Development
143
,
1099
-
1107
.
Wang
,
S.
,
Sengel
,
C.
,
Emerson
,
M. M.
and
Cepko
,
C. L.
(
2014
).
A gene regulatory network controls the binary fate decision of rod and bipolar cells in the vertebrate retina
.
Dev. Cell
30
,
513
-
527
.
Webber
,
A. L.
,
Hodor
,
P.
,
Thut
,
C. J.
,
Vogt
,
T. F.
,
Zhang
,
T.
,
Holder
,
D. J.
and
Petrukhin
,
K.
(
2008
).
Dual role of Nr2e3 in photoreceptor development and maintenance
.
Exp. Eye Res.
87
,
35
-
48
.
Wehman
,
A. M.
,
Staub
,
W.
,
Meyers
,
J. R.
,
Raymond
,
P. A.
and
Baier
,
H.
(
2005
).
Genetic dissection of the zebrafish retinal stem-cell compartment
.
Dev. Biol.
281
,
53
-
65
.
Wetts
,
R.
and
Fraser
,
S. E.
(
1988
).
Multipotent precursors can give rise to all major cell types of the frog retina
.
Science
239
,
1142
-
1145
.
Wilson
,
S. W.
and
Houart
,
C.
(
2004
).
Early steps in the development of the forebrain
.
Dev. Cell
6
,
167
-
181
.
Wu
,
F.
,
Li
,
R.
,
Umino
,
Y.
,
Kaczynski
,
T. J.
,
Sapkota
,
D.
,
Li
,
S.
,
Xiang
,
M.
,
Fliesler
,
S. J.
,
Sherry
,
D. M.
,
Gannon
,
M.
, et al. 
(
2013
).
Onecut1 is essential for horizontal cell genesis and retinal integrity
.
J. Neurosci.
33
,
13053
-
13065
,
13065a
.
Wu
,
R. S.
,
Lam
,
I. I.
,
Clay
,
H.
,
Duong
,
D. N.
,
Deo
,
R. C.
and
Coughlin
,
S. R.
(
2018
).
A rapid method for directed gene knockout for screening in G0 zebrafish
.
Dev. Cell
46
,
112
-
125.e114
.
Yamada
,
R.
,
Mizutani-Koseki
,
Y.
,
Koseki
,
H.
and
Takahashi
,
N.
(
2004
).
Requirement for Mab21l2 during development of murine retina and ventral body wall
.
Dev. Biol.
274
,
295
-
307
.
Yu
,
G.
,
Wang
,
L.-G.
,
Han
,
Y.
and
He
,
Q.-Y.
(
2012
).
clusterProfiler: an R package for comparing biological themes among gene clusters
.
Omics
16
,
284
-
287
.

Competing interests

The authors declare no competing or financial interests.

Supplementary information