A comprehensive series of temporal transcription factors in the fly visual system

The brain consists of thousands of different neuronal types that are generated through multiple divisions of neuronal stem cells. These stem cells have the capacity to generate different neuronal types at different stages of their development. In Drosophila, this temporal patterning is driven by the successive expression of temporal transcription factors (tTFs). While a number of tTFs are known in different animals and across various parts of the nervous system, these have been mostly identified by informed guesses and antibody availability. We used single-cell mRNA sequencing to identify the complete series of tTFs that specify most Drosophila medulla neurons in the optic lobe. We tested the genetic interactions among these tTFs. While we verify the general principle that tTFs regulate the progression of the series by activating the next tTFs in the series and repressing the previous ones, we also identify more complex regulations. Two of the tTFs, Eyeless and Dichaete, act as hubs integrating the input of several upstream tTFs before allowing the series to progress and in turn regulating the expression of several downstream tTFs. Moreover, we show that tTFs not only specify neuronal identity by controlling the expression of cell type-specific genes. Finally, we describe the very first steps of neuronal differentiation and find that terminal differentiation genes, such as neurotransmitter-related genes, are present as transcripts, but not as proteins, in immature larval neurons days before they are being used in functioning neurons; we show that these mechanisms are conserved in humans. Our results offer a comprehensive description of a temporal series of tTFs in a neuronal system, offering mechanistic insights into the regulation of the progression of the series and the regulation of neuronal diversity. This represents a proof-of-principle for the use of single-cell mRNA sequencing for the comparison of temporal patterning across phyla that can lead to an understanding of how the human brain develops and how it has evolved.

sequencing for the comparison of temporal patterning across phyla that can lead to an understanding of how the human brain develops and how it has evolved. The brain is the most complex organ of the animal body: the human brain consists of over 80 billion neurons 1 that belong to thousands of different neuronal types and form ~10 14 synapses. Understanding the generation of this complexity in humans is an almost insurmountable problem. Researching the brains of simpler genetic model organisms as diverse as mice and flies provides windows into the underlying molecular mechanisms of how the production of neuronal diversity is achieved. This has highlighted the importance of two essential factors: the spatial location and age of neuronal progenitors. Spatial information has long been acknowledged for its importance in patterning the dorsoventral and anteroposterior axes of animal bodies. Classic examples include patterning of the Drosophila embryo as well as the vertebrate spinal cord 2,3 .
Morphogenetic gradients are converted into discrete spatial domains (e.g., the French Flag model 2 ); in the nervous system, these domains express different transcription factors (spatial transcription factors -sTFs) and give rise to fate-restricted groups of neural stem cells, each of which can generate a unique subset of neuronal types 4,5 .
There have been reports of orthologous transcription factors (TFs) acting as sTFs in mice and flies: for example, Drop (also called Msh -muscle segment homeobox)/Intermediate neuroblasts defective (Ind)/Ventral nervous system defective (Vnd) are expressed along the DV axis of the Drosophila neuroectoderm, and their mouse orthologs Msx1/Gsx1/Nkx2-2 are expressed in progenitor cells along the DV axis of the mouse neural tube at embryonic day 11.5 6 .
Temporal patterning describes the developmental trajectory neural stem cells follow that allows them to generate different neuronal types as they age 7,8 . It is a powerful mechanism that allows neural stem cells to produce different neuronal types in the correct order and stoichiometry. The first mechanism of temporal patterning in neuronal systems was described in the Drosophila ventral nerve cord (VNC) where a cascade of temporal transcription factors (tTFs) is expressed in embryonic neural stem cells (neuroblasts) as they divide and age [9][10][11] . It was later suggested that tTFs also contribute to the generation of neuronal diversity in different mammalian neuronal tissues, such as the retina [12][13][14][15] and the cortex 16 . However, only a few tTFs have been discovered to play a role in both insects and vertebrates, such as Ikzf1 (a mouse ortholog of Drosophila Hunchback that is the first tTF in the fly VNC) 11,17 that specifies young neural stem cells in the mouse cortex and retina 12,16 , as well as Pou2f1/Pou2f2 15  cortical radial glia, as well as their immediate descendants, at different time points during development, offering an excellent resource to identify mouse tTFs. Although many genes were found to be dynamically expressed as the apical cortical progenitors age, a series of tTFs was not reported. Nonetheless, they discovered that neurogenesis of cortical excitatory neurons is governed by two orthogonal (i.e. independent) processes, specification and differentiation, where different neuronal identities are specified through time in radial glia while neuronal differentiation follows a precast and highly similar program.
The Drosophila optic lobe is an ideal system to address how neuronal diversity is generated and how neurons proceed to differentiate. It is an experimentally manageable, albeit complex structure, for which we have a very comprehensive catalogue of neuronal cell types. It consists of four neuropils, the lamina, medulla, lobula, and lobula plate.
Meticulous work from the last decades have identified ~100 cell types in the optic lobes based solely on morphological characters 19 . Recent work took advantage of elaborate molecular genetic tools, as well as scRNASeq, to expand the number of neuronal cell types to ~200, based on both morphology and molecular identity [20][21][22] . Importantly, because the optic lobe processes visual information generated in each of the 800 unit eyes (ommatidia), it is formed of 800 similar circuits running in parallel, with many of the optic lobe cell types present in multiple copies in the brain (ranging from ~5 to 800). This makes clustering of scRNASeq data and cluster annotation easier. Moreover, the neuroblasts that generate the medulla, which is the largest neuropil of the optic lobe, are formed by a wave of neurogenesis over a period of days, and they all progress through the same tTF temporal series 23 . This means that at any given developmental stage from mid third larval stage (L3) to the beginning of pupation, the neurogenic region contains neuroblasts at all stages of their development when each must express one of the tTFs ( Figure 1A).
We performed scRNASeq of the developing Drosophila optic lobe at the time when neural stem cells (neuroblasts) divide to generate the ~200 neuron types that compose this brain structure. We focused on the neural stem cell that form the main structure of the optic lobe, the medulla. This allowed us to identify most, if not all, tTFs that are expressed in medulla neural stem cells, as well as to describe the genetic interactions among them that allow the tTF series to progress. Although we can define a general rule that a given tTF activates the next tTF and represses the previous one, we uncovered an unexpected degree of complexity, in which different tTFs assume different roles in the timing of the progression of the series and the generation of neuronal diversity. We also showed that tTFs control neuronal identity by regulating the expression of downstream cell type-specific TFs. Finally, we describe the very first steps that a neuron takes to differentiate; we find that all neurons express terminal differentiation genes as early as L3 and they follow a similar route for differentiation, independent of their identity, and that this progression is conserved in the human brain.

Single-cell transcriptomes recapitulate the structure of the developing optic lobe
The adult Drosophila optic lobes start developing at L3 from the lateral parts of the larval central nervous system. The medulla part of the optic lobe, which is by far the largest and most complex optic neuropil, is formed by medulla neuroblasts that are generated by a neurogenic epithelium called the outer proliferation center (OPC) 23 . Over a period of two days, the OPC is progressively converted into neuroblasts by a neurogenic wave that initiates medially and continues laterally until the entire epithelium is consumed. This process results in the generation of seemingly identical neuroblasts that produce neuronal types throughout optic lobe development, meaning that at any single point in time, there are medulla neuroblasts of different ages (young to old) ( Figure 1A). This characteristic of medulla development provides a distinct opportunity to study neuroblast and neuronal trajectories in unparalleled detail since all developmental stages coexist in the same brain. To achieve that, we performed scRNASeq on optic lobes microdissected from the central brain using the Chromium system (10x Genomics). We obtained 49,893 single-cell transcriptomes from 40 L3 optic lobes. We used the Seurat v3 integration pipeline 24 to remove batch effects between the ten different libraries that were generated (Extended Data Figure 1).  Figure 1B -Extended Data Figure 2A).
The OPC neuroepithelium generates two of the optic lobe structures: it is progressively converted from the medial side to give rise to medulla neuroblasts while its lateral side gives rise to the lamina precursor cells that will form the lamina 23 ( Figure 1A). Notably, medulla neuroepithelium, neuroblasts, GMCs and neurons were arranged in the UMAP following a progression that resembles the in vivo differentiation process ("medulla differentiation trajectory" ( Figure 1B). Similarly, lamina neuroepithelium, lamina precursor cells, and lamina monopolar neurons were also arranged following a similar differentiation trajectory ("lamina differentiation trajectory") but in the opposite orientation of that of the medulla, highlighting the similarities of the UMAP trajectories with the actual differentiation process in the brain ( Figure 1B). Emanating from the GMCs, different neuronal branches emerge that appear to represent developmental trajectories of different neurons ( Figure 1B). Lobula plate neurons are generated from the neuroblasts of the inner proliferation center (IPC) that are distinct from those of the OPC that generate the lamina and medulla. These neuroblasts and the neurons that are generated from the IPC follow a different trajectory in the UMAP plot ("IPC neuron differentiation trajectory", Figure 1B) and will not be discussed further.
To verify whether these trajectories retain the information of time, as suggested by the progression of neuroblasts to GMCs to neurons, we merged the larval single-cell dataset with the early pupal stage 15 (P15) dataset that we had previously generated 21 . These P15 neurons fell at the tip of each of the neuronal branches, confirming that these branches indeed represent neuronal trajectories ( Figure 1C). Importantly, since the P15 dataset is annotated, this allowed us to identify the neuronal types at the tip of the trajectories. We could in fact identify neurons from all the neuropils of the optic lobe (lamina, medulla, lobula, and lobula plate) in the larval UMAP, which together accounted for 85% of the dataset. The annotation of the different neuropils was confirmed by looking at known markers of the lamina (dac, gcm, so, eya, sim 33 ) and lobula plate (D, tll, acj6, dac 34,35 ) (Extended Data Figure 2B). The remaining cells included a small number of central brain neurons and neuroblasts that were retained when cutting off the optic lobe ( Figure 1D).
We then looked at the expression of the known spatial and temporal TFs in the neuroepithelium and neuroblasts, respectively. The neuroepithelium is divided into three broad domains by the expression of three spatial factors (Vsx, Optix, and Rx) 5 . These spatial factors were expressed in largely non-overlapping subsets of the neuroepithelial cells ( Figure 1E, Extended Data Figure 2C). We clustered the neuroepithelial cells and used Vsx1, Optix, and Rx expression in each cluster to assign them to a spatial domain ( Figure 1E'). The number of neuroepithelial cells corresponding to the different domain correlated with their size: Optix represented the largest spatial domain (spanning ~65% of the epithelium), followed by Rx (23%) and Vsx (12%). However, as previously shown by immunostainings, their expression was not maintained in neuroblasts and neurons (Extended Data Figure 2C).
The medulla neuroblasts express a series of five tTFs (Hth, Ey, Slp, D, and Tll) in a temporal manner 36 . tTFs showed a very distinct pattern in the UMAP plot: not only were they expressed in subsets of neuroblasts, but neuroblasts and GMCs were organized based on their age, with progenitors expressing a tTF positioned between those expressing the previous and the next one in the temporal cascade ( Figure 1F): Hth was present in the bottom row of the cluster, followed by Ey, Slp, D, and Tll with partial overlap among them, similar to what is observed with immunostaining in vivo 36 .
Interestingly, neuroblasts positioned between Hth and Ey were not expressing any of the known tTFs, as was expected from in vivo stainings 36 ( Figure 1F -arrow, discussed later) In general, we observed that the UMAP plot recapitulated remarkably well what is happening in the developing tissue: there were two different axes of time in the UMAP, a vertical axis that represents neuroblasts progressing through their temporal series of tTFs and a horizontal axis that represents cell state and differentiation status (i.e., neuroblast to GMC to immature neuron to mature neurons). We observed two bottlenecks along the developmental axis ( Figure 1B, red arrows), one when the neuroepithelium is converted into neuroblasts and one when the GMCs with different temporal identities converge transcriptionally before they diverge again towards separate neuronal trajectories. This might be due to the fact that gene transcription during epithelial to mesenchymal transition in the first case and before the terminal division of the GMCs in the second case are obscuring more specific identity features.

A comprehensive temporal series of transcription factors in the developing medulla
Although they cover almost the entire life of neuroblasts, the existing tTFs were discovered from educated guesses and screening of available antibodies. There is also clear evidence that there are additional tTFs, as the existing TFs are not able to explain the entire neuronal diversity in the optic lobe and there are neuroblasts between the Hth and Ey temporal windows that do not express any of the five tTFs in vivo 36 . We could confirm that there are cells in the UMAP plot that express none of the known tTFs ( Figure 1F -arrow). The larval scRNAseq dataset gave us the opportunity to look for all potential tTFs in an unbiased and comprehensive way. We isolated the cluster of medulla neuroblasts from the scRNASeq data and used Monocle 37 to reconstruct their developmental trajectory. To confirm the accuracy of the trajectory, we looked at the expression of the known tTFs: Hth, Ey, Slp, D, and Tll. Indeed, as was already clear from the UMAP plot, these tTFs were expressed in the correct temporal order along the trajectory ( Figure 2A -TFs in purple). We then examined the expression dynamics of all Drosophila TFs and identified 39 candidate tTFs that exhibited expression restricted to a temporal window. These fell into two distinct categories: 14 of them were expressed at relatively high levels and included the 6 known tTFs (Slp includes two genes, Slp1 and Slp2) (Extended Data Figure 3A), while 25 of them were expressed at lower and fluctuating levels along the trajectory (Extended Data Figure 3B). We tested the expression pattern of four of the 25 lowly expressed candidates (ap, cut, gcm, and gem) in the developing optic lobes using antibodies against Ap and the cut-Gal4, gcm-Gal4, and gem-Gal4 lines, but none were expressed in a temporal manner; therefore, we decided not to pursue these candidates further as their fluctuations likely represent noise. Moreover, although Klumpfuss (Klu) was previously suggested to be a tTF 38 , Klu mRNA was found to be continuously expressed throughout neuroblast life in our dataset.
We then tested the protein expression of the eight newly discovered candidate tTFs in medulla neuroblasts (surface view - Figure 2B). Using antibodies against these potential tTFs and the already known ones, we verified that their protein expression was limited to a specific temporal window ( Figure 2B' and 2C). They are described below in their order of expression in neuroblasts.
Opa (Odd-paired) is expressed in two waves: it is first expressed in young neuroblasts immediately after and partially overlapping with the Hth temporal window ( Figure 2D arrow). Then its expression ceases before reappearing just before Slp ( Figure 2D' and Extended Data Figure 4A). Erm (Earmuff) immediately follows Hth ( Figure 2E and Extended Data Figure 4B), partially overlaps with Opa and precedes (partially overlapping) Ey expression (Extended Data Figure 4C). Esg (Escargot) is expressed within the Ey temporal window, albeit in a salt-and-pepper manner, indicating that it is likely not a bona fide tTF as it is not expressed in all neuroblasts ( Figure 2F and Extended Data Figure 4D). Hbn (Homeobrain) expression almost completely overlaps Ey in neuroblasts (Extended Data Figure 4E), right before Slp1 ( Figure 2G and Extended Data Figure 4F). Scro (Scarecrow) expression starts immediately after Ey ( Figure 2H and Extended Data Figure 4G), but it remains expressed until the end of the neuroblast divisions ( Figure 2B'). BarH1 is expressed after D ( Figure 2I and Extended Data Figure   4H) and before Tll (Extended Data Figure 4I), partially overlapping both. Finally, in the absence of a functioning antibody, we tested the expression of Oaz using an Oaz-Gal4 line driving UAS-GFP. While based on the bioinformatic analysis it was expected to be expressed in young neuroblasts up to the Slp temporal window, it showed expression in all medulla NBs, potentially due to the perdurance of Gal4 in older neuroblasts (Extended Data Figure 4J).
Therefore, we could confirm that the predicted medulla tTF proteins (except potentially for Oaz) are expressed temporally in the developing optic lobe, defining new temporal windows as the neuroblasts progress through divisions ( Figure 2C).

Different tTFs assume different roles in the progression of the series
The known tTFs (except Hth) contribute to the progression of the series by activating the next tTF in the series and repressing the previous one 36 . While the TFs discovered above were expressed temporally, this does not imply that they actively participate in the progression of the temporal series. To test which of the newly identified tTFs were involved in the progression of the temporal series, we used mutant clones. We also  Figure 5A).

Early unit
Hth begins the temporal series: importantly, hth transcripts are present in the very first neuroblasts as well as in the neuroepithelium that has not yet been transformed into neuroblasts, indicating that its activation is likely regulated by upstream patterning events in the neuroepithelium. We had previously shown that Hth is not required for the progression of the temporal series as the next known tTF, Ey, is expressed normally in hth mutants 36 . Therefore, another overlapping factor must be responsible for activating Ey. In fact, we identified two factors that regulate the expression of Ey in different manners, Erm and Opa. Erm acts like the known tTFs as it is required to activate its next tTF, Ey and to inhibit the previous, Hth: In erm mutant clones, Ey is not expressed and Hth expression is expanded ( Figure 3B). At the same time, Opa, which is co-expressed in the last Hth NBs, is required for the activation of Ey at the correct time: opa mutant neuroblasts have strongly delayed expression of Ey ( Figure 3C), which leads to a delayed expression of the temporal series after Ey (Extended Data Figure 5E); Hth and Erm are unaffected in opa mutant neuroblasts (Extended Data Figure 5C). Once the Ey temporal window is initiated at the correct time by the combined action of Erm and Opa, Ey represses the expression of its activators, Opa and Erm: in ey mutant clones, both Erm ( Figure 3D) and Opa ( Figure 3E) are expanded to later temporal windows.
Therefore, Erm is a tTF essential for the progression of the cascade while Opa contributes to the correct timing of expression of the next tTFs.

Middle unit
We had shown that Ey activates and is then inhibited by Slp 36 . However, the developmental trajectory of neuroblasts uncovered a much more complex situation.
First, we found that Ey also activates Scro and is inhibited by it: in ey mutant clones Scro expression was completely lost ( Figure 3F), while when Scro was knocked down by RNAi, Ey remained expressed until the last division of the neuroblasts ( Figure 3G).
Moreover, Ey also activated Hbn and was inhibited by it: in ey mutant clones, Hbn expression was lost ( Figure 3H), while in hbn mutant clones, Ey was extended to later temporal windows ( Figure 3I). Then, Hbn allows the temporal series to progress by activating Slp as hbn mutant clones lacked Slp expression ( Figure 3I). This suggests that the activation of Slp by Ey is mediated by Hbn. Slp then inhibits Ey 36 and Hbn: slp mutant clones showed extension of Hbn into later temporal windows ( Figure 3J). Finally, Hbn activates the second temporal window of Opa, as in hbn mutant clones, the second wave of Opa expression was absent ( Figure 3K).
The complex genetic interactions that involve the activation of ey (temporal regulation by opa and regulation of expression by erm), as well as the fact that two genes (Scro and slp) are required for it to be repressed, indicate that Ey plays the role of a hub factor in the initiation and progression of the temporal series, as it integrates several signals before activating the expression of several downstream tTFs.
Late unit: Finally, D requires both Slp and Scro to be expressed. We had previously shown that in slp mutant clones, D is not expressed 36 . Similarly, when Scro was knocked down by RNAi, D was not expressed in the neuroblasts ( Figure 3L and Extended Data Figure   5G). Scro is therefore important for the progression of the series, as it is activated by Ey, which it then inhibits together with Slp; it then activates the expression of D and remains expressed until the end the neuroblast's life. Once D is activated, it inhibits Slp and activates BarH1: in D mutant clones, BarH1 expression was lost ( Figure 3M). At the same time, D activates Tll and Tll is sufficient to inhibit D, as was shown before 36 .
We have thus been able to identify most, if not all, temporally expressed TFs in a developing neuronal system and to show that these tTFs participate in the progression of the temporal series. By exhaustively examining the genetic interactions between the new and old tTFs ( Figure 3 and Extended Data Figure 5), we show that the temporal series is more complex than previously described and that not all tTFs play a similar role. While some tTFs directly inhibit the previous factor and/or activate the next (Erm, Hbn, Slp, and Tll), one does not participate in the progression of the temporal series

Temporal transcription factors often remain expressed in neurons and regulate neuronal identity
Besides their participation in the progression of the temporal series, tTFs have been shown to also regulate neuronal identity either by being expressed in the neuronal subsets that are generated during their temporal window and acting as effector TFs (i.e. activating effector genes) 36 , or by activating the expression of downstream transcription factors 36,38 , which then regulate the expression of effector genes in the absence of the tTF.
We first looked at the expression of tTFs in the neuronal progeny in the scRNASeq data.
Hth, Erm, Hbn, Ey, and D were expressed in largely non-overlapping neuronal clusters We then asked which neuronal types are generated from each temporal window; for this purpose, we used the expression of tTFs in neurons and GMCs described earlier to assign each neuronal type to a specific temporal window (Table S1). Except for the neurons described above that co-express Erm/Ey and D/Hbn, most tTFs are only expressed in neurons that come from their respective temporal window. We evaluated the assignments of the neurons to specific temporal windows by looking at neurons whose temporal window is already known, such as Mi1 (Hth temporal window), Tm3 (Ey temporal window), and Tm20 (Slp temporal window). We found that proximal medulla Otd ( Figure 4D-D') were no longer found. We also tested whether Opa was required for the generation of TfAP-2 positive neurons that are generated during both Opa temporal windows. In opa mutant clones, we found a significant reduction of TfAP-2 positive neurons compared to the adjacent wild-type tissue ( Figure 4E). This shows that Hbn and Opa are regulating neuronal diversity not only by allowing the temporal series to progress (by activating the expression of Slp and timing the expression of Ey, respectively), but also by regulating neuronal transcription factor expression.

The implementation of neuronal identity occurs very early in neuronal life and follows a fixed trajectory
Drosophila neurons are already specified at birth, based on the spatial and temporal identity of the neuroblasts from which they were born, as well as the Notch status of the neuron. However, it remains unclear how neurons proceed to differentiate. To study the very first steps of neuronal differentiation after specification, we initially focused on a well-studied and easily identifiable neuronal type, Mi1 (Mi1 is the only adult neuronal type in the medulla that expresses Bsh). We selected from the scRNASeq dataset of the L3 developing optic lobes the Bsh-positive clusters that correspond to Mi1 neurons at different levels of maturity, as well as the GMCs most closely linked to them in the UMAP plot (Extended Data Figure 6A-A'). We merged this dataset with the Mi1 clusters from pupal stages P15, P30, P40, P50, and P70, and used Monocle3 42 to reconstruct their differentiation trajectory ( Figure 5A). We used the expression of Ase in GMCs and Bsh in neurons to mark the beginning of the elongated trajectory at L3 and P15 (Extended Data Figure 6B quickly reaching a plateau that was maintained until P15. Their expression then started to increase again until adulthood when the products of these genes support neuronal function ( Figure 5B). This indicates that not only is neuronal identity specified during the first hours of neuronal development, but it is already implemented through genes involved in neuronal function very early (at least at the transcript level).
To evaluate whether this applies to all neurons or if each neuron follows a distinct differentiation trajectory, we aggregated and plotted the genes that belonged to the GO terms that were differentially activated during the Mi1 trajectory on the UMAP plot and observed that all optic lobe neurons followed the same differentiation trajectories, as indicated by the expression of the GO terms in all neuronal branches of the UMAP plot ( Figure 5C). This differentiation trajectory was not only observed in the medulla, but also in the lamina and lobula plate neurons ( Figure 5C).
As the expression of genes that should only be required very late was unexpected, we asked whether the transcript expression observed as early as late L3 was also translated into protein expression. We focused on neurotransmitter-related genes. We first verified that the correct neurotransmitter identity was indeed established as early as L3. We looked for the expression of ChAT, VGlut, and Gad1 (markers for cholinergic, glutamatergic, and GABAergic cell types, respectively) in the scRNASeq data and observed that they were expressed in the medulla in non-overlapping neuronal sets ( Figure 5D). We confirmed using the adult scRNASeq data that we published recently 21 that the neurotransmitter identity of each cell type at L3 was the same as in the adult (Table S1). We then stained late L3 brains with antibodies for choline acetyltransferase (ChAT) and the vesicular glutamate transporter (VGlut) in late L3. While we observed their expression in the mature neurons of the larval ventral nerve cord, we saw no expression in the developing optic lobe ( Figure 5E). This suggests that transcribing the loci might be a means towards commitment to a specific neurotransmitter identity but that other factors prevent translation of these mRNAs.
In conclusion, we show that the differentiation process of Drosophila optic lobe neurons is fixed and independent of neuronal identity: Acquisition and implementation of identity are two consecutive processes, where the temporal and spatial information inherited from the neuroblasts specify the genes that are expressed, while the differentiation trajectory decides the timing of their expression. This agrees with recent data from the mouse cortex, where specification and differentiation were proposed to be two independent processes that occur mainly in different cell types (stem cells vs neurons) and where differentiation follows a precast path in all neurons, independent of their identities 18 .

Common differentiation trajectory in Drosophila optic lobe and human cortical neurons
Understanding how neuronal differentiation occurs in human cortical neurons in vivo is necessary for the development of accurate in vitro differentiation protocols that can be used for neuronal replacement therapies 43 . We therefore wondered whether the differentiation trajectory we described in Drosophila optic lobe neurons was also implemented during human neuronal differentiation. We generated single-nuclear RNAsequencing data from the developing human fetal cortical plate at gestational week 19.
We used Monocle 3 to reconstruct their developmental trajectory (Extended Data Figure   6C). We could see a trajectory from radial glia progenitors to intermediate progenitors Despite this difference, these results show that the neurons follow a similar differentiation trajectory in Drosophila and humans that can be either attributed to convergence or to their common origin.

Drosophila tTF expression in mouse cortical radial glia
While similarities in spatial patterning between vertebrates and flies has long been noted 3,6 , it is not clear how temporal patterning has evolved 48,49 . Sporadic evidence suggests that temporal factors identified in the Drosophila ventral nerve cord neuroblasts pattern stem cells that generate the mouse retina and even the cortex [12][13][14][15][16] . To address the similarities and differences between Drosophila neuroblasts and mouse cortical radial glia, we probed for the expression of known Drosophila tTFs a recently published scRNASeq dataset from the mouse cortex, where radial glia at different stages of development (E12-E15) were sequenced 18 . We first looked for the orthologs of the known tTFs of the Drosophila optic lobes described earlier, as well as in this study. None of the medulla neuroblast tTFs were expressed in strict temporal windows in ageing radial glia between day 12, where they produce deep cortical layers, and day 15, when more superficial layers are generated (Extended Data Figure 7A), with the exception of the Ey ortholog, Pax6, which was enriched in older progenitors. Foxg1, an ortholog of Slp was previously described to be enriched in young radial glia 49,50 ; while it is slightly reduced in older radial glia in this dataset, this reduction was not statistically significant.
We also looked at the Drosophila orthologs of the mouse TFs that were described to be expressed in a temporal manner in the mouse CNS 51 ; none of them were found to be tTFs in our trajectory analysis. The lack of conservation of a common temporal series of transcription factors sensu stricto between flies and mice shows that the acquisition of the specific temporal series occurred later and independently in each phylum, which is consistent with the several cascades of tTFs observed in various brain structures, such as the ventral nerve cord 11 , Type II central brain neuroblasts 52 or optic lobe neuroblasts.

Discussion
We present here a comprehensive series of transcription factors that temporally pattern a developing neural structure. We took advantage of the unique structure of the developing Drosophila optic lobe and generated detailed trajectories of neural stem cells starting from the time they are born to the time they terminally differentiate. All known tTFs were confirmed using in this approach; moreover, all the candidate tTFs that were identified computationally were verified experimentally, thus providing a proof-ofprinciple for the combination of scRNASeq and trajectory inference that can be applicable to most other neuronal tissues in different animals that lack the genetic toolkit of Drosophila. We show that most tTFs are expressed in overlapping windows creating a combinatorial code that differentiates neural stem cells of different ages and therefore provide them with the ability to generate diverse neurons after every division. We conservatively assigned them into 12 distinct temporal windows ( Figure 2C), which, when integrated with spatial patterning (5 spatial domains) and the Notch binary cell fate decision (2 outcomes -Notch ON and Notch OFF ), can explain the generation of ~120 cell types (12 times 5 times 2), which is close to the entire neuronal type diversity of the Drosophila medulla.
We also identified the regulatory interactions between these tTFs. Importantly, we show that not all tTFs function in the same way for the progression of the temporal cascade: while several tTFs directly control the progression by activating the next tTF and which is a unique resource to study these neurons (Table S1).
We also provide a detailed transcriptomic description of the first steps in the differentiation trajectory of a neuron. We find that all optic lobe neurons follow a similar differentiation program once they become postmitotic consisting of four main steps: Finally, we probed a scRNASeq dataset of mouse radial glia for the expression of the optic lobe tTFs. None of the Drosophila tTFs are expressed in strict temporal windows; only Pax6 was found to be expressed in a gradient being enriched in older radial glia, while Foxg1 was slightly elevated in younger ones. Notably, the regulatory interaction between Ey/Pax6 and Slp/Foxg1 is potentially conserved between flies and mice 54 .
Interestingly, the RNA binding protein Imp/Igf2bp1 is expressed in a gradient, being enriched in young progenitors in both Drosophila central brain neuroblasts and in the mouse cortical radial glia 55,56 (Extended Data Figure 7B). Moreover, Bach2 (which has been suggested to be an ortholog of Chinmo 57 ) is expressed in very young radial glia, but more importantly in neurons that come from the high-Igf2bp1 radial glia (Extended Data Figure 7B), This is reminiscent of the expression of its reported ortholog, Chinmo, that is also expressed in young neurons of the fly central brain and mushroom body under the regulation of Imp. This suggests that a similar temporal program may exist between radial glia and Drosophila neuroblasts that regulates neuronal identity.
The absence of strict temporal windows in mouse cortical apical progenitors was also illustrated by a recent paper that identified shared tTFs in the mouse CNS.    (D) Using this annotation, we were able to identify cell types that belong to the four optic lobe neuropils (lamina, medulla, lobula, and lobula plate), as well as central brain cells that were not removed during the dissections, and glial cells.   (M) In negatively marked D mutant clones (GFP-), Bar-H1 expression is lost but Scro expression is not affected, consistent with our observations that Scro is expressed until the very end of the temporal series (see Fig. 2). Previously published results showed that D activates Tll while Tll is sufficient but not necessary to inhibit D.
Scale bar: 10um  Library preparation and sequencing

Single-nucleus RNA seq
Human cortical plate sample preparation Tissue was collected from de-identified prenatal autopsy specimens without neuropathological abnormalities under approved IRB protocol. The cortical plate was dissected fresh from the anterior frontal lobe of anatomically intact brain specimens with postmortem time interval less than 24 hours, and immediate fresh-frozen on dry ice.
Isolation and fluorescence-activated nuclear sorting (FANS) with hashing.
All buffers were supplemented with RNAse inhibitors (Takara). 25mg of frozen postmortem human brain tissue was homogenized in cold lysis buffer (0.32M Sucrose, 5 and 83.6%, respectively. We therefore chose a dimensionality of 150 for the larval dataset. The dataset was then clustered with a resolution of 2. Notably, in this developing structure, cells are clustered both by identity and by differentiation stage. For example, Mi1 cells fall into 2 clusters, an immature (cluster 23) and a mature cluster (cluster 53).
Larval and pupal datasets were merged using default parameters.150 PCs were used subsequently for generating the UMAP to remain consistent with the integration of the larval dataset.

Spatial patterning analysis
To focus on the heterogeneity within the neuroepithelial cells, the larval dataset is further Trajectory analysis: identification of candidate tTFs To study temporal patterning in neuroblasts, we first identified the cluster that corresponded to the medulla neuroblasts (cluster 9) based on the expression of Dpn, as well as the expression of the known temporal factors. We extracted the counts from these cells and inputted them into Monocle. We used default parameters to order the cells in pseudotime. We used the DDRTree method for dimensionality reduction. The cells were then ordered in pseudotime and the beginning and end of the trajectory were defined based on the expression of the known tTFs (i.e. Hth marked the beginning of the trajectory and Tll marked the end). We then looked at the expression along the pseudotime of 629 genes annotated as transcription factors in FlyBase to identify the candidate tTFs.
Merging of larval and pupal Mi1 and DE analysis over pseudotime Larval and pupal (P15, P30, P40, P50, and P70) datasets were merged after cells were batch corrected for each stage separately. Standard Monocle workflow was followed to generate trajectories. The L3 and P15 trajectories were ordered manually.
Based on the way the developing optic lobe develops, there are cells at the same differentiation stage in the L3 and P15 datasets. We, therefore, decided to align these two datasets in order to get a continuum of expression. We tested different genes and ended up using "Ggamma30A" as a reference gene. Ggamma30A starts increasing in the middle of the L3 trajectory and continues all the way to P15 in linear manner. We adjusted the expression of Ggamma30A in P15 using linear regression, which was then applied to all genes of P15. This does not change the dynamics of expression, just the relative levels, and serves the purpose of aligning the trajectories over pseudotime of L3 and P15.
To identity differentially expressed genes along the differentiation trajectory from L3 to P70, we used two methods: "principal graph" and "knn". We selected genes that were identified as differentially expressed with at least one of the two methods. We then used the find_gene_modules function to group the differentially expressed genes into modules of genes that co-vary. These genes were then used for GO analysis.

GO enrichment analysis
We performed GO enrichment analysis and calculated enrichment for 'Biological Process' using The Gene Ontology Resource (http://geneontology.org/) using a Fisher's exact test to calculate p-value. Multiple testing correction was performed by calculating the False Discovery Rate.
To find the expression of GO terms over time, we added and normalized the expression of all genes that belong to a specific GO term and plotted it over pseudotime or on the UMAP.

Analysis of human data
We mapped the sequenced libraries to the H. sapiens genome assembly GRCh38 We used the find_gene_modules function to group genes into 6 modules of genes that co-vary. These modules were then used for GO analysis.

Analysis of mouse cortical data
The dataset that was generated by Telley et al. 18 was downloaded from GEO (GSE118953). The raw counts were inputted into Seurat and the standard workflow was followed (log-normalization, followed by clustering and UMAP using 25 PCs, and clustering was done with a resolution of 2). The radial glia clusters (clusters 2 and 3) were identified based on the expression of known radial glia markers, such as SOX2 and PAX6. Radial glia from different embryonic days 12, 13, 14, and 15 were used to generate the violin plots of Extended Data Figure 7A.