Temporal analysis of enhancers during mouse cerebellar development reveals dynamic and novel regulatory functions

  1. Miguel Ramirez
  2. Yuliya Badayeva
  3. Joanna Yeung
  4. Joshua Wu
  5. Ayasha Abdalla-Wyse
  6. Erin Yang
  7. FANTOM 5 Consortium
  8. Brett Trost
  9. Stephen W Scherer
  10. Daniel Goldowitz  Is a corresponding author
  1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Canada
  2. University of British Columbia, Canada
  3. RIKEN, Japan
  4. The Centre for Applied Genomics, The Hospital for Sick Children, Canada

Abstract

We have identified active enhancers in the mouse cerebellum at embryonic and postnatal stages which provides a view of novel enhancers active during cerebellar development. The majority of cerebellar enhancers have dynamic activity between embryonic and postnatal development. Cerebellar enhancers were enriched for neural transcription factor binding sites with temporally specific expression. Putative gene targets displayed spatially restricted expression patterns, indicating cell-type specific expression regulation. Functional analysis of target genes indicated that enhancers regulate processes spanning several developmental epochs such as specification, differentiation and maturation. We use these analyses to discover one novel regulator and one novel marker of cerebellar development: Bhlhe22 and Pax3, respectively. We identified an enrichment of de novo mutations and variants associated with autism spectrum disorder in cerebellar enhancers. Furthermore, by comparing our data with relevant brain development ENCODE histone profiles and cerebellar single-cell datasets we have been able to generalize and expand on the presented analyses, respectively. We have made the results of our analyses available online in the Developing Mouse Cerebellum Enhancer Atlas, where our dataset can be efficiently queried, curated and exported by the scientific community to facilitate future research efforts. Our study provides a valuable resource for studying the dynamics of gene expression regulation by enhancers in the developing cerebellum and delivers a rich dataset of novel gene-enhancer associations providing a basis for future in-depth studies in the cerebellum.

Editor's evaluation

This manuscript is a valuable study to understand regulatory elements in the early developing mouse cerebellum. The authors have generated important ChIP-seq datasets from embryonic and early postnatal mouse cerebellum. The authors use these data to convincingly identify enhancers in the developing cerebellum. The authors have also made their data and analyses available online in the "Developing Mouse Cerebellum Enhancer Atlas."

https://doi.org/10.7554/eLife.74207.sa0

Introduction

Neuronal development is a complex and dynamic process that involves the coordinated generation and maturation of countless cell types. For cerebellar granule cells, the most numerous neurons in the brain, neuronal differentiation consists of several steps beginning with the commitment of neural stem cells to become specified neural precursors, followed by multiple migratory stages to reach and mature at its final destination (Consalez et al., 2020). Underpinning these events is the expression of gene regulatory networks that drive dynamic molecular processes required for proper brain formation (Ziats et al., 2015). However, the transcriptional mechanisms that precisely regulate these gene expression programs have not been fully described.

Gene expression is typically activated when transcription factors (TFs) bind to non-coding regulatory elements and recruit the necessary components to begin transcription. Among the several classes of non-coding sequences that regulate gene expression, enhancers are the most common, with thousands predicted to coordinate transcriptional regulation during development (Heinz et al., 2015). Enhancers are regulatory sequences that serve as binding sites for TFs and activate distal target gene expression. In the brain, enhancers help to ensure that gene expression is spatially- and temporally-specific, defining what genes will be active during distinct stages of development (Nord and West, 2020). Transcriptional regulation by enhancers has been shown to be critical for cellular identity, maturation during central nervous system (CNS) development, and activity-dependent responses in mature neurons (Frank et al., 2015; Pattabiraman et al., 2014). A detailed understanding of the enhancers that govern changes in gene expression during embryonic and early postnatal cerebellum development remains limited. Profiling genome-wide enhancer activity at different time points and identifying their gene regulatory targets can provide insight into developmental processes regulated by enhancer elements.

Several molecular properties have been associated with enhancer activity, and the advancement of sequencing technology has facilitated their identification genome-wide in several developing brain structures (Carullo and Day, 2019). Enhancers are marked with histone post-translational modifications H3K4me1 and H3K27ac, both of which contribute to opening chromatin for TF binding (Calo and Wysocka, 2013). H3K27ac delineates active from poised elements and has been a reliable marker for enhancer activity genome-wide (Creyghton et al., 2010). Analysis of these marks, in conjunction with transcriptomic and epigenomic datasets, has revealed that the vast majority of non-coding variants associated with neurological and psychiatric disorders are found within these regulatory elements, highlighting their importance in functional readout in the brain (Barešić et al., 2020). Thus, profiling enhancer-associated histone modifications in the brain across time provides a comprehensive understanding of gene-regulatory principles, disease-associated variants, and the genetics of brain development (Nott et al., 2019).

The cerebellum has been a long-standing model to study the developmental genetics of the brain. This is, in part, due to the limited number of cell types, well-defined epochs of development for these cell types and a simple trilaminar structure in which these cells are organized, making for an enhanced resolution of events in time and space (Wang and Zoghbi, 2001). Our lab has previously generated transcriptomic time course datasets spanning embryonic and postnatal development using microarrays and Cap Analysis of Gene Expression followed by sequencing (CAGE-seq; Forrest et al., 2014; Ha et al., 2019), which led to the discovery of novel TFs critical for proper development. More recently, the developing cerebellum has served as an ideal setting for pioneering mouse and human single-cell RNA-seq time courses (Aldinger et al., 2021; Carter et al., 2018; Peng et al., 2019; Wizeman et al., 2019). These studies have provided an unbiased classification of cerebellar populations using their associated gene expression profiles and insight into the cellular origins of distinct postmitotic cell subtypes from neural progenitors.

When considering the non-coding regulatory elements, such as enhancers, that fine tune the spatial and temporal expression of these genes, there are only a few studies that have profiled these sequences genome-wide and specific to the developing cerebellum. Enhancer-associated histone marks (H3K27ac and H3K4me1) have been identified in postnatal and adult mouse cerebella (Frank et al., 2015). During embryonic development, these marks have been profiled in the mouse hindbrain, but not specifically in the cerebellum (Gorkin et al., 2020). Recently, the open chromatin landscape has been examined during embryonic and postnatal cerebellum development using single-nuclear ATAC-seq (snATAC-seq) (Sarropoulos et al., 2021). While this analysis has provided a comprehensive atlas of predicted regulatory sequences at cell-type resolution, additional signals associated with enhancer activity, such as H3K4me1 and H3K27ac, have yet to be quantified during embryonic and early postnatal development. Profiling these chromatin marks can provide further evidence of regulatory activity and is an important step in establishing a catalog of active enhancers involved in cerebellar development. Identifying enhancers and predicting their target genes would also provide a valuable resource for the research community and would facilitate the discovery of novel genetic drivers of the precisely-timed and cell-specific molecular events in the developing cerebellum.

We utilize chromatin immunoprecipitation followed by sequencing (ChIP-seq) of enhancer associated histone marks H3K4me1 and H3K27ac at 3 stages of embryonic and early postnatal cerebellar development. We identify temporally specific enhancers using a differential peak analysis comparing postnatal and embryonic timepoints. TF motif enrichment and prediction of gene targets led to the elucidation of molecular processes regulated by enhancers during these stages. As examples of the use of this data for discovery, we identify two novel regulators of cerebellar development, Pax3 and Bhlhe22: a novel marker of GABAergic progenitors and a regulator of postnatal granule cell migration, respectively. Finally, we also demonstrate how this data may be used to explore the role of enhancers in neurodevelopmental disabilities by identifying an enrichment of autism spectrum disorder (ASD) associated SNPs and de novo variants found in ASD-affected individuals in cerebellar enhancers. We have made the results of our analyses available online in the Developing Mouse Cerebellum Enhancer Atlas which can be easily queried, curated and exported by the research community. This provides a rich resource that can be utilized to discover novel genetic regulators of cerebellar development and contribute to understanding the impact of genetic variation associated with neurodevelopmental disorders.

Results

Enhancer identification during cerebellar development

To identify enhancers active during embryonic and postnatal cerebellar development, we generated genome-wide H3K27ac and H3K4me1 ChIP-seq profiles from mouse cerebella dissected at embryonic day 12 (E12), postnatal day 0 (P0) and postnatal day 9 (P9) (Figure 1A). These stages represent a window of time that spans important neurodevelopmental events such as cell specification, emergence from the cell cycle, differentiation, migration, and maturation for all of the major neuronal types of the cerebellum including the cerebellar granule cells, Purkinje cells, cerebellar interneurons, and cerebellar nuclear neurons (Goldowitz and Hamre, 1998; Hatten and Heintz, 1995; Wang and Zoghbi, 2001). H3K27ac and H3K4me1 signals were reproducible between biological replicates as exemplified in a region on chromosome 14 (Figure 1B). There was a high correlation between replicates for both marks at each age (Figure 1—figure supplement 1A). Therefore, we had confidence in using our H3K27ac and H3K4me1 data in downstream analyses. Robust cerebellar enhancers were identified by the presence of overlapping peaks between the two enhancer-associated histone marks at each age. This highlighted a total of 9,622 peaks; 5,859, 474, and 3,289 peaks that were in both the H3K27ac and H3K4me1 datasets at E12, P0, and P9, respectively (Figure 1C). Duplicate peaks between ages were removed, producing a list of 7024 active cerebellar enhancers derived from overlapping H3K27ac and H3K4me1 signals (Supplementary file 1).

Figure 1 with 1 supplement see all
Enhancer identification during cerebellar development.

(A) An overview of the stages of cerebellar development profiled in this study. The datasets collected at these ages and the downstream analyses are shown in the flow chart. Labels: NE: Neuroepithelium, RL: Rhombic lip, EGL: External granular layer, PL: Purkinje layer, IGL: Inner granular layer, CN: Cerebellar nuclei (B) A region of the mouse genome chr14:122,715,876–122,899,964 (mm9) in the Integrative Genomics Viewer (IGV) showing H3K27ac and H3K4me1 profiles across biological replicates of E12, P0, P9 cerebella. Active cerebellar enhancers are highlighted (gray box). (C) Venn diagrams displaying overlap between H3K27ac and H3K4me1 peaks at E12, P0 and P9. (D–E) An example of a cerebellar enhancer identified from the E12 cerebellum. Shown is normalized H3K27ac and H3K4me1 signal at the enhancer (gray box), as well as (E) normalized CAGE-seq expression of the nearest gene, Ascl1, across developmental time, at E12, P0, P9. TPM, Transcripts Per Million. (F–G) An example of a cerebellar enhancer identified from the P9 cerebellum. Shown is normalized H3K27ac and H3K4me1 signal at the enhancer (gray box), as well as (G) normalized (TPM) CAGE-seq expression of the nearest gene, Pax6, across developmental time, at E12, P0, P9.

The relationship between enhancer activity and genes relevant to cerebellar development is shown in genomic regions flanking Ascl1 and Pax6, two genes critical to cerebellar development (Kim et al., 2008; Yeung et al., 2016). We identified an enhancer active at E12 located in close proximity to Ascl1 (Figure 1D). A decrease in the H3K27ac ChIP-seq signal at this enhancer corresponded to a decrease in Ascl1 gene expression (Figure 1E). We identified two active enhancers at P9 located near Pax6 (Figure 1F). H3K27ac ChIP-seq signal also showed a pattern of activity similar to Pax6 expression, increasing from embryonic to postnatal ages (Figure 1G). These results provide validation for the enhancers identified in our dataset relative to genes critical to cerebellar development.

We compared our list of robust cerebellar enhancers to four previously published cerebellar enhancer datasets. First, P7 H3K27ac ChIP-seq and DNase-seq profiles previously generated by Frank et al., 2015 were overlapped with robust cerebellar enhancers. Greater than 90% of our reported robust cerebellar enhancers are replicated by H3K27ac and DNAse-seq peaks from this study (Figure 1—figure supplement 1B-C). Second, we compared the number of robust cerebellar enhancers and our histone profiles identified at P9 to H3K27ac and H3K4me1 ChIP-seq datasets from the adult cerebellum (Gorkin et al., 2020). We found that 73% (7037/9545) of H3K27ac peaks, 59% (20674/35010) H3K4me1peaks and 60% (1974/3289) of robust cerebellar enhancers overlapped with peaks identified in the adult cerebellum (Figure 1—figure supplement 1D). This indicates that a subset of the enhancers identified at P9 are active specifically during postnatal development when compared to adult stages. Third, enhancers retrieved from the enhancer database EnhancerAtlas 2.0, reporting enhancer activity in the mouse cerebellum at P0-P14 (Gao and Qian, 2020), were compared to robust cerebellar enhancers. We found that 73%, and 80% of our enhancers overlapped with the postnatal cerebellum enhancer dataset at P0, and P9, respectively (Figure 1—figure supplement 1E). Fourth, mouse enhancers that had experimentally validated hindbrain activity at E11.5 from the VISTA Enhancer Browser (Visel et al., 2007) were compared to the cerebellar enhancers reported here. We found that 56% of VISTA enhancers overlap with our cerebellar enhancer sequences at E12 (Figure 1—figure supplement 1F). Additionally, we found that a higher proportion of H3K27ac peaks overlapped with VISTA enhancers compared to open chromatin regions identified in the developing cerebellum (Sarropoulos et al., 2021) and the developing hindbrain (Gorkin et al., 2020; Figure 1—figure supplement 1G). This indicates that our identified H3K27ac peaks may be more predictive of active enhancers than previously generated chromatin accessibility datasets. Taken together, these confirmative findings indicate our approach was effective in capturing active cerebellar enhancers.

Enhancer dynamics during cerebellar development

The dynamics of enhancer activity over cerebellar development were examined through a differential peak analysis of H3K27ac signal. The majority, 89% (6238/7023), of cerebellar enhancers had significant differences in peak signal (adjusted P-value ≤0.05) throughout cerebellar development (Figure 2A, Supplementary file 2). At P9, 1273 cerebellar enhancers were significantly active compared to either P0, E12 or both. At E12, 4432 active enhancers were differentially active compared to either P9, P0 or both. At P0, in contrast, only a small number of enhancers with differential signal was identified (403 and 154 showed significant changes when compared to E12 and P9, respectively). However, none of these P0 cerebellar enhancers were differentially active when compared to both E12 and P9, indicating that enhancer activity did not spike at birth. Quality control metrics evaluating sensitivity of the H3K27ac and H3K4me1 ChIPs (Supplementary file 3), such as the fraction of reads in peaks (FRiP) and relative strand correlation (RSC), all met ENCODE standards (Landt et al., 2012). The P0 samples had slightly lower scores, on average, in both metrics than the E12 and P9 samples, suggesting a minor decrease in sensitivity compared to the other stages. We acknowledge that the difference in sensitivity may have masked P0-specific enhancer activity in our analysis and that future experiments may reveal temporally-specific enhancers at this stage. Taken together, this analysis highlights two temporally specific windows of enhancer activity at Early (embryonic) and Late (postnatal) stages (Figure 2B).

Figure 2 with 4 supplements see all
Enhancer activity is dynamic throughout cerebellar development.

(A) Volcano plots showing robust cerebellar enhancers with differential H3K27ac peak signal for three comparisons: E12 vs P9, E12 vs P0, and P0 vs P9. Differential signal strength was identified for 4433 and 4355 robust cerebellar enhancers when comparing E12 to P9 and to P0, respectively. At P9, 1275 and 403 robust cerebellar enhancers had differential signal when compared to E12 and P0, respectively. Enhancers with significant differential activity are colored at a cutoff of an adjusted P-value <0.05. Displayed on the y-axis is the negative log10 adjusted p-value and on the x-axis is the difference in ChIP-seq signal between to the ages for a given peak. (B) Diagram displays how Early and Late active cerebellar enhancers are classified based on differential peak analysis results. (C) Boxplot shows mean ChIP-seq signal (y-axis) for all Early (upper) and Late (lower) active enhancers. Error bars represent the standard error of the mean. (D) Mean profile and heatmaps of H3K27ac signal at the midpoint of our predicted cerebellar enhancers (rows ±3 kb) in Early and Late groups at E12, P0, and P9.

Distinct patterns of enhancer activity were observed for temporally classified enhancers. For Early active enhancers, there was a loss of mean H3K27ac signal over time, with a steep decline after E12 (Figure 2C). Late active enhancers exhibited a gain in activity over time, with mean H3K27ac signal increasing steadily through development. These patterns are seen when looking at the changes in signal flanking the summits of our peaks across time (Figure 2D). These results indicate that the majority of cerebellar enhancers are dynamic throughout time and exhibit temporally specific activity.

A subset of cerebellar enhancers is active in other developing regions of the brain

To examine cerebellar enhancer activity in the context of other developing brain regions, we overlapped our histone profiles and robust cerebellar enhancer coordinates with H3K4me1 and H3K27ac profiles previously generated from either the developing mouse hindbrain, midbrain or forebrain at E12 and P0 (Gorkin et al., 2020). At E12, we found that the majority of cerebellar H3K27ac and H3K4me1 peaks overlapped with profiles generated in either the hindbrain, midbrain, or forebrain (Figure 2—figure supplement 1A). Similar results were observed at P0, with the majority of cerebellar H3K27ac peaks identified in other brain regions; however, only a modest overlap was found for H3K4me1 profiles at this stage (Figure 2—figure supplement 1B). We found that most robust cerebellar enhancers overlapped with H3K27ac peaks detected in the hindbrain (89.1%), midbrain (78.6%), or forebrain (78.6%; Figure 2—figure supplement 1C). When looking at multiple brain regions, 69.8% of robust cerebellar enhancers were active in all three developing brain regions (Figure 2—figure supplement 2). Interestingly, the majority of the enhancers active in these brain regions (73.9%) were classified as Early active enhancers in the cerebellum which display temporally specific activity to embryonic development. Importantly, we identified 467 (6.7%) robust cerebellar enhancers specifically active in the cerebellum and not the other brain regions, indicating our analysis has uncovered novel enhancer sequences potentially critical for cerebellar formation (Figure 2—figure supplement 2). Over half of the cerebellar specific enhancers (55.03%) were classified as Late active enhancers with activity peaking during postnatal development, while only 15.92% of these enhancers were Early active enhancers and the remaining 29.05% of enhancers were not differentially active (labeled as ‘Neither’). We then examined whether Early and Late and enhancers displayed similar patterns of activity in other brain regions from E12 to P0. Early enhancers displayed a decrease in average activity over time in the hindbrain, midbrain and forebrain (Figure 2—figure supplement 3). Late enhancers showed an increase in hindbrain samples, but minimal changes in the developing midbrain and forebrain (Figure 2—figure supplement 3). These collective results indicate that Early enhancers may also be active in other developing regions of the brain while Late active enhancers may be more likely to be active specifically in the cerebellum.

Comparison of robust cerebellar enhancers with chromatin accessibility identified in single cells

To gain a more granular view of the spatial activity of robust cerebellar enhancers, we overlapped their coordinates with cis-regulatory elements (CREs) active in the developing cerebellum previously identified by single-nuclei ATAC-seq (snATAC-seq; Sarropoulos et al., 2021). We found that 6342/7024 (90.1%) robust cerebellar enhancers overlapped with CREs with open chromatin conformation. Cell-types were then assigned to robust cerebellar enhancers based on the activity of an overlapping CRE. CREs were previously aggregated into 26 clusters based on activity using an iterative clustering procedure (Sarropoulos et al., 2021). Using these clusters, we were able to assign a predicted cell type to the 6342 robust cerebellar enhancers that overlapped with CREs (Figure 2—figure supplement 4A). The predicted cell types with the largest proportion of robust cerebellar enhancers were granule cells (19.32%), progenitor cells (18.35%), Purkinje cells (11.48%) and multiple early-born neuron types (9.69%; Figure 2—figure supplement 4A). When splitting cerebellar enhancers into Early and Late active groups, we found that the majority of Early enhancers were predicted to be active in progenitor cells, and multiple early born neuron types such as Purkinje cell precursors while Late enhancers were predicted to be active in developing granule cells, which is the dominant population of cells produced during postnatal stages, and other late-born cell types (Figure 2—figure supplement 4B).

We then examined the average open chromatin signal at our established Early and Late robust cerebellar enhancers in progenitor cells and the predominant neuron types: granule cells, Purkinje cells, and interneurons. Early enhancers peaked in accessibility at embryonic stages and steadily declined throughout developmental time for all cell types (Figure 2—figure supplement 4C). Late enhancers gradually increased in accessibility throughout time peaking at postnatal stages for progenitor cells, granule cells and interneurons. In Purkinje cells, Late active enhancers showed a minimal change in accessibility throughout time (Figure 2—figure supplement 4C). These findings indicate that Early and Late enhancers display their respective activity patterns in the context of chromatin accessibility in individual cell-types.

Cerebellar enhancers are enriched for neural transcription factor binding sites in an age-dependent manner

We then sought to identify transcription factors whose activity is dictated by the availability of robust cerebellar enhancers, as many neural lineage-defining factors drive cell commitment in the developing brain through enhancer binding (Elsen et al., 2018; Lindtner et al., 2019). We used HOMER to search for enriched motifs (adjusted p-value <1E-11) in Late and Early active cerebellar enhancers; which were then matched to known transcription factor motifs in the JASPAR database (Heinz et al., 2010). This analysis revealed a distinct set of significantly enriched motifs for Early and Late enhancers (Figure 3A). The TFs with the best matching DNA binding motifs belonged to protein families with known regulatory roles in cerebellar development, serving as validation for our analysis (Figure 3A). Since TF motifs are similar between protein family members, it is possible that Early and Late cerebellar enhancers are bound by TFs in the same protein family as the predicted best match indicated in Supplementary file 4. TFs enriched in the Early active enhancer group show a decrease in expression over time while TFs enriched in the Late active enhancer group show an increase in expression over time. This correspondence between enriched TF expression and enhancer activity provides validation for our findings and indicates the timing of enhancer activity may be dictated by the expression and binding of these enriched TFs.

Figure 3 with 2 supplements see all
Neural transcription factors with known and novel function in the developing cerebellum are enriched in dynamic cerebellar enhancers.

(A) Dot plot displaying significantly enriched (adjusted p-value <1E-11) motifs and the predicted matching transcription factor (TF). Displayed are the results for Early (top) and Late (bottom) active enhancers. TFs with an unknown functional role in cerebellar development are indicated with a red arrow. Size of the dots indicate the negative log10 adjusted p-value for a given motif and the color scale displays the z-score normalized expression throughout the cerebellar developmental time course. (B) Top: Immunofluorescent staining of Pax3 in the mouse cerebellum at E12, E15, P0, and P3. Bottom: Pax3 and Ptf1a immunofluorescent co-staining of the E12 mouse cerebellum. Immunofluorescent co-staining of Pax3 and Pax2 in the mouse cerebellum at E15 and P0. Labels: CP: Cerebellar parenchyma, EGL: External granular layer, NGL: Nascent granular layer, RL: Rhombic lip, VZ: Ventricular zone, Scalebars = 100 µm.

The top three enriched TF motifs for Early active enhancers were Ascl1, Meis2, and Atoh1 (Figure 3A). These TFs have established roles in cerebellar development, acting as markers of GABAergic or glutamatergic cell types and regulators of differentiation (Ben-Arie et al., 1997; Kim et al., 2008; Wizeman et al., 2019). Importantly, many of the motifs enriched in the Early group matched with TFs which have received little to no attention in the cerebellum, including Sox4, Lhx2, Rfx4, Pou3f1, and Pax3 (Figure 3A). These TFs have been previously associated, however, with the development of other brain areas (Frantz et al., 1994; Porter et al., 1997; Su et al., 2016; Zhang et al., 2006). In contrast, the TFs matching the motifs enriched in the Late active enhancers have a previously identified role in cerebellar development; but not necessarily in the same cellular processes (Figure 3A). For example, the top 3 enriched motifs matched with Neurod1, Nfia/b/x, and NF1, which have all been associated with granule cell differentiation (Miyata et al., 1999; Sanchez-Ortiz et al., 2014; Wang et al., 2007). However, two other TFs with enriched binding sites, Pax6 and Smad4 have been found to be critical for granule cell precursor proliferation, a process preceding differentiation (Swanson and Goldowitz, 2011). These results suggest a dynamic role for the majority of our Early and Late active enhancers, driven by TFs involved in distinct stages of neuron development.

Early active enhancers are enriched for Pax3 binding sites, a novel marker for GABAergic cells

The TF motif enrichment analysis of Early enhancers led to the discovery of several novel TFs in the context of embryonic cerebellar development; potentially involved in seminal aspects of development such as cellular specification or commitment. As a case study, we focused on Pax3, as other members of the Pax protein family have been shown to play key roles in the developing cerebellum (Leto et al., 2009; Urbánek et al., 1997; Yeung et al., 2016). Immunofluorescent staining was conducted to profile Pax3 expression in the embryonic and early postnatal mouse cerebellum. We observed robust expression in the ventricular zone (VZ); a neural progenitor region for GABAergic cells in the cerebellum (Leto et al., 2006; Figure 3B). A molecular marker of the VZ is Ptf1a, a GABAergic lineage-defining molecule in the cerebellum (Hoshino et al., 2005), and we examined if there was co-expression between Pax3 and Ptf1a. Colocalization between 50% of Pax3 positive cells in the VZ and Ptf1a (Figure 3—figure supplement 1A), (Hoshino et al., 2005) confirmed Pax3 expression within GABAergic neural progenitors. At E15, Pax3 + cells are seen in the region just dorsal to the VZ, which consist of post-proliferative cells such as Purkinje cells and interneurons (Hoshino et al., 2005; Leto et al., 2006; Figure 4B). We examined Pax3 co-labeling with markers for these cell types using Foxp2 (marker of Purkinje cells) and Pax2 (a general marker of GABAergic interneurons; Fujita et al., 2008; Maricich and Herrup, 1999). While colocalization between 91% of Pax3 + cells and Pax2 was found (Figure 3B, Figure 3—figure supplement 1A), no co-staining between Pax3 and Foxp2 was observed (Figure 3—figure supplement 1B). These results extend to P0 and P3, where Pax3 + cells are found in the inner granule cell layer (IGL) as well as the cerebellar parenchyma (Figure 4B, Figure 3—figure supplement 1A). To identify if Pax3 was expressed in interneuron precursors at these stages, co-labeling was conducted with Pax2. We found that 76% of Pax3 + cells co-localized with Pax2 and P0 and that 79% of Pax3 + cells co-localized with Pax2 at P3 (Figure 4B, Figure 3—figure supplement 1A). Pax3 did not colocalize with Calbindin, a Purkinje cell marker, at P0 (Figure 3—figure supplement 1B). At P9, Pax3 exhibited no observable expression (Figure 3—figure supplement 1C).

Figure 4 with 2 supplements see all
Correlated Early target genes are expressed in spatially distinct areas and have diverse roles in cerebellar development.

(A) Line plot and heatmap showing mean z-score expression for Early target genes throughout the cerebellar time course. (B) Line graph representation of expression pattern throughout time for each cluster. (C) Known cerebellar genes in each cluster and in situ hybridization (ISH) images showing spatial expression at peak expression ages. ISH images were taken from the Developing Mouse Atlas at E13.5 for clusters 1 and 2, and E11.5 for clusters 3 and 4. (D) Gene Ontology (GO) enrichment analysis of target genes from each cluster, displaying the top enriched GO terms. Size of the dots indicates the observed vs expected fold enrichment of genes within that GO category (FoldEnrichment). Color scale indicates the adjusted p-value for each GO term.

We then examined Pax3 expression throughout developmental time using a bulk tissue CAGE-seq dataset generated previously (Ha et al., 2019). Pax3 shows high expression during embryonic stages, peaking at E12, and decreasing steadily after E12 (Figure 3—figure supplement 2A). To further validate the cell types in which Pax3 are expressed, previously generated single-cell RNA-seq datasets quantified in the developing mouse and human cerebellum were examined (Carter et al., 2018; Aldinger et al., 2021). Pax3 showed similar expression profiles in the developing mouse and human cerebellum and is highly expressed in Lhx1/5+GABAergic progenitors and Pax2 +GABAergic interneurons supporting our immunofluorescent analysis (Figure 3—figure supplement 2B-C). Pax3 expression was also observed in other cell lineages, such as rhombic lip progenitors, UBC/Granule cell progenitors and astrocytes. We do not observe Pax3 expression through immunofluorescent staining in these cell types, which suggests a misalignment with single-cell data. This highlights the need to confirm the fidelity between histological and single-cell datasets. Taken together, our results indicate that Pax3 is a novel marker for GABAergic progenitors and interneuron precursors in the developing cerebellum.

Co-expressed putative target genes are expressed in spatially distinct areas of the developing cerebellum

We next investigated the molecular processes regulated by robust cerebellar enhancers through predicting their downstream targets (Osterwalder et al., 2018; Yao et al., 2015). This was done by calculating the correlation between H3K27ac signal and gene expression at E12, P0, and P9 (Ha et al., 2019) for genes located within the same conserved topological associating domain (TAD) identified previously (Dixon et al., 2012) (See Materials and methods). Overall, at least one positively correlated target gene was identified for 5815/7023 (70.61%) cerebellar enhancers with an average Pearson correlation coefficient of 0.86 (Supplementary file 5). In total, we identified 2261 target genes. Using the Mouse Genome Informatics (MGI) database, we identified 98 target genes that when knocked out result in a cerebellar phenotype. To evaluate whether cerebellar genes were enriched in putative target genes, we conducted a permutation analysis by generating 10,000 permutations of 2261 genes randomly selected from all genes expressed in the cerebellum and assessing the number of cerebellar genes in each permutation. We found that cerebellar genes were indeed enriched in target genes (p-value = 0.0405, Figure 4—figure supplement 1A) demonstrating the validity of our high-throughput approach.

An unbiased k-means clustering was then conducted for Early and Late target genes to delineate them into the various co-expression programs coordinating molecular events during development. For this analysis, the target gene expression time course was expanded to 12 different timepoints during cerebellar development, quantified previously by CAGE-seq (Ha et al., 2019). To determine the k-value, we conducted an Elbow analysis, identifying that 4 clusters were optimal for Early and Late active enhancers (Figure 4—figure supplement 1B).

For Early active enhancers, 4 Clusters of co-expressed target genes were identified (Figure 4A). Genes in these clusters had decreasing expression over time, similar to their corresponding enhancer activity. However, a distinct mean expression profile was observed for each Cluster (Figure 4B). Interestingly, genes with known function during cerebellar development showed distinct spatial expression patterns, observed using ISH data from the Developing Mouse Brain Atlas (Thompson et al., 2014; Figure 4C).

For example, in Cluster 3, cerebellar genes Ascl1 and Neurog2 are expressed exclusively in the ventricular zone at E11.5 while Cluster 4 contains Lhx9 and Meis2 which are expressed in the Nuclear Transitory Zone (neurons destined for the cerebellar nuclei). A Gene Ontology (GO) enrichment analysis revealed that each cluster is enriched for molecular processes known to be regulated by cerebellar genes within the cluster (Figure 4D). For example, Cluster 1 is enriched for axonogenesis (GO:0007409, p-value: 3.31E-4), neuron migration (GO:0001764, p-value: 3.3E-4) and Purkinje layer development (GO:0021691, p-value: 0.01) and also contains Lhx1 and Lhx5 which are expressed in migrating Purkinje cells in cerebellar parenchyma and has previously been associated with the regulation of Purkinje cell differentiation during embryonic cerebellar development (Zhao et al., 2007). Together, these findings support the notion that Early active enhancers regulate their targets in a spatially specific manner, regulating distinct processes in their respective cell types.

We then used the results of our putative target gene analysis of Early active cerebellar enhancers to better understand the potential function of Pax3 during embryonic cerebellar development. We determined the putative target of genes of robust cerebellar enhancers containing Pax3 DNA binding motifs. We found that 1111 Early active robust cerebellar enhancers contained a Pax3 binding motif and that a putative target gene was identified for 923 Pax3-motif enhancers (Supplementary file 6). These target genes were enriched for GO terms pertaining to neural progenitor function and specification (Figure 4—figure supplement 2). Several putative Pax3 target genes have been previously associated with GABAergic neural progenitors in the cerebellar ventricular zone and the specification and differentiation of GABAergic cell types such as Ascl1 and Neurog2 (Florio et al., 2012; Kim et al., 2008). These genes also display a similar spatial expression pattern to Pax3 at E11, with high expression in the cerebellar ventricular zone (Figure 4—figure supplement 2).

For the Late active enhancers, 4 Clusters of co-expressed target genes were identified (Figure 5A). We observed relatively distinct expression patterns in each of the 4 Clusters with a gradual rise in mean expression over time (Figure 5B). Genes with known function during cerebellar development also show distinct spatial expression patterns, identified using the Developmental Mouse Atlas (Figure 5C). For example, Clusters 1 and 3 contained known cerebellar genes critical for granule cell development, such as Neurod1 and Zic1, while Clusters 2 and 4 contained cerebellar genes important for in Purkinje cell development, such as Atxn1 and Hcn1 (Figure 5C; Aruga and Millen, 2018; Ebner et al., 2013; Miyata et al., 1999; Rinaldi et al., 2013).

Correlated Late target genes are expressed in developing granule cells or Purkinje cells with common roles in cerebellar development.

(A) Line plot and heatmap showing mean z-score expression throughout the cerebellar time course. (B) Line plot representation of expression pattern throughout time for each cluster. (C) Known cerebellar genes in each cluster and in situ hybridization showing spatial expression at peak expression ages. ISH images were taken from the Developing Mouse Atlas provided by the Allen Brain Atlas conducted at P4.5 for all clusters. (D) Gene Ontology (GO) enrichment analysis of all target genes of Late active enhancers, displaying the top enriched GO terms. Size of the dots indicates the observed vs expected fold enrichment of genes within that GO category (FoldEnrichment) and the x-axis represents the -log10 adjusted p-value for each GO term. Color scale indicates the adjusted p-value for each GO term.

A GO enrichment analysis was conducted for each Cluster; with no significantly enriched Cluster-specific GO terms. However, if all Late enhancer target genes were combined, several enriched GO terms emerged including ones involved in postnatal neuronal development, such as neuron death (GO:0070997, p-value: 0.003), neurotransmitter transport (GO:0006836, p-value: 0.006) and regulation of synaptic vesicle exocytosis (GO:2000300, p-value: 0.005) (Figure 5D). Overall, this analysis provides a working framework for the placement of hundreds of genes into the overall structure of embryonic or postnatal cerebellar development.

Bhlhe22 is a novel regulator of granule cell development

To demonstrate the utility of our results, we sought to elucidate target genes not previously identified in cerebellar development. We focused on Late Cluster 1, which contained target genes expressed in granule cells. We hypothesized that genes within this cluster regulated postnatal granule cell differentiation. To identify genes in this cluster regulating granule cell development, we filtered these genes relative to their interaction with Atoh1, the lineage defining molecule for granule cells and other glutamatergic neurons in the developing cerebellum (Ben-Arie et al., 1997). The genes were filtered using the following criteria: (1) Atoh1 is bound to the predicted enhancer during postnatal development (Klisch et al., 2011) and (2) the genes are differentially expressed in the Atoh1-null mouse (Klisch et al., 2011). The resulting list of genes were sorted based on differential expression in the Atoh1-null mouse and filtered for TFs. These criteria filtered gene candidates for validation from Late Cluster 1 from 254 genes to 26 genes. Among the top 15 genes in the filtered list, we identified 4 novel genes and 11 genes that have previously been implicated in postnatal granule cell development (Table 1).

Table 1
A list of enhancer-regulated target genes from Late Cluster 1 found to be significantly differentially expressed in the conditional Atoh1 knockout mouse.

The second and third column contain the observed P-value and fold change from the differential expression analysis, respectively. The fourth and fifth columns indicate whether the gene has previously been implicated in cerebellar development and the corresponding reference PubMed ID.

Genep-value (Atoh1-null)Fold Change (Atoh1-null)Cerebellar DevelopmentReference (PMID)
Neurod19.196E-2290.2X19609565
Nfix1.1991E-430.53X21800304
Zic11.3082E-370.35X21307096
Barhl12.2505E-350.22X9412507
Zic22.0968E-200.34X11756505
Insm15.4814E-200.25X18231642
Tcf49.9139E-200.69X30830316
Nfia1.8675E-160.6X17553984
Bhlhe224.7662E-100.53
Purb2.5878E-090.53
Neurod24.424E-090.37X11356028
Klf131.7938E-060.72
Zfp5213.5899E-060.8X24676388
Sox183.7168E-050.71
Nfib0.000210090.63X17553984

The known genes provided validation for our approach. The novel genes included Bhlhe22 (also known as Bhlhb5), Purb, Klf13, and Sox18. We focused on Bhlhe22 as it has previously been implicated in the differentiation of neurons in the cortex (Joshi et al., 2008). An enhancer ~2 kb upstream of the Bhlhe22 transcriptional start site was identified and is bound by Atoh1 during postnatal development (Figure 6—figure supplement 1A). This enhancer displayed H3K27ac activity highly correlated (Pearson correlation coefficient = 0.98) to Bhlhe22 expression, which consistently rises throughout cerebellar development and peaks at P9.5 (Figure 6—figure supplement 1B).

To attain a cellular resolution of the expression pattern for Bhlhe22 during postnatal cerebellar development cerebellar development, immunofluorescent staining was conducted. Bhlhe22 expression was observed in cells within the inner external granule layer (EGL), molecular layer (ML) and in the internal granule layer (IGL) of the postnatal cerebellum (Figure 6A). Cells in the ML are Bhlhe22 +may be of two phenotypes: (1) small cell-bodied radially oriented migrating granule cells and (2) larger cell-bodied rounded interneurons.

Figure 6 with 2 supplements see all
Bhlhe22 is expressed in differentiating granule cells in postnatal cerebellar development.

(A) Bhlhe22 (green) and Neurod1 (red) immunofluorescent co-staining at P9.5 of taken from a posterior lobe IX. (B) Bhlhe22 (green) and NeuN (red) immunofluorescence co-staining at P6 taken from posterior lobe IX. (C) Bhlhe22 (red) and Dcx (green) immunofluorescent co-staining at P6 showing the posterior lobe IX; Labels: EGL = external granular layer IGL = internal granular layer, ML = molecular layer, Scalebars = 100 µm.

To identify whether Bhlhe22 is expressed in differentiating granule cells, co-staining experiments were performed with Neurod1 and NeuN which mark differentiating and more mature granule cells, respectively (Miyata et al., 1999; Weyer and Schilling, 2003). At P6.5, colocalization between 98% of Bhlhe22 + cells and Neurod1 was observed, indicating expression in differentiating and migrating granule cells (Figure 6A, Figure 6—figure supplement 1C). Co-staining between 82% of Bhlhe22 + cells and NeuN expression was also observed, indicating expression in maturing granule cells found within the IGL (Figure 6B, Figure 6—figure supplement 1C). To confirm whether the Bhlhe22-positive cells within the molecular layer were migrating granule cells, we performed a double labelling experiment for a neuronal migration marker Doublecortin (Takács et al., 2008). Colocalization of Doublecortin in 95% of Bhlhe22 + cells was observed within the inner EGL and the molecular layer, confirming Bhlhe22 expression in migrating granule cells (Figure 6C, Figure 6—figure supplement 1C). To assess Bhlhe22 expression at a single-cell resolution, we examined Bhlhe22 expression in single-cell RNA-seq datasets previously generated in the developing mouse and human cerebellum (Carter et al., 2018; Aldinger et al., 2021). In developing mouse and human cerebellum, Bhlhe22 is highly expressed in granule cells and granule cell progenitors, supporting the results of our immunofluorescent analysis (Figure 6—figure supplement 2). Bhlhe22 is also expressed in other cell types in the glutamatergic lineage, such as unipolar brush cells and excitatory cerebellar nuclei as well as high expression in inhibitory GABAergic interneurons neurons indicating importance for both neuronal lineages in the cerebellum.

We then investigated the role that Bhlhe22 plays in postnatal granule cell development using a well-established in vitro system (Lee et al., 2009). Three sets of experiments were performed using P6.5 isolated granule cells transfected with siRNA targeting Bhlhe22 transcripts (Figure 7A). First, to determine if Bhlhe22 expression was diminished, changes in gene expression were assessed after 3 days in vitro (DIV) using reverse transcriptase quantitative PCR (RT-qPCR). A 50% reduction of Bhlhe22 expression, on average, was found in treated cultures compared to controls (Figure 7B).

Knockdown of Bhlhe22 reduces migration of cultured cerebellar granule cells.

(A) Workflow for dissociated and reaggregate postnatal granule cell cultures. (B) RT-qPCR analysis of Bhlhe22 gene expression in dissociated postnatal granule cell cultures after treatment with Bhlhe22 siRNA. Gene expression was normalized relative to the expression of the co-transfected EGFP protein to account for transfection variability between cultures. Data are presented as mean ± SD (n=3). (C) Image of cultured cerebellar granule cell reaggregates treated with control and Bhlhe22 siRNA. Shown are EGFP-positive cells indicating successful transfection. Scalebars = 100 µm. (D) Box plot displaying mean distance of granule cell migration from the aggregate. Value above indicates a statistical difference between control cultures and those treated with Bhlhe22 siRNA (p-value = 0.0013). (E) Bar plot showing the percentage of cells migrated at different distances from the aggregate for control and Bhlhe22 siRNA-treated cerebellar granule cell cultures. (F) RT-qPCR analysis of gene expression of cell adhesion molecules in dissociated postnatal granule cell cultures after treatment with Bhlhe22 siRNA. Gene expression was normalized relative to the expression of the co-transfected EGFP protein to account for transfection variability between cultures. Data are presented as mean ± SD (n=3). Symbols: *: p ≤ 0.05, **: p ≤ 0.01, ***: p ≤ 0.001, which indicate statistical differences observed between Bhlhe22 siRNA-treated samples and controls. All error bars represent the standard error of the mean,.

Second, phenotypes of the transfected cells were examined: neuritic outgrowths from the aggregate and the migration of granule cells from the aggregates, within the first 24 hr of plating (Gärtner et al., 2006). Neuritic outgrowth was unaffected. In contrast, a marked reduction in migration was found (Figure 7C). Bhlhe22 siRNA transfected cells travelled on average 54.2 µm from the edge of the aggregate, a 50% reduction compared to control samples (Figure 7D). Furthermore, examining the distribution of migrated cells from the edge of the aggregate identified a higher percentage of Bhlhe22 siRNA transfected granule cells migrating less than 50 µm, while the majority of the cells in control samples migrated 100 µm and beyond (Figure 7E).

Third, changes in the expression of cell adhesion molecules that are known to be involved in granule cell development were assessed (Consalez et al., 2020; Wang et al., 2007). A significant reduction of Efnb1, Efnb2, Tag1, Cdh2, and Astn2 was observed in Bhlhe22 knockdown granule cell cultures compared to controls (Figure 7F). In addition to these genes, we also found a significant reduction in Doublecortin (Dcx) expression. Taken together, these in vitro knockdown experiments reveal a novel function for Bhlhe22, a gene that was identified by our temporal enhancer-target gene analysis and was predicted to have a critical role in postnatal granule cell development.

Active cerebellar enhancers are enriched for common and de novo genetic variants associated with autism spectrum disorder

Given the emerging importance of the cerebellum in the etiology of autism spectrum disorder (ASD) (Limperopoulos et al., 2014; Stoodley and Limperopoulos, 2016) and to demonstrate the utility of our dataset in potentially functionally annotating variation associated with neurodevelopmental disorders, we tested whether ASD-associated variants and de novo mutations are enriched in cerebellar enhancers. Analysis of genome wide association studies (GWAS) have revealed that the majority of variants associated with neurodevelopmental diseases are found within non-coding regulatory sequences, particularly enhancers (Visel et al., 2009). The software tool GREGOR (Genomic Regulatory Elements and Gwas Overlap algoRithm) was used to evaluate the enrichment of common genetic variants associated with ASD in robust cerebellar enhancers (Schmidt et al., 2015). Robust cerebellar enhancers were converted from mouse (mm9) to human (hg38) genomic coordinates (6630/7024, 94.4% converted). The majority (89.6%) of robust cerebellar enhancers with a putative target gene are located at orthologous regions distal to the same gene promoters indicating conservation between mouse and humans and bringing credence to our analysis (Supplementary file 7). ASD-associated SNPs were retrieved from the GWAS Catalog (Buniello et al., 2019) and a stringent filter was applied to identify SNPs associated with the ASD (see Materials and methods). We examined 174 ASD-associated SNPs with a maximum p-value of 9E-06 (Buniello et al., 2019). ASD-associated SNPs were enriched in robust cerebellar enhancers (p-value = 2.34E-03) and in H3K27ac peaks at E12, P0, and P9 (p-values of 1.29E-03, 1.05E-02 and 1.42E-04, respectively) (Figure 8A). As a negative control, we conducted the same analysis with SNPs associated with chronic kidney disease and found no enrichment in cerebellar enhancers. For the 13 cerebellar enhancers containing ASD-associated SNPs, we identified 12 predicted target genes (Supplementary file 8). Among these, three (PAX6, TCF4, and ZMIZ1) are ASD risk genes according to the Simons Foundation Autism Research Initiative (SFARI) gene database (Banerjee-Basu and Packer, 2010).

Figure 8 with 3 supplements see all
Cerebellar enhancers are enriched for GWAS SNPs and DNMs associated with ASD.

(A) Enrichment analysis of ASD-associated and chronic kidney disease associated (negative control) GWAS variants in cerebellar enhancers and H3K27ac peaks called from E12, P0, and P9 samples. (B) Enrichment of de novo single nucleotide variants and indels in ASD-affected individuals compared with their unaffected siblings. Counts are not equal to the sum of the four enhancer types because some enhancers are categorized as more than one type. (C) Gene targets for enhancers overlapped by de novo CNVs in the SSC cohort. (D) Left: Line graph representing Cdc424bpb normalized expression in the developing mouse cerebellum from E11 to P9. TPM = transcripts per million. Right: In situ hybridization showing Cdc42bpb expression in the lateral (left) and medial (right) adult mouse cerebellum (Developing Mouse Brain Atlas). Note expression is found in granule cells, particularly those of the lateral cerebellum.

De novo mutations (DNMs) (variants present in the genome of a child but not his or her parents) have been found to play a significant role in the etiology of ASD, including those found in non-coding regions of the genome (Grove et al., 2019; Yuen et al., 2016). We hypothesized that DNMs within cerebellar enhancers would be more prevalent in ASD-affected individuals compared with their unaffected siblings. We used whole-genome sequencing data from 2603 ASD-affected individuals and 164 unaffected siblings from the MSSNG cohort (C Yuen et al., 2017), as well as 2340 ASD-affected individuals and 1898 unaffected siblings from the Simons Simplex Collection (SSC) (Fischbach and Lord, 2010) to analyze the prevalence of DNMs in ASD-affected individuals compared with their unaffected siblings.

We found that DNMs (specifically de novo single nucleotide variants and indels) in cerebellar enhancers and H3K27ac peaks from E12, P0, and P9 were enriched in ASD-affected individuals, with odds ratios ranging from 1.04 to 1.10 (Figure 8B). While these differences were not statistically significant for cerebellar enhancers and peak coordinates individually, statistical significance was achieved when combined (odds ratio = 1.06; P-value = 0.0043). We also identified de novo CNVs overlapping cerebellar enhancers. Since the number of such CNVs was too small to perform statistical enrichment tests, we selected a subset of seven of these CNVs (four deletions and three duplications) for further characterization to identify candidates for association with ASD (Figure 8C). The most promising candidate was an ~11 kb deletion overlapping an enhancer predicted to regulate CDC42BPB (803 kb upstream of the TSS), which has previously been implicated in neurodevelopmental phenotypes (Chilton et al., 2020). By visual validation in Integrative Genomics Viewer (Robinson et al., 2011), we verified that this deletion was truly de novo (Figure 8—figure supplement 1). Within these coordinates, there is a proportion of sequence alignment between mouse and human genomes located within a robust cerebellar enhancer on chromosome 12, and this enhancer is located in cis a similar distance from the Cdc42bpb TSS as in humans (833 kb) indicating that enhancer activity occurs at orthologous regions (Figure 8—figure supplement 2A). Cdc42bpb is expressed steadily throughout mouse cerebellar development and is expressed in the granule cell layer of the adult mouse cerebellum (Figure 8D). Expression in granule cells is found in the lateral aspects of the adult cerebellum but not the medial adult cerebellum. In the context of cerebellar development, scRNA-seq data from mice and humans showed relatively consistent expression levels across most cell types, with highest expression in glutamatergic cerebellar nuclei and developing Purkinje cells (Figure 8—figure supplement 3A).

Discussion

In our study on the cerebellum, we performed a novel assessment of enhancer activity through genome-wide profiling of H3K4me1 and H3K27ac at three time points during embryonic and early postnatal stages. These datasets were utilized to define functional enhancer elements with temporally specific activity during these developmental ages. The biological processes under enhancer regulation were described through motif enrichment analysis and target gene prediction, identifying temporally and spatially specific regulatory functions. As a result, our study provides a novel dataset for the developmental biology and neuroscience communities, especially those interested in functionally annotating enhancers in the context of brain development. We have created an online resource that can be used to access, curate and export our dataset. In addition to genomic coordinates of our cerebellar enhancers, this resource provides the activity patterns (Early or Late) for each enhancer throughout developmental time, as well as putative target genes. Predicted cell-types in which a robust cerebellar enhancer may be active as well as activity patterns in other brain regions are also included in our resource based upon comparisons with previously generated datasets (ENCODE and single-cell chromatin accessibility). As demonstrated in our analysis of Pax3 and Bhlhe22 in the context of the cerebellum, this resource can facilitate the discovery of novel genetic regulators of cerebellar development.

Cerebellar enhancers regulate gene expression important for distinct stages of neuronal development

Identification of enriched TF binding sites and putative target genes indicated that cerebellar enhancers likely play a regulatory role in various phases of neuronal development. In agreement with our results, previous examinations of active non-coding regulatory sequences revealed that neural progenitor cells and mature neurons exhibit distinct signatures of enhancer-associated histone profiles, DNA methylation, chromatin conformation and enhancer-promoter interactions (Bonev et al., 2017; Torre-Ubieta et al., 2018; Lister et al., 2013; Whyte et al., 2012). Bonev et al., 2017 examined changes in enhancer-promoter interactions between transgenic cell lines that were FACS sorted for embryonic stem cells, neural progenitors and mature neurons and identified that changes in enhancer-promoter contacts are cell-state specific and correlate with changes in gene expression (Bonev et al., 2017). A global shift in regulatory sequence usage was observed between neuro-progenitors and mature neurons, indicated by dynamic changes in enhancer-promoter interactions. These changes were also reflected at the level of TF binding, as interactions at Pax6-bound sites, a TF marking neural progenitors, were stronger in neural progenitors than in neurons, while NeuroD2-bound sites, a TF marking mature neurons, were stronger in neurons than NPCs (Bonev et al., 2017). This shift in enhancer usage throughout cortical development is also reflected in DNA methylation profiles, where fetal enhancers are hypermethylated and decommissioned in the adult brain, while enhancers regulating adult gene expression were hypomethylated (Lister et al., 2013). Hypermethylation was accompanied by a decrease in H3K4me1, H3K27ac and DNase hypersensitivity while the increase was observed after hypomethylation (Lister et al., 2013). Our study supports the importance of temporally-specific activity during different stages of neuron development in vivo and details the processes driven by enhancer-regulated expression during embryonic and early postnatal brain development.

Expression analysis of two genes novel to cerebellar development, Pax3 and Bhlhe22, supported the notion that enhancer profiles are specific to developmental stage. TF enrichment analysis identified Pax3 preferentially enriched in Early active enhancers and robust expression of Pax3 was localized to GABAergic interneuron progenitor cells. Pax3 has previously been associated with neural tube and neural crest development (Epstein et al., 1991; Olaopa et al., 2011); however, the study of it’s function in cerebellar development is limited. Pax3 has been shown to be upregulated by BDNF and NGF in cerebellar in vitro cultures (Kioussi and Gruss, 1994). In the context of human disorders and disease with cerebellar phenotypes, it has been found to be downregulated Group 3 medulloblastomas (Zagozewski et al., 2020) and is located downstream of a deletion on chromosome 2 in Dandy-Walker Syndrome patients (Jalali et al., 2008). Analysis of predicted gene targets of Late active enhancers identified Bhlhe22 as a novel gene expressed in postnatal differentiated granule cells, and in vitro knockdown experiments in primary granule cells indicated Bhlhe22 regulates granule cell migration potentially through regulation of cell adhesion molecule expression. These results are supported by findings in the developing cortex, where Bhlhe22 has been shown to regulate post-mitotic acquisition of area identity in layers II-V of the somatosensory and caudal motor cortices (Joshi et al., 2008). The contrasting expression profiles of Pax3 and Bhlhe22 highlight the wide-ranging developmental impact of enhancer-mediated gene expression regulation.

Co-expressed gene targets of cerebellar enhancers display cell-type-specific expression patterns

In addition to being temporally-specific, recent evidence indicates that enhancer activity is cell type specific in the brain (Blankvoort et al., 2018). In the context of cerebellar development, a recent study examining open chromatin regions genome-wide using snATAC-seq identified a catalog of cerebellar cis-regulatory elements (CREs) and has revealed that the cell types of the cerebellum have unique chromatin signatures throughout embryonic and postnatal development (Sarropoulos et al., 2021). We found that the majority (~90%) of the robust cerebellar enhancers identified in our study overlap these CREs. Additionally, open chromatin signals at these sequences closely resemble the changes in chromatin marks for Early and Late cerebellar enhancers in the major neuronal cell types in the cerebellum. This overlap further validates our findings and provides additional evidence for enhancer activity at these genomic locations. The comparison of our results with CREs identified in single cells allowed the prediction of the cell types in which our robust cerebellar enhancers may be active. These predictions will help guide future studies looking to validate enhancer activity and regulatory function.

Spatial specificity of enhancer activity is also highlighted in the identification of cerebellar enhancer target gene clusters for Early and Late active enhancers with cell specific patterns of expression (Figures 4 and 5). For example, distinct boundaries can be seen in gene expression from Early Clusters 3 and 4 at E11.5 between cells in the subpial stream (Cluster 4) and neuroepithelium (Cluster 3) where neural precursors of two separate lineages, the glutamatergic cerebellar nuclei and GABAergic cerebellar nuclei and Purkinje cells, are found, respectively. These sharp borders are reminiscent of the small domains of distinct enhancer activity identified in neural progenitors in the telencephalon, which were found to fate-map to specific prefrontal cortex subdivisions (Pattabiraman et al., 2014). We see a similar pattern in the more developed postnatal cerebellum, observing a spatial distinction between Late Clusters 1/3 and 2/4 delineating expression in granule cells and Purkinje cells, respectively. This cell-type-specific enhancer usage is demonstrated in the adult brain. Blankvoort et al., 2018 used ChIP-seq analysis of microdissected subregions of the adult mouse cortex to reveal unique enhancer profiles pertaining to each region. Additionally, Nott et al., 2019 identified enhancer-promoter interactome maps specific to the major cell types in the cortex, which included neurons, microglia, oligodendrocyte, and astrocytes. Enriched GO terms for each cerebellar target gene clusters were cell-type and temporally specific, highlighting enhancer specificity. Functionally annotating their respective clusters provides a working hypothesis for hundreds of genes, which can be used as a jumping point for future in-depth studies in the cerebellum. Collectively, these findings support the notion that the cell types in the cerebellum have specific enhancer signatures which are reflected by the expression and functions of their target genes.

The results from our analyses revealed that the majority of the identified cerebellar enhancers are predicted to be active in the most abundant cell types at each developmental stage. For example, most Early active enhancers were predicted to be active in progenitors while Late active enhancers were predicted to be active in the developing granule cells and Purkinje cells. We attribute this apparent bias to our bulk tissue approach, which inherently identifies active enhancers based on signal abundance (ie. peak calling) making it more likely to capture signals specific to more abundant cell types compared to other less abundant and rare cell types. In order to assess histone marks in less abundant cell types, a more granular approach should be implemented such as isolation of certain cell types or single-cell analysis. Recent technological developments have allowed the examination of histone modifications genome-wide in single cells, which has seen preliminary success in mouse brain tissue (Zhu et al., 2021, Bartosovic et al., 2021, Rang et al., 2022).

Cerebellar enhancers display activity in other developing brain regions, especially during embryonic development

The comparison between active cerebellar enhancers with histone profiles generated in other developing brain regions revealed that the majority of Early enhancers were active in multiple brain regions. Previous genome-wide and locus-specific functional studies have shown evidence of enhancers with activity in multiple developmental contexts (Singh and Yi, 2021; Preger-Ben Noon et al., 2018; Lonfat et al., 2014; Andersson et al., 2014; Nord et al., 2013; Hiller et al., 2012). For example, using H3K27ac as a predictor of enhancer activity across multiple tissues, Nord et al., 2013 identified that 52% of predicted enhancers are active in more than one organ, of which 31% are active in two organs and 21% are active in three organs. Additionally, these studies have also revealed that pleiotropic enhancers may serve as binding sites for TFs that may be reused in different contexts. In Drosophila, Preger-Ben Noon et al., 2018 found that enhancer elements regulating the shavenbaby gene drive expression in multiple tissues and developmental stages. In one of these enhancers, the same TF binding site is used during both embryonic and pupal expression, while another enhancer utilizes different binding sites. In the context of brain development, a study that compared the binding of p300 at enhancers identified pleiotropic activity in both the embryonic mouse forebrain and midbrain at E11.5 (Hiller et al., 2012). Overall, a more detailed examination of enhancer activity across brain sub-regions and validation of the bound TFs and downstream target genes of enhancers with pleiotropic activity in the brain in figure studies will provide insight into the functional role these enhancers serve. Single-cell analysis of enhancer activity, as mentioned above, may provide further insight on how often enhancers function across cell types in the developing brain and the timing in which they are active.

GWAS SNPs and DNMs associated with ASD are enriched in cerebellar enhancers

Having established and characterized enhancer sequences in the cerebellum, we sought to elucidate the potential involvement of these regions in the etiology of neurological disorders and demonstrate how our dataset may be utilized to functionally annotate variation associated with neurodevelopmental disorders; imaging and quantitative data show consistent cerebellar abnormalities, particularly in cases of individuals with autism (Limperopoulos et al., 2014; Stoodley and Limperopoulos, 2016). Our results indicate that cerebellar enhancer sequences are significantly enriched for GWAS variants and DNMs associated with ASD, suggesting an important role for enhancers in contributing to the condition. PAX6 was among 12 target genes of cerebellar enhancers containing ASD-associated variants and is classified in the SFARI database as a ASD risk gene. The deletion of Pax6 in the murine cerebellum results in aberrant development of the glutamatergic cells in the cerebellum: the cerebellar nuclei, unipolar brush cells, and granule cells (Yeung et al., 2016). Behavioral analysis of Pax6 animal models has also indicated a possible link between this gene and autistic-like behavior (Umeda et al., 2010). Additionally, Pax6 has been linked with WAGR (Wilm’s tumor, Aniridia, Genitourinary malformations, and mental Retardation syndrome) which is co-morbid for ASD. Our analysis invites future investigation these target genes and how perturbation of their expression may lead to ASD phenotypes.

Of the target genes of the enhancers that overlapped de novo CNVs in the SSC cohort, none have been previously associated specifically with cerebellar development. Interestingly, one of these target genes, CDC42BPB, has recently been associated with neurodevelopmental disorders including ASD (Chilton et al., 2020). This gene is a serine/threonine protein kinase and codes for MRCKβ (myotonic dystrophy-related Cdc42-binding kinase beta), a regulator of cell cytoskeletal reorganization and cell migration (Pichaud et al., 2019). Of note, the CNV associated with this gene deletes the entire enhancer. CDC42BPB shows expression in the granule cell layer of the lateral adult cerebellum, which has been associated with cognitive functions (Koziol et al., 2014).

Together, our dataset provides an atlas of active enhancers during embryonic and postnatal cerebellar development. Not only have we demonstrated the utility of our dataset through the discovery of TFs with novel functions in the developing cerebellum, we also provide an invaluable resource for future studies. The online Developmental Cerebellum Enhancer Atlas can be used by the scientific community as a confirmatory tool for ongoing studies and to identify novel candidate genes involved in cerebellar development for future studies.

Materials and methods

Mouse strains and husbandry

Request a detailed protocol

C57BL/6  J mice were originally purchased from JAX laboratory and maintained and bred in our pathogen-free animal facility with 12/12  hour light/dark cycle and a controlled environment. Embryonic ages utilized in these experiments were confirmed based upon the appearance of a vaginal plug. The morning that a vaginal plug was detected was designated as E0.5. Pregnant females were cervically dislocated and embryos were harvested from the uterus. Postnatal ages were determined based upon the date of birth with the morning of the observation of newborn pups considered as P0.5. All studies were conducted according to the protocols approved by the Institutional Animal Care and Use Committee and the Canadian Council on Animal Care at the University of British Columbia.

Tissue preparation for chromatin immunoprecipitation

Request a detailed protocol

C57BL/6  J mice (male and female) at E12.5, P0.5 and P9.5 (referred to as E12, P0, and P9) were decapitated for dissection of cerebella. Cerebella were dissected and collected in ice cold Dulbecco’s PBS (DPBS) without magnesium or calcium and subsequently washed 2 x for 5 minutes. Samples from each litter were pooled and trypsinized in DPBS containing 0.25% trypsin for 10, 15, and 30 min at room temperature for E12, P0, and P9, respectively. Following three washes with fresh DPBS, the tissue was triturated with three progressively smaller (1, 0.5, 0.1 mm) bore polished and sterile pipettes in DPBS containing 250 U/ml DNase, 0.25% glucose, and 8 mg/ml BSA. The triturated cells were diluted 1:4 with cold DPBS and passed through a cell strainer (40 μm mesh) to remove large cellular debris. The cells were collected by mild centrifugation, washed in fresh DPBS and counted. The cells were split into 100,000 cell aliquots, pelleted and snap frozen using liquid nitrogen. Cell pellets were stored at –80 °C.

Histone chromatin immunoprecipitation

Request a detailed protocol

We performed native chromatin immunoprecipitation (ChIP) using validated antibodies against H3K4me1 and H3K27ac according to previously established protocols by the International Human Epigenomics Consortium (IHEC) (Lorzadeh et al., 2017). Briefly, cells were lysed in mild non-ionic detergents (0.1% Triton X-100 and Deoxycholate) and protease inhibitor cocktail (Calbiochem) to preserve the integrity of histones harbouring epitopes of interest during cell lysis. Cells were digested by Microccocal nuclease (MNase) at room temperature for 5 min and 0.25 mM EDTA was used to stop the reaction. Antibodies to H3K4me1 (Diagenode: Catalogue #C15410037, lot A1657D) and H3K27ac (Hiroshi Kimura, Cell Biology Unit, Tokyo Institute of Technology) were incubated with anti-IgA magnetic beads (Dynabeads from Invitrogen) for 2 hr. Digested chromatin was incubated with magnetic beads alone for 1.5 hr. Digested chromatin was separated from the beads and incubated with antibody-bead complex overnight in immunoprecipitation buffer (20 mM Tris-HCl pH 7.5, 2 mM EDTA, 150 mM NaCl, 0.1% Triton X-100, 0.1% Deoxycholate). The resulting immunoprecipitations were washed twice by low salt (20 mM Tris-HCl pH 8.0, 2 mM EDTA, 150 mM NaCl, 1% Triton X-100, 0.1% SDS) and then high salt (20 mM Tris-HCl pH 8.0, 2 mM EDTA, 500 mM NaCl, 1% Triton X-100, 0.1% SDS) wash buffers. Immunoprecipitations were eluted in an elution buffer (1% SDS, 100 mM Sodium Bicarbonate) for 1.5 hr at 65 °C. Remaining histones were digested by Protease (Invitrogen) for 30 min at 50 °C and DNA fragments were purified using Ampure XP beads (Beckman Coulter). The library preparation was conducted by Diagenode ChIP-seq/ChIP-qPCR Profiling service (Diagenode Cat# G02010000) using the MicroPlex Library Preparation Kit v2 (Diagenode Cat. C05010013). 50 bp single-end sequencing was performed on all libraries by Diagenode (Belgium) on an Illumina HiSeq 3000 platform. Two independent biological replicates were performed for each antibody and developmental time point.

Histone modification ChIP-seq data processing

Request a detailed protocol

The sequencing data were uploaded to the Galaxy web platform (usegalaxy.org) for analyses (Afgan et al., 2016). 50 bp single-end ChIP-seq reads were aligned to the NCBI37/mm9 reference genome and converted to binary alignment/map (BAM) format by Bowtie2 v.2.3.4 (Langmead et al., 2009) with default parameters. Duplicate reads were marked using Picard v.1.52. Quality control metrics after alignment, which were used to evaluate immunoprecipitation sensitivity, were calculated using ChIPQC v.4.2 (Carroll et al., 2014). Peak enrichment was computed using MACS v.2.1.1 (Zhang et al., 2008) with a false discovery rate (FDR) cutoff of 0.01 (P-value <1E-5) using input samples as a control for each replicate. bigWigs were generated and normalized by the total number of mapped reads using the BamCompare and profiles were generated from these bigWigs by calculating average coverage in 50 bp bins using Deeptools v.3.3 (Ramírez et al., 2016) for downstream analysis and visualization.

Identification of active cerebellar enhancers

Request a detailed protocol

We first determined consensus peaks between replicates for both H3K27ac and H3K4me1 peaks collected at E12, P0, and P9 using the intersect function from Bedtools v.2.28 (Quinlan and Hall, 2010). Robust active cerebellar enhancers were identified by overlapping replicated H3K27ac and H3K4me1 peaks called for E12, P0, and P9 samples. The genomic coordinates of the H3K27ac peaks that overlapped with H3K4me1 enriched regions at the same age were used for our list of robust active cerebellar enhancers. We then removed peaks found within promoter sequences by eliminating any peaks found 500 bp upstream or downstream of transcription start sites (TSSs) in the developing cerebellum as determined previously (Zhang et al., 2018). The resulting list of robust active cerebellar enhancer sequences at E12, P0, and P9 were used for downstream analysis.

For the comparative analysis with cerebellar postnatal enhancers previously published by Frank et al., 2015, H3K27ac and DNase-seq peak coordinates were downloaded from Gene Expression Omnibus (GEO) (GSE60731). The following sequences were downloaded from public enhancer databases: (1) enhancers downloaded from the VISTA Enhancer Browser (https://enhancer.lbl.gov/) (Visel et al., 2007) with hindbrain activity were filtered using the ‘Advanced Search’ tool, selecting ‘hindbrain’ under Expression Pattern and retrieving only mouse sequences with positive signal and (2) mouse cerebellar neonate enhancer coordinates were downloaded from the Enhancer Atlas 2.0 repository (http://www.enhanceratlas.org/downloadv2.php) (Gao and Qian, 2020). For the comparisons, sequences were overlapped with our robust cerebellar enhancer peaks and H3K27ac peaks at E12, P0 and/or P9 using Bedtools v.2.28 (Quinlan and Hall, 2010).

Differential binding analysis

Request a detailed protocol

Aligned read counts (BAM file format) from our H3K27ac ChIP-seq experiments mapped to our robust active cerebellar enhancers for E12, P0, and P9 samples were used as input to the package DiffBind v1.4.2 (Stark and Brown, 2011). Read counting at each genomic location was conducted, which was subsequently normalized by experimental input samples. The result of counting is a binding affinity matrix containing normalized read counts for every sample at each robust active cerebellar enhancer. For differential binding affinity analysis, three contrasts were set up in DiffBind: E12 vs P0, E12 vs P9, and P0 vs P9. Differential binding was determined by DiffBind using a negative binomial test at an FDR <0.05 threshold. The FDRs and normalized signal difference for each contrast were plotted using the EnhancedVolcano package in R (Blighe et al., 2020).

Temporal classification of cerebellar enhancers

Request a detailed protocol

To determine cerebellar enhancers with embryonic-specific activity, H3K27ac signal at E12 was compared to P0 and P9. Enhancers with significantly higher signal at E12 for either contrast were considered ‘Early’ active enhancers. A region found to be enriched for both contrasts was counted as one enhancer. To determine cerebellar enhancers with postnatal-specific activity, H3K27ac signal at P9 was compared to E12 and P0. Enhancers with significantly higher signal at P9 for either contrast were considered ‘Late’ active enhancers. A region found to be enriched for both contrasts was counted as one enhancer. To determine cerebellar enhancers with activity specific to birth, H3K27ac signal at P0 was compared to P9 and E12. Enhancers with significantly higher signal at P0 in both contrasts would identify enhancers that peaked in activity at P0. We did not identify any enhancers that peaked in activity at P0 and conducted the remaining analysis for only Early and Late enhancers.

Comparison with ENCODE developing mouse brain datasets

Request a detailed protocol

Cerebellar H3K4me1 and H3K27ac peaks and robust cerebellar enhancer coordinates were compared to ENCODE histone profiles quantified in the developing mouse hindbrain, midbrain and forebrain (Gorkin et al., 2020). ENCODE H3K4me1 and H3K27ac Bigwig and BED files for E12, P0 and Adult samples were downloaded from the ENCODE data Collection and Coordination website (https://www.encodedcc.org) using the experiment IDs listed in the supplementary information provided by the study. Cerebellar histone peaks were overlapped with hindbrain, midbrain and forebrain peaks using the intersect function from Bedtools v2.30.0 (Quinlan and Hall, 2010). H3K27ac ChIP-seq signal was then examined for Early and Late robust cerebellar enhancers in hindbrain, midbrain and forebrain samples. This was calculated using the ‘multiBigwigSummary’ function from DeepTools v.3.3 (Ramírez et al., 2016).

Comparative analysis with developing cerebellar cis-regulatory sequences (CREs) identified by snATAC-seq

Request a detailed protocol

Robust cerebellar enhancers were compared to CREs identified previously in the developing cerebellum by snATAC-seq (Sarropoulos et al., 2021). Datasets containing CRE coordinates, global CRE clusters and associated information about peaks and cells were downloaded from https://apps.kaessmannlab.org/mouse_cereb_atac/. The CRE coordinates and global CRE clusters were isolated from the file ‘Mouse_Cereb_ATAC_CRE_info.txt’. CRE locations were overlapped with robust cerebellar enhancer genomic coordinates using the ‘intersect’ function from Bedtools v2.30.0 (Quinlan and Hall, 2010). Robust cerebellar enhancers were assigned a cell-type based on the cell-type in which an overlapping CRE was detected. CREs were previously clustered into 26 clusters based on activity using an iterative clustering procedure (Sarropoulos et al., 2021). Average chromatin accessibility signal for each cell type at all available developmental stages were isolated from the file ‘Mouse_Cereb_ATAC_CREs_SE.rds’ using a custom script (Source Code File 1) provided by the corresponding authors. The CREs were filtered for those that overlapped with Early and Late active robust cerebellar enhancers and ATAC-seq signal was plotted for progenitors, granule cells, Purkinje cells, and interneurons.

Transcription factor motif enrichment analysis

Request a detailed protocol

Transcription factor motif enrichment was calculated using the software HOMER using the script FindMotifsGenome.pl with default parameters (Heinz et al., 2010). Analyses for Early and Late active enhancers were conducted separately. Motif enrichment was statistically analyzed using a cumulative binomial distribution. Enriched motifs were aligned with known transcription factor binding sites to determine the best matches. Top known motif matches were filtered based on expression within the developing cerebellum at E12 for ‘Early’ active enhancers and P9 for ‘Late’ active enhancers.

Cerebellar enhancer target gene prediction and co-expression analysis

Request a detailed protocol

To identify possible gene targets of our robust cerebellar enhancers, the correlation between H3K27ac signal and mRNA expression of genes located in cis at E12, P0, and P9 was calculated. For a given enhancer, a gene located in cis was considered a possible target if it was positively correlated with H3K27ac signal throughout time. These genes were then filtered based on location using conserved topologically associating domains (TADs), which are areas of the genome that preferentially interact (Dixon et al., 2012). These TADs are conserved between different cell types and even across species and were established using Hi-C data generated, previously. Gene target candidates for a given enhancer were curated for those located within the same TAD. Predicted gene targets were then ranked based on their Pearson correlation coefficient value. For the predicted gene targets of Early and Late active enhancers, we conducted k-means clustering of predicted gene targets separately. Input for this analysis was gene expression captured from cerebella at 12 embryonic and postnatal time points (Ha et al., 2019). Briefly, gene expression was quantified using Cap Analysis of Gene Expression followed by sequencing (CAGE-seq) for mouse cerebellar samples dissected every 24 hr from E11-P0 and then every 72 hr until P9 (12 in vivo time points in total). The number of clusters for the k-means clustering was determined using the Elbow analysis for each classified group of enhancers: Early (n=4) and Late (n=4).

A permutation analysis was conducted to evaluate whether cerebellar genes were enriched in robust cerebellar enhancer putative target genes. We generated 10,000 permutations of 2261 (total number of putative targets) randomly selected from all expressed genes in the developing cerebellum (~30,000 genes) (Ha et al., 2019). The number of cerebellar genes were then assessed in each permutation. One-sided p-values were calculated by dividing the number permutations containing a given number cerebellar genes by the total number of permutations (10,000).

Expression analysis of mouse and human single-cell RNAseq (scRNA-seq) datasets

Request a detailed protocol

Pax3 and Bhlhe22 expression values were analyzed in scRNA-seq datasets generated in the developing mouse and human developing cerebellum (Carter et al., 2018; Aldinger et al., 2021). The mouse cerebellar dataset was received as an ‘h5’ file from the corresponding authors, containing normalized expression values and annotations for all cells analyzed. Average normalized expression values were then calculated for 18 cell-type clusters previously determined using cell-type-specific markers (Carter et al., 2018). For the human cerebellar scRNA-seq dataset, we received a dataset from the corresponding author containing average normalized expression values for the 21 cell-type clusters previously determined using cell-type-specific markers (Aldinger et al., 2021). To evaluate Pax3, Bhlhe22, and Cdc42bpb expression in the cell types of the mouse and human developing cerebellum, normalized expression values for each cell-type cluster were plotted in bar and radar plots.

Tissue preparation for histology

Request a detailed protocol

Embryos harvested between E12.5 to E15.5 were fixed by immersion in 4% paraformaldehyde in 0.1 M phosphate buffer (PB, pH 7.4) for 1 hr at 4 °C. Postnatal mice between P0.5 to P9.5 were perfused through the heart with a saline solution followed by 4% paraformaldehyde/0.1 M PBS. The brain tissues were isolated and further fixed in 4% paraformaldehyde in 0.1 M PB for 1 hr at room temperature. Fixed tissues were rinsed with PBS, followed by cryoprotection with 30% sucrose/PBS overnight at 4 °C before embedding in the Optimal Cutting Temperature compound (Tissue-Tek). Tissues were sectioned at 12 µm for immunofluorescence experiments and cryosections were mounted on Superfrost slides (Thermo Fisher Scientific), air dried at room temperature, and stored at −80 °C until used. Sagittal sections were cut from one side of the cerebellum to the other (left to right, or vice versa). In all cases, observations were based on a minimum of 3 embryos per genotype per experiment.

Cerebellar immunostaining

Request a detailed protocol

Tissue sections were first rehydrated in PBS (3x5 minute washes) followed by a phosphate buffered saline with Triton X-100 (PBS-T) rinse. Sections were then incubated at room temperature for 1 hr with blocking solution (0.3% BSA, 10% normal goat serum, 0.02% sodium azide in PBS-T). Following the blocking step, the slides were incubated with primary antibody in incubation buffer (0.3% BSA, 5% normal goat serum, 0.02% sodium azide in PBS-T) at room temperature overnight in a humid chamber. Following the overnight incubation, the slides were rinsed in 3x10 min PBS-T washes. The sections were then incubated with the appropriate secondary antibody at room temperature for 1 hr, followed by three 0.1 M PB washes and one 0.01 M PB wash. Coverslips were applied to the slides using FluorSave mounting medium (345789, Calbiochem). The primary antibodies used were: rabbit anti-Bhlhe22 (1:1000, a gift from Dr. Michael Greenburg, Harvard University), mouse anti-Neurod1 (1:500, Abcam, ab60704), mouse anti-Pax3 (1:500, R&D systems, MAB2457), rabbit anti-Pax2 (1:200, Invitrogen, 71–6000), mouse anti-NeuN (1:100, Millipore, MAB377), rabbit anti-Calbindin (1:1000, Millipore, AB1778), rabbit anti-Foxp2 (1:2000, Novus Biologicals NB100-55411), chicken anti-Doublecortin: (1:100, Abcam ab153668). For immunofluorescence, secondary antibodies (Invitrogen) labeled with fluorochrome were used to recognize the primary antibodies.

Granule cell culture

Request a detailed protocol

Granule cells were isolated and cultured as previously described (Lee et al., 2009). Briefly, cerebella from litters of P6 mice were pooled and digested at 37 °C for 20 min in 10 U ml−1 papain (Worthington), and 250 U ml−1 DNase in EBSS using the Papain Dissociation Kit (Worthington, Cat #:LK003150). The tissue was mechanically triturated and suspended cells were isolated and resuspended in EBSS with albumin-ovomucoid inhibitor solution. Cell debris was removed using a discontinuous density gradient and cells were resuspended in HBSS, glucose and DNase. The cell suspension was then passed through a 40 µm cell strainer (Falcon 2340), layered on a step gradient of 35% and 65% Percoll (Sigma), and centrifuged at 2500 rpm for 12 min at 25 °C. Granule cells were harvested from the 35/65% interface and washed in HBSS-glucose. Granule cells were then resuspended in Neurobasal medium and 10% FBS and pre-plated on lightly coated poly-D-lysine-coated dishes for 20 min. This step allows any heavier cells to drop and adhere to the coated surface while the granule cells are retained in the media. Granule cells in the media were then collected, washed and counted using a Hemocytomoter. The cells were then plated on 25 mm or 12 mm poly-D-lysine (Sigma), laminin-coated coverglasses placed in six-well plates with Neurobasal medium containing B-27 serum-free supplement, 2 mm l-glutamine, 100 U/ml penicillin, and 100 μg/ml streptomycin (pen-strep; Invitrogen, Grand Island, NY) and 0.45% d-glucose (Gibco). Granule cells were incubated at 37 °C at 5% CO2 were incubated for 3 days in vitro (DIV).

For aggregate cultures, aggregates were allowed to form by incubating purified granule cells for 20 hr on uncoated tissue culture dishes in DMEM containing 10% FBS, 0.45% D-glucose, Pen-strep, 2 mM L-glutamine at 4E6 cells/ml. Aggregates were then washed and cultured in Neurobasal/B27 medium on poly-d-lysine/laminin-treated chamber slides at 37 °C/5% CO2. Neuronal processes extend from aggregates and most form neurite bundles. After several hours, small bipolar granule cells migrate unidirectionally away from the cell clusters along these neurites and neurite bundles by extension of processes, followed by translocation of cell bodies outside of the aggregate cell cluster margin. For immunofluorescence experiments, cells were fixed in 4% paraformaldehyde for 10 min and washed with calcium and magnesium-free PBS.

RNA interference

Request a detailed protocol

For the knockdown of Bhlhe22, we purchased ON-TARGETplus SMARTPool Mouse Bhlhe22 siRNA from Horizon Discovery (Cat ID: L-063262–01). Control samples were transfected with ON-TARGETplus Non-targeting Control Pool (Cat ID: D-001810–10). siRNA molecules were electroporated into isolated postnatal cerebellar granule cells using the Nucleofector 2b Device (Lonza, AAB-1001) as previously described (Gärtner et al., 2006). Briefly, after cells were isolated (described above), 6-7E6 cells were resuspended in nucleofection solution and mixed with 3 µg of pCAG-EGFP plasmid (Addgene, 89684) and 600 nM of siRNA. Cuvettes were loaded with cellular solution and nucleofected using program O-03. After electroporation, cells were allowed to recover in DMEM media in a humidified 37 °C/5% CO2 incubator for 90 min. Cells were washed and resuspended in either culture media for plating (dissociated cultures) or DMEM media for overnight incubation (aggregate cultures).

RNA isolation and reverse transcription followed by quantitative PCR (RT-qPCR)

Request a detailed protocol

RNA was collected from cultured granule cells using the Monarch Total RNA Miniprep Kit (NEB). Then cDNA was reverse transcribed using SuperScript IV First-Strand Synthesis System (Invitrogen) using random hexamers. Quantitative PCR was conducted using the Applied Biosystems Fast SYBR Green Master Mix reagent and Applied Biosystems 7500 Real-time PCR system. PCR conditions were as follows: 95  °C for 20 seconds, 40 cycles of 95  °C for 3  s, and 60  °C for 30 s followed by 95  °C for 15  s, 60  °C for 1  min, 95  °C for 15  s and 60  °C for 15  s. Three biological replicates were analyzed for each target gene. Amplification of eGFP was used as a reference gene to normalize the relative amounts of successfully transfected cells between treated and control experiments. Gene specific primers are listed in Supplementary file 9. Expression profiles for each gene were calculated using the average relative quantity of the sample using the deltaCT method (Livak and Schmittgen, 2001). For comparisons between siRNA treated and control samples, means were compared using a two-tailed t-test. Results were expressed as the average ± SE, and p-values <0.05 were considered significant.

Image analysis and microscopy

Request a detailed protocol

Analysis and photomicroscopy were performed using a Zeiss Axiovert 200 M microscope with the Axiocam/Axiovision hardware-software components (Carl Zeiss) and downstream image analysis was conducted using the AxioVision software v.4.9.1 (Carl Zeiss).

For cerebellar granule cell aggregate cultures, aggregate size determined using the tracing tool and all aggregates analyzed were within 1000 µm2 of each other. Transfected cells were identified by examining eGFP expression and for each biological replicate/experimental treatment, 20 aggregates were examined. Granule cell migration was measured by calculating the distance of migrated cells from the edge of the aggregate on captured images. Mean migration distance was calculated for each aggregate, and the average of all 20 aggregates was used for statistical analysis. The distribution of migrated cells from the aggregate was calculated for the following ranges:<25 µm, 25–50 µm, 50–75 µm, 75–100 µm, >100 µm. For each range, the average percentage was calculated for 20 aggregates per replicate. For comparisons between siRNA treated and control samples, means were compared using a two-tailed t-test. Results were expressed as the average ±SE, and p-values <0.05 were considered significant.

For the quantification of Pax3 and Bhlhe22 expressing cells, images captured for immunofluorescent staining were used for cell counting. Three randomly selected 100 µm square portions of each sagittal section of the developing cerebellum were analyzed at the indicated developmental stages. The percentage Pax3 and Bhlhe22 expressing cells that colocalized with cell type markers was determined by counting the number of co-stained cells dividing the total by the number of Pax3 or Bhlhe22 expressing cells in a square. Three square portions were counted for three sagittal sections for each of the three biological replicates. Each data point on the respective bar graphs represent an average calculated for each sagittal section.

Plots and statistical methods

Request a detailed protocol

All plots and correlation analysis were generated in R version 3.2.3 and figures were produced using the package ggplot2. Bedtools v.2.28 (Quinlan and Hall, 2010) was used for comparing and overlapping the genomic coordinates of peaks and existing genomic features described in the manuscript. Gene Ontology enrichment analyses were conducted using the clusterProfiler package in R (Yu et al., 2012). Boxplots represent the median (centre line), first and third quartiles (top and bottom of box, respectively) and confidence intervals (95%; black lines). Genome browser screenshots were taken from the IGV genome browser (Robinson et al., 2011). Bar plots results were expressed as the average and the corresponding error bars represent standard error.

GWAS SNP enrichment analysis

Request a detailed protocol

The NCBI LiftOver tool was used to convert the coordinates of cerebellar enhancers from mm9 to hg19 to hg38, and BEDTools (Quinlan and Hall, 2010). 6630/7024 robust cerebellar enhancers were converted from mm9 to hg38. Single-nucleotide polymorphisms (SNPs) were retrieved from the GWAS Catalog (Buniello et al., 2019) downloaded on March 8th, 2020. The SNPs were then filtered by their associated traits. Traits containing the word ‘autism’ were selected and from this list any traits containing the word ‘or’ were excluded. This resulted in a final list of 8 traits (Supplementary file 10) and the associated SNPs were used as input for our analysis. The software Genomic Regulatory Elements and Gwas Overlap algoRithm (GREGOR) v.1.4.0 (Schmidt et al., 2015), a tool to test for enrichment of an input list of trait-associated index SNPs in experimentally annotated regulatory domains, was used to identify enrichment of trait-specific disease variants within enhancers. An underlying hypothesis of GREGOR is that both trait-associated SNPs and variants in strong linkage disequilibrium (LD) may be deemed as causal. For this, we used the European population reference file (EUR; LD window size = 1 Mb; LD r2 ≥0.7) from 1000 G data (Release date: May 21, 2011). The probability of an overlap of either a SNP or at least one of its LD proxies with our enhancers relative to a set of matched control variants was used to evaluate significance of overlap. The enrichment p-value is the probability that the overlap of control variants with our enhancers is greater than or equal to the overlap of the GWAS variants with converted robust cerebellar enhancers.

De novo mutation analysis

Request a detailed protocol

De novo mutations were detected using whole-genome sequencing data from the MSSNG (Yuen et al., 2016) and Simons Simplex Collection (SSC) (Isoda et al., 2017) cohorts using a pipeline involving DeNovoGear (Ramu et al., 2013) as previously described (C Yuen et al., 2017). To maximize data homogeneity, we included only individuals sequenced on the Illumina HiSeq X platform. Individuals having a total DNM count more than three standard deviations above the mean of the cohort were excluded. The NCBI LiftOver tool was used to convert the coordinates of cerebellar enhancers from mm9 to hg19 to hg38, and BEDTools (Quinlan and Hall, 2010) was used to identify DNMs overlapping these coordinates. Contingency tables (2x2) were generated containing counts of the number of DNMs in ASD-affected individuals and unaffected siblings either overlapping or not overlapping each dataset (cerebellar enhancer or H3K27ac peak coordinates). Fisher’s exact test was used to determine statistical significance. Copy number variants (CNVs) ≥ 1000 bp were detected from the MSSNG and SSC WGS data using a pipeline involving the algorithms ERDS (Zhu et al., 2012) and CNVnator (Abyzov et al., 2011) as previously described (Trost et al., 2018). A CNV was deemed to be de novo if it was detected by both ERDS and CNVnator in the child but by neither algorithm in both parents. We then used BEDtools (Quinlan and Hall, 2010) to identify de novo CNVs overlapping our cerebellar enhancers.

Data availability

Sequencing data have been deposited in GEO under accession code: GSE183697.

The following data sets were generated
    1. Ramirez M
    2. Badayeva Y
    3. Yeung J
    4. Wu J
    5. Yang E
    (2021) NCBI Gene Expression Omnibus
    ID GSE183697. Temporal analysis of enhancers during mouse brain development reveals dynamic regulatory function and identifies novel regulators of cerebellar development.
The following previously published data sets were used
    1. Frank CL
    2. Liu F
    3. Wijayatunge R
    4. Song L
    5. Biegler MT
    6. Yang MG
    7. Vockley CM
    8. Safi A
    9. Gersbach CA
    10. Crawford GE
    11. West AE
    (2015) NCBI Gene Expression Omnibus
    ID GSE60731. Regulation of chromatin accessibility and Zic binding at enhancers in the developing cerebellum.
    1. Aldinger KA
    (2021) Human Cell Atlas
    ID aldinger20. Spatial and cell type transcriptional landscape of human cerebellar development.
    1. Carter RA
    (2018) European Nucleotide Archive
    ID PRJEB23051. A Single-Cell Transcriptional Atlas of the Developing Murine Cerebellum.
    1. Sarropoulos I
    (2021) Developing mouse cerebellum snATAC-seq atlas
    ID kaessmannlab. Developmental and evolutionary dynamics of cis-regulatory elements in mouse cerebellar cells.

References

  1. Software
    1. Blighe K
    2. Rana S
    3. Lewis M
    (2020)
    EnhancedVolcano: publication-ready volcano plots with enhanced colouring and labeling
    EnhancedVolcano.
    1. Forrest ARR
    2. Kawaji H
    3. Rehli M
    4. Baillie JK
    5. de Hoon MJL
    6. Haberle V
    7. Lassmann T
    8. Kulakovskiy IV
    9. Lizio M
    10. Itoh M
    11. Andersson R
    12. Mungall CJ
    13. Meehan TF
    14. Schmeier S
    15. Bertin N
    16. Jørgensen M
    17. Dimont E
    18. Arner E
    19. Schmidl C
    20. Schaefer U
    21. Medvedeva YA
    22. Plessy C
    23. Vitezic M
    24. Severin J
    25. Semple CA
    26. Ishizu Y
    27. Young RS
    28. Francescatto M
    29. Alam I
    30. Albanese D
    31. Altschuler GM
    32. Arakawa T
    33. Archer JAC
    34. Arner P
    35. Babina M
    36. Rennie S
    37. Balwierz PJ
    38. Beckhouse AG
    39. Pradhan-Bhatt S
    40. Blake JA
    41. Blumenthal A
    42. Bodega B
    43. Bonetti A
    44. Briggs J
    45. Brombacher F
    46. Burroughs AM
    47. Califano A
    48. Cannistraci CV
    49. Carbajo D
    50. Chen Y
    51. Chierici M
    52. Ciani Y
    53. Clevers HC
    54. Dalla E
    55. Davis CA
    56. Detmar M
    57. Diehl AD
    58. Dohi T
    59. Drabløs F
    60. Edge ASB
    61. Edinger M
    62. Ekwall K
    63. Endoh M
    64. Enomoto H
    65. Fagiolini M
    66. Fairbairn L
    67. Fang H
    68. Farach-Carson MC
    69. Faulkner GJ
    70. Favorov AV
    71. Fisher ME
    72. Frith MC
    73. Fujita R
    74. Fukuda S
    75. Furlanello C
    76. Furino M
    77. Furusawa J
    78. Geijtenbeek TB
    79. Gibson AP
    80. Gingeras T
    81. Goldowitz D
    82. Gough J
    83. Guhl S
    84. Guler R
    85. Gustincich S
    86. Ha TJ
    87. Hamaguchi M
    88. Hara M
    89. Harbers M
    90. Harshbarger J
    91. Hasegawa A
    92. Hasegawa Y
    93. Hashimoto T
    94. Herlyn M
    95. Hitchens KJ
    96. Ho Sui SJ
    97. Hofmann OM
    98. Hoof I
    99. Hori F
    100. Huminiecki L
    101. Iida K
    102. Ikawa T
    103. Jankovic BR
    104. Jia H
    105. Joshi A
    106. Jurman G
    107. Kaczkowski B
    108. Kai C
    109. Kaida K
    110. Kaiho A
    111. Kajiyama K
    112. Kanamori-Katayama M
    113. Kasianov AS
    114. Kasukawa T
    115. Katayama S
    116. Kato S
    117. Kawaguchi S
    118. Kawamoto H
    119. Kawamura YI
    120. Kawashima T
    121. Kempfle JS
    122. Kenna TJ
    123. Kere J
    124. Khachigian LM
    125. Kitamura T
    126. Klinken SP
    127. Knox AJ
    128. Kojima M
    129. Kojima S
    130. Kondo N
    131. Koseki H
    132. Koyasu S
    133. Krampitz S
    134. Kubosaki A
    135. Kwon AT
    136. Laros JFJ
    137. Lee W
    138. Lennartsson A
    139. Li K
    140. Lilje B
    141. Lipovich L
    142. Mackay-Sim A
    143. Manabe R
    144. Mar JC
    145. Marchand B
    146. Mathelier A
    147. Mejhert N
    148. Meynert A
    149. Mizuno Y
    150. de Lima Morais DA
    151. Morikawa H
    152. Morimoto M
    153. Moro K
    154. Motakis E
    155. Motohashi H
    156. Mummery CL
    157. Murata M
    158. Nagao-Sato S
    159. Nakachi Y
    160. Nakahara F
    161. Nakamura T
    162. Nakamura Y
    163. Nakazato K
    164. van Nimwegen E
    165. Ninomiya N
    166. Nishiyori H
    167. Noma S
    168. Noma S
    169. Noazaki T
    170. Ogishima S
    171. Ohkura N
    172. Ohimiya H
    173. Ohno H
    174. Ohshima M
    175. Okada-Hatakeyama M
    176. Okazaki Y
    177. Orlando V
    178. Ovchinnikov DA
    179. Pain A
    180. Passier R
    181. Patrikakis M
    182. Persson H
    183. Piazza S
    184. Prendergast JGD
    185. Rackham OJL
    186. Ramilowski JA
    187. Rashid M
    188. Ravasi T
    189. Rizzu P
    190. Roncador M
    191. Roy S
    192. Rye MB
    193. Saijyo E
    194. Sajantila A
    195. Saka A
    196. Sakaguchi S
    197. Sakai M
    198. Sato H
    199. Savvi S
    200. Saxena A
    201. Schneider C
    202. Schultes EA
    203. Schulze-Tanzil GG
    204. Schwegmann A
    205. Sengstag T
    206. Sheng G
    207. Shimoji H
    208. Shimoni Y
    209. Shin JW
    210. Simon C
    211. Sugiyama D
    212. Sugiyama T
    213. Suzuki M
    214. Suzuki N
    215. Swoboda RK
    216. ’t Hoen PAC
    217. Tagami M
    218. Takahashi N
    219. Takai J
    220. Tanaka H
    221. Tatsukawa H
    222. Tatum Z
    223. Thompson M
    224. Toyodo H
    225. Toyoda T
    226. Valen E
    227. van de Wetering M
    228. van den Berg LM
    229. Verado R
    230. Vijayan D
    231. Vorontsov IE
    232. Wasserman WW
    233. Watanabe S
    234. Wells CA
    235. Winteringham LN
    236. Wolvetang E
    237. Wood EJ
    238. Yamaguchi Y
    239. Yamamoto M
    240. Yoneda M
    241. Yonekura Y
    242. Yoshida S
    243. Zabierowski SE
    244. Zhang PG
    245. Zhao X
    246. Zucchelli S
    247. Summers KM
    248. Suzuki H
    249. Daub CO
    250. Kawai J
    251. Heutink P
    252. Hide W
    253. Freeman TC
    254. Lenhard B
    255. Bajic VB
    256. Taylor MS
    257. Makeev VJ
    258. Sandelin A
    259. Hume DA
    260. Carninci P
    261. Hayashizaki Y
    262. FANTOM Consortium and the RIKEN PMI and CLST (DGT)
    (2014) A promoter-level mammalian expression atlas
    Nature 507:462–470.
    https://doi.org/10.1038/nature13182
    1. Frantz GD
    2. Bohner AP
    3. Akers RM
    4. McConnell SK
    (1994)
    Regulation of the POU domain gene SCIP during cerebral cortical development
    The Journal of Neuroscience 14:472–485.
  2. Software
    1. Stark R
    2. Brown G
    (2011)
    DiffBind: differential binding analysis of chip-seq peak data
    DiffBind.

Decision letter

  1. Genevieve Konopka
    Reviewing Editor; University of Texas Southwestern Medical Center, United States
  2. Kathryn Song Eng Cheah
    Senior Editor; University of Hong Kong, Hong Kong
  3. Genevieve Konopka
    Reviewer; University of Texas Southwestern Medical Center, United States
  4. Alex S Nord
    Reviewer; University of California, Davis, United States

In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.

Decision letter after peer review:

Thank you for resubmitting your work entitled "Temporal analysis of enhancers during mouse brain development reveals dynamic regulatory function and identifies novel regulators of cerebellar development" for further consideration by eLife. Your revised article has been evaluated by Kathryn Cheah (Senior Editor) and a Reviewing Editor.

The manuscript has been improved but there are some remaining issues that need to be addressed, as outlined below:

Importantly, we think this manuscript would fit a "Resource" category better than its current status as a research article. However, even to meet that bar we believe there are substantial revisions that need to be made. The reviewers felt that while the data might constitute a valuable resource for specific researchers that the novelty and impact of the dataset overall was not sufficient for a research article. However, if the authors can integrate previously published datasets more effectively and make your dataset more easily digestible/accessible by a broader neuroscience audience, we believe this could be a valuable dataset for the community.

Reviewer #1:

This manuscript aims to identify and understand regulatory elements in the early developing mouse cerebellum. To this end, the authors generate datasets of H3K4me1 and H3K27ac ChIP-seq in embryonic and early postnatal mouse cerebellum. The authors use these data to identify potential enhancers and carry out a number of confirmatory experiments that overlap their data with previously published ones as well as showing expression patterns of target genes that fit with the enhancer data.

The strengths of the manuscript include:

1) the generation of embryonic data that have been missing in the literature.

2) the integration with previously published datasets.

3) the IHC/ISH experiments that correlate expression with enhancer activity.

4) the functional characterization of Bhlhe22, which had not previously been described in the literature.

The weaknesses of the manuscript are:

1) The bulk approach that could potentially mask patterns driven by rare cell types.

2) Lack of mutational testing of specific enhancers.

In general, I think the manuscript is well written and easy to understand. The authors make a case for the need and novelty of their datasets. They also try to bolster their findings using a number of additional wet-bench approaches as well as integration with published datasets.

The biggest weakness is the lack of cell-type resolution presented using this bulk approach. The authors do a reasonable job trying to address this with the data in hand and mention in the discussion that finer resolution is needed. I'm just not sure that given the current state of what is possible technically that the current approach provides sufficiently important data (although contribution of embryonic data is always needed).

Reviewer #2:

The authors here sought to identify enhancers that are active in developing cerebellum and distinguish early (e.g. embryonic) from late (e.g. early postnatal) activity. To accomplish this, the authors performed new ChIP-seq experiments on two histone modifications that differentiate poised/latent from active enhancers, H3K4me1 and H3K27ac respectively, interrogating E12, P0, and P9 mouse cerebellum. The authors performed differential binding analysis and intersected these epigenomic datasets with published transcriptomic and autism human genetics datasets. While the scope is limited, the datasets and analysis are of sufficient rigor and reproducibility and the authors validate two genes with novel cerebellar function, thus the results and resulting resource of annotated enhancers should be sound.

The main weakness with this study is limited novelty and impact due to studies that have generated similar datasets and general findings. Other studies have generated epigenomics datasets for developing mouse cerebellum, including recently with single cell resolution using snATAC-seq, and more detailed time courses and complete chromatin signatures. The datasets here are largely redundant with work from the mouse ENCODE project, which profiled hindbrain across eight developmental stages that overlap with the two early timepoints here, as well as adult cerebellum. Similarly, the general concept of enhancers anchoring dynamic developmental processes, including in brain and specifically cerebellum, is now firmly established. Nonetheless, focused studies such as this one have value in that the results and datasets are digestible, curated, and directed toward a specific audience that are likely to appreciate the findings and use the results as a resource. The value of the work is the focus on identification and characterization of enhancers and their target genes with differential early developmental versus post-natal activity in the cerebellum.

While this focused study should have an audience, it would be stronger if the authors interfaced their data and results with other larger and more granular (e.g. single cell or more chromatin marks) datasets. If this is done effectively, this manuscript will have clear stand alone value as a complement to the larger, but not curated or cerebellum focused ENCODE data. And would link to the emergine single cell datasets that also focus on dynamic transcriptional and epigenomic changes. This would put the results from this work in the context of the larger and centralized efforts, while maintaining the useful focus.

Analysis and results:

1. The IHC in the manuscript lacks any quantitation and co-labeling/co-expression patterns are not clear with the lower magnification images. Figure 3 and 6 (and supp fig) analysis should have a quantitative summary of staining and of co-labeling. For example, Figure 3 is used to argue for presence or absence of co-expression/co-labeling between Pax3 and Ptf1a, Pax2, and Foxp2. I would like to see neuroanatomical annotation lines as overlay on sections for non-experts. Finally, I think higher magnification insets that clearly show co-labeling (or lack thereof) would make these figures stronger.

2. The P0 enhancer set is much smaller than the E12 or P9, raising questions about technical sensitivity. The authors make the conclusion that there are no P0 enriched enhancers, but this could also simply be a lack of sensitivity. The authors could address this by presenting information regarding the signal to noise and comparing their data to published ENCODE P0 hindbrain data.

3. The authors should compare their datasets to the mouse ENCODE developmental hindbrain datasets. Ideally, they would verify the patterns they found, showing that their datasets are high quality and correlated with the ENCODE data and then also review the reproducibility of the early vs. late enhancers. The authors can also use the P9 data to compare to the adult cerebellum in mouse ENCODE to see how similar these datasets are – it would be of interest if P9 is capturing novel enhancers. The authors could also look at dynamic chromatin accessibility to see if the H3K27ac + H3K4me1 signature here is more sensitive to active enhancers than accessibility alone (which would be expected and would further justify the value of the data here).

4. It would be of interest to see how the enhancer sets identified here map to the recently published snATAC-seq datasets from mouse cerebellum across development. This publication was referenced in the manuscript, but not described – Sarropoulos Science 2021. More specifically, can the early and late enhancers here be assigned to shared or distinct cerebellar cell types? Also as above, it would add value to this study if the H3K27ac+H4K3me1 activity signatures produces improved sampling of active enhancers vs. the snATAC.

5. The authors could intersect their target early and late genes with recent snRNA-seq data from developing human cerebellum (Aldinger Nat Neuro 2021).

6. Some general issues with the target gene set analsyis:. For the k-means clustering of target genes presented in Figures 4 and 5, the authors should provide rationale for using the k value selected. K-means is a great algorithm for this due to simplicity and wide use, but can require optimization and nested approaches. What other values of k were used? Did a nested approach reveal any subclusters within the four major groups at either timepoint? It is fine to use some arbitrary/post-hoc selection of best clusters, but the rationale should be described. Figure 5 should include N for genes per cluster. Finally I'm not sure what the GeneRatio measure is supposed to capture. IT would be more straightforward to use an observed vs. expected fold enrichment, as counts can skew depending on number of genes and fold enrichment is standard in GO analysis.

7. Motif analysis: Motifs are often similar for families, so motif to TF matching is imprecise at best. Fine to select best guess TF as authors do, but good to also show motif logo in Figure 3A and include mention that the actual TF could be a different family member unless orthogonal evidence exists.

8. Intersection with ASD variants is flimsy. The authors should use some presumed negative for non-brain and other brain GWAS to show that the ASD SNP findings are relevant e.g. have stronger signal than presumed negative or less relevant phenotypes. The overall enrichment with de novo mutations in enhancesr does not demonistrate any specific cerebellar relevance unless these enhancers are not active in other brain regions. Other studies have showed general enrichment of de novo mutations in brain enhancers. Not sure whether the rare CNV with the cerebellar gene is useful or convincing of relevance as the N is 1 and the gene expression outside the cerebellum is bot shown or discussed. All this analysis is fine, but shouldn't be treated as a major finding. I think the other parts of this study are much stronger, so de-emphasizing this would not hurt the paper.

9. The argument that MGI mouse cerebellum phenotype is evidence for validity needs to be compared to some random set in order to evaluate how strong the 98/2261 finding is- e.g. via permutation analysis selecting random sets of 2261 genes to establish background.

10. Figure 2B is not needed – just put numbers in legend or overlay on volcano plots in 2A.

11. Figure 2C would be more informative with longitudinal data showing change over time.

Text

1. The authors highlight their study as the first to characterize developmentally active enhancers, both in the abstract and discussion. This is an over-statement at best and the authors would do better to both remove this language and to directly deal with some of the other datasets, as described below.

2. The introduction should include brief description of the other relevant studies, particularly the mouse ENCODE data and recent snATAC-seq and explain the novelty of approach here. The discussion should devote space to comparing the results from this study and the mouse ENCOE and snATAC-seq specifically.

3. The last paragraph of the introduction is essentially a detailed summary of the findings. I would prefer this be the first paragraph of the discussion to avoid redundancy and present results before making these claims. A brief general summary statement is fine.

Reviewer #3:

Strengths of the manuscript:

The authors have carried out a stringent analysis of genomic enhancers that are dynamically active in specific cell populations over a period of cerebellar development. Their data represents a highly useful resource for the field of cerebellar research.

An enrichment analysis of the enhancers for specific transcription factor binding sites identified potential transcription factors previously not well-known to regulate cerebellar development.

The authors demonstrate how their dataset could be exploited to identify novel genes in cerebellar development (e.g., Bhlhe22) as well as potential disease genes and mechanisms (e.g. in autism spectrum disorder).

Weaknesses:

The authors select three timepoints during cerebellar development for their ChIP-seq analysis. The rational for selecting these specific timepoints is not very clear, especially as they do not cover the entire period of cerebellar development in the mouse.

While the resource itself is very valuable, the presentation of example genes that the authors extract from their analyses (Pax3, Bhlhe22, and Cdc42bpb) is rather short and lacking in detail. For example, the authors could have further strengthened their findings by utilising recent single-cell datasets, in particular the single-cell transcriptional atlas of the developing cerebellum published by Carter et al.

It won't be easy for the scientific community to interact with the generated data.

The authors use in silico tools to identify transcription factor binding sites enriched in the identified active enhancers. No experimental validation is provided for this analysis. This would have been particularly relevant for transcription factors that have not been previously well-described in cerebellar development including Pax3, which the authors present as a novel marker for GABAergic progenitors. Furthermore, it would have been helpful to go into further detail and to provide some predicted candidate genes that might be regulated by Pax3 in GABAergic progenitors.

Figure 3: Very limited expression data is shown for Pax3. It would have been helpful to show additional immunostainings at least for P9 (same timepoint as used for the ChIP-seq analysis) if not for further timepoints during cerebellar development. Moreover, it would have been really interesting to compare this to the data available on Pax3 from Carter et al. The authors mention this in their discussion; it would be helpful to include this data in the relevant figure and also to comment on how expression of Pax3 changes during the course of cerebellar development.

With regards to Pax3, there are some interesting papers that link Pax3 to cerebellar development and cerebellar disorders and that the authors could comment on in their discussion; for example, Kioussi and Gruss (1994) and Jalali et al., (2008), and Zagozewski et al., (2020).

I am not sure that I agree with the authors about the fact that 98 out of 2261 identified target genes result in a cerebellar phenotype when knocked out particularly demonstrates the validity of the high throughput approach (page 12). Is this a significant number?

The use of CAGE-sequencing data in this manuscript is very unclear. The authors refer to a previous manuscript (Zhang et al., 2018); however, in this study reports transcriptomic data from only three timepoints (E13, E15, E18) and not 12 timepoints as mentioned in by the authors? The authors should better describe which data they used. Also, given that the CAGE-Seq data is from bulk-sequencing, it would have been much more informative to look at gene expression of their candidate transcription factors in the Carter et al., single-cell dataset. Furthermore, this would have allowed the authors to look beyond the P9 timepoint and thus for a more complete gene expression analysis over cerebellar development. This is also true for Figures 1E,G.

The approach of using the enhancer data set to identify novel target genes not previously identified in cerebellar development is unclear to me. Using their described filtering approach, the authors describe genes that differentially regulated by Atoh1. Could these not have been identified from the original dataset by Klisch et al., 2001? What does the initial identification of the enhancers add here?

Figure 6: As above, the authors should include analysis of Bhlhe22 expression in the Carter et al., single-cell dataset.

With regards to Bhlhe22, the authors should comment on the protein that is encoded by this gene and its potential function during cerebellar development.

Figure 7 could be much more clearly presented. The font size is really small and the figures not very clear. It would be nice to see successful knockdown of Bhlhe22 by immunostaining.

The authors normalize gene expression to co-transfected GFP. Can the authors demonstrate that cells that are transfected with GFP were also transfected with the siRNA construct?

Can the authors comment on the evolutionary conservation of their identified enhancer elements? This is particularly relevant for their analysis of variants associated with ASD. Are the loci of interest conserved between mouse and human and thus, do enhancer activity and transcription factor occupancy occur at orthologous regions distal to the same gene promoters?

For the identification of CDC42BPB, it would have been nice to link this back to expression data over time and to provide more detail on the potential transcription factors involved.

Finally, given that this is a very important dataset for the scientific community, I wish it would be easier for others to interact with the data rather than just depositing data files into GEO. It would have also been really helpful to provide supplementary tables listing enhancer elements and associated transcription factors as well as nearby genes that the authors have identified. Currently, it is not transparent which lists of genes etc were used for the various analyses shown.

https://doi.org/10.7554/eLife.74207.sa1

Author response

Reviewer #1:

The biggest weakness is the lack of cell-type resolution presented using this bulk approach. The authors do a reasonable job trying to address this with the data in hand and mention in the discussion that finer resolution is needed. I'm just not sure that given the current state of what is possible technically that the current approach provides sufficiently important data (although contribution of embryonic data is always needed).

The bulk tissue dataset that we present provides consensus (i.e. data from many cells) and high confidence signals that point us to coordinates that will serve as benchmarks for validation for past and future single-cell assays. We address the issue of cell-type resolution by comparing our data with recent single cell/nuclei data (Sarrapolous et al., 2021; Carter et al., 2018; Aldinger et al., 2021). The outcome is that we find a significant overlap between our data and the single cell/nuclei data. These analyses and results can be found in a new section in our manuscript entitled ‘Comparison of robust cerebellar enhancers with chromatin accessibility identified in single cells’ (pp. 11) and in the sections that report our findings Pax3 (pp.14-15) and Bhlhe22 (pp.19).

Reviewer #2:

The main weakness with this study is limited novelty and impact due to studies that have generated similar datasets and general findings. Other studies have generated epigenomics datasets for developing mouse cerebellum, including recently with single cell resolution using snATAC-seq, and more detailed time courses and complete chromatin signatures. The datasets here are largely redundant with work from the mouse ENCODE project, which profiled hindbrain across eight developmental stages that overlap with the two early timepoints here, as well as adult cerebellum. Similarly, the general concept of enhancers anchoring dynamic developmental processes, including in brain and specifically cerebellum, is now firmly established. Nonetheless, focused studies such as this one have value in that the results and datasets are digestible, curated, and directed toward a specific audience that are likely to appreciate the findings and use the results as a resource. The value of the work is the focus on identification and characterization of enhancers and their target genes with differential early developmental versus post-natal activity in the cerebellum.

We believe that our dataset is not redundant with the chromatin profiles generated from the entire rhombencephalon as reported in the ENCODE work of Gorkin et al., (2020). The rhombencephalon develops into the pontine nuclei, the medulla and other brainstem regions in addition to the cerebellum. Three points are notable here: (1) Our analysis, is specific to chromatin profiling of the developing cerebellum which is a novel dataset. (2) In our analysis, we examined the P9 timepoint which was not examined in the ENCODE dataset, and this is a key watershed in cerebellar development; at this timepoint we find a large subset of enhancers that are specific to postnatal cerebellar development. (3) When comparing our cerebellar histone profiles to those generated in the hindbrain, we find that ~25% of our H3K27ac peaks at E12 are not found in the ENCODE hindbrain data (Figure 2 —figure supplement 1A). Additionally, we identified 765 robust cerebellar enhancers active in the cerebellum and not in the hindbrain. Thus, our study has led to the discovery of novel enhancer elements, potentially important specifically for cerebellar development. These findings are presented in the Results section ‘A subset of cerebellar enhancers is active in other developing regions of the brain’ (pp. 10) in the revised manuscript. We also utilize our cerebellar histone profiles to uncover novel regulators of cerebellar development and characterize the expression and function for transcription factors Pax3 and Bhlhe22, which have not previously been described in the literature. In comparison to the ENCODE analysis previously conducted in the mouse hindbrain, our study provides results and datasets that are digestible and curated, which can be utilized to discover temporally-specific regulators of embryonic and postnatal cerebellar development.

While this focused study should have an audience, it would be stronger if the authors interfaced their data and results with other larger and more granular (e.g. single cell or more chromatin marks) datasets. If this is done effectively, this manuscript will have clear stand alone value as a complement to the larger, but not curated or cerebellum focused ENCODE data. And would link to the emergine single cell datasets that also focus on dynamic transcriptional and epigenomic changes. This would put the results from this work in the context of the larger and centralized efforts, while maintaining the useful focus.

We agree with this reviewer that the integration of the previously generated data will provide our study with stand-alone value as a complement to these larger datasets. We have addressed this reviewer’s suggest to interface our data with larger and more granular datasets (see our Cover Letter). Overall, by the addition of relevant comparisons, our study provides a valuable resource of active enhancers focused on the developing cerebellum and can be used as a means to discover novel genetic regulators of cerebellar development.

Analysis and results:

1. The IHC in the manuscript lacks any quantitation and co-labeling/co-expression patterns are not clear with the lower magnification images. Figure 3 and 6 (and supp fig) analysis should have a quantitative summary of staining and of co-labeling. For example, Figure 3 is used to argue for presence or absence of co-expression/co-labeling between Pax3 and Ptf1a, Pax2, and Foxp2. I would like to see neuroanatomical annotation lines as overlay on sections for non-experts. Finally, I think higher magnification insets that clearly show co-labeling (or lack thereof) would make these figures stronger.

We agree with this reviewer that a quantitative assessment of co-labeling will add further value to our findings and that higher magnification images of co-staining experiments will bring more clarity to our results of our co-staining experiments. To these ends, in the revised manuscript, we have now included a quantitative analysis of co-staining relative to the protein of interest. For Pax3 we quantified the percentage of Pax3 cells that colocalized with Ptf1a, Pax2, Foxp2, and Calbindin. These changes for Pax3 are detailed in Figure 3 and Figure 3 —figure supplement 1. In the revised manuscript, this analysis can be found in the section ‘Early active enhancers are enriched for Pax3 binding sites, a novel marker for GABAergic cells.’ (pp. 14).

For Bhlhe22 and markers for differentiating granule cells and migrating cells, we quantified the percentage of Bhlhe22 cells that colocalized with Neurod1, NeuN and Dcx. These data for Bhlhe22 are shown in Figure 6—figure supplement 1. In the revised manuscript, this analysis can be found in the section ‘Bhlhe22 is a novel regulator of granule cell development’ (pp. 20-21).

Additionally, we ensure all the images used to assess immunofluorescent co-staining with cell-type specific markers were taken at 20X magnification. We have also drawn neuroanatomical annotation lines to identify regions of interest in the developing cerebellum.

2. The P0 enhancer set is much smaller than the E12 or P9, raising questions about technical sensitivity. The authors make the conclusion that there are no P0 enriched enhancers, but this could also simply be a lack of sensitivity. The authors could address this by presenting information regarding the signal to noise and comparing their data to published ENCODE P0 hindbrain data.

We thank the reviewer for this suggestion and in the revised manuscript we present quality control metrics recommended by ENCODE to evaluate the success and sensitivity of our H3K27ac and H3K4me1 immunoprecipitations: fraction of reads in peaks (FRiP) and relative strand correlation (RSC) (Landt et al., 2012). These metrics are explained in greater detail in Landt et al., 2012 and briefly in the reviewed Methods section. All the H3K27ac and H3K4me1 samples generated in our study met ENCODE standards for a successful immunoprecipitation, with FRiP score higher than 1% and an RSC value greater than 0.8 (Landt et al., 2012). The P0 samples had slightly lower scores, on average, in both metrics than the E12 and P9 samples, suggesting that the lower number of P0 enriched enhancers could have been due to a lack of sensitivity compared to the other stages. Additionally, the number of peaks identified in both biological replicates for H3K27ac was greater in the P0 ENCODE hindbrain compared to our cerebellar dataset, also indicating issues in sensitivity. However, 91.7% of the H3K27ac peaks identified in the cerebellum were replicated in the hindbrain, indicating robust signals. The comparisons between our dataset and the ENCODE hindbrain data are presented in the Results section ‘Enhancer dynamics during cerebellar development’ (pp. 8).

3. The authors should compare their datasets to the mouse ENCODE developmental hindbrain datasets. Ideally, they would verify the patterns they found, showing that their datasets are high quality and correlated with the ENCODE data and then also review the reproducibility of the early vs. late enhancers.

We thank the reviewer for this suggestion and agree that identifying whether our findings are aligned with the ENCODE hindbrain datasets would support our enhancer predictions and provide greater confidence that our datasets are of high quality. As noted in Response #1 to Reviewer #1, our study focused on enhancer activity that is specific to the developing cerebellum while ENCODE looked at enhancer activity throughout the hindbrain. As suggested by this Reviewer (2), we compared chromatin profiles generated in our study to those in the hindbrain and identified that the majority of H3K27ac and H3K4me1 peaks overlapped between the two datasets. There is a subset of peaks specific to the cerebellum (~25% of peaks H3K27ac peaks at each stage) that highlights the identification of cerebellarspecific enhancers, bolstering the import of our analysis. These findings are detailed in a new section of the manuscript entitled “A subset of robust cerebellar enhancers is active in other developing brain regions” (pp. 10) and in Figure 2 —figure supplement 1.

In the revised manuscript, we also assessed whether the activity patterns/histone ChIP-seq signal of Early and Late active enhancers could be reproduced in mouse hindbrain samples. We found that in the hindbrain, Early and Late enhancers showed similar trends as we find in the cerebellum with a slight increase and decrease in activity, respectively. These findings are detailed in the new section entitled ‘A subset of robust cerebellar enhancers is active in other developing regions of the brain’ (pp. 10) and in Figure 2—figure supplement 2 and 3.

The authors can also use the P9 data to compare to the adult cerebellum in mouse ENCODE to see how similar these datasets are – it would be of interest if P9 is capturing novel enhancers.

This is an interesting point and we agree that it is relevant to this study to identify whether cerebellar enhancers demonstrate postnatally specific activity compared to the adult. The postnatal day (P)9 cerebellum is far from developmentally static and previous postnatal findings by Frank et al., (2015) indicate that there may be enhancers with activity specific to postnatal development compared to adult stages. To examine this in our data, we have now included a comparison of mouse cerebellar P9 histone profiles with the ENCODE adult profiles. We find that 30-40% of peaks were unique to P9. This indicates that a subset of the enhancers identified in our study may have activity specific to postnatal development. These findings can be found in the revised Results section ‘Enhancer identification during cerebellar development’ (pp. 7).

The authors could also look at dynamic chromatin accessibility to see if the H3K27ac + H3K4me1 signature here is more sensitive to active enhancers than accessibility alone (which would be expected and would further justify the value of the data here).

In our revised manuscript, we address this comment by looking at the overlap of our predicted enhancers with experimentally validated enhancers in the VISTA enhancer database. We find that 0.05% H3K27ac peaks overlapped with experimentally validated enhancers with activity in the hindbrain. This compares very favorably with the same overlap of cerebellar snATAC-seq and hindbrain ATAC-seq peaks, where only 0.005% and 0.006% of peaks overlapped with VISTA enhancers, respectively. This comparison is detailed in the revised Results section ‘Enhancer identification during cerebellar development’ (pp. 7) and in Figure 1—figure supplement 1. This analysis further justifies the value of our data as a resource for the identification of regulatory sequences that are likely to be active enhancers in the cerebellum.

4. It would be of interest to see how the enhancer sets identified here map to the recently published snATAC-seq datasets from mouse cerebellum across development. This publication was referenced in the manuscript, but not described – Sarropoulos Science 2021. More specifically, can the early and late enhancers here be assigned to shared or distinct cerebellar cell types?

We thank the reviewer for this suggestion and is aligned with this revision’s status as a resource paper.

We now compare our cerebellar histone profiles to single-cell profiles of chromatin accessibility (Sarrapolous et al., 2021). Mapping our cerebellar enhancers with this recently published snATAC-seq dataset we find that 90.1% of robust cerebellar enhancers identified in our bulk approach overlapped with CREs identified using snATAC-seq. We then predicted the cell type in which our robust cerebellar enhancers were active based on the cell-type annotation of the overlapping CREs (see new Figure 2 —figure supplement 4). Specifically, as suggested by this Reviewer, we used this dataset to quantify Early and Late active robust cerebellar enhancer accessibility in distinct cell types. The comparison of our findings with the snATAC-seq dataset provides cell-type resolution of chromatin accessibility for the majority of cerebellar enhancers. Importantly, our dataset provides novel evidence of enhancer activity by identifying enhancer-associated histone marks at thousands of CREs identified by snATAC-seq. These sequences are more likely to function as enhancers and are prime candidates for validation. The analysis and results can be found in a new Results section in our revised manuscript entitled ‘Comparison of robust cerebellar enhancers with chromatin accessibility identified in single cells’ (pp. 11).

5. The authors could intersect their target early and late genes with recent snRNA-seq data from developing human cerebellum (Aldinger Nat Neuro 2021).

The suggestion made by this reviewer is excellent and in the revised manuscript we now quantify the expression of Pax3 and Bhlhe22 in various cell-types in the developing human cerebellum using the snRNA-seq dataset previously generated by Aldinger et al., (2021). Pax3 and Bhlhe22 are Early and Late genes, respectively, discovered by our analysis with novel expression patterns and functions in the context of the developing cerebellum. The results of this analysis support the spatial expression patterns identified in our immunofluorescent staining for both Pax3 and Bhlhe22. These findings are presented in the revised Results sections ‘Early active enhancers are enriched for Pax3 binding sites, a novel marker for GABAergic cells’ (pp. 14-15) for Pax3 and ‘Bhlhe22 is a novel regulator of granule cell development’ (pp. 21) for Bhlhe22.

6. Some general issues with the target gene set analsyis:. For the k-means clustering of target genes presented in Figures 4 and 5, the authors should provide rationale for using the k value selected.

We agree that the rationale for choosing the k-value should be more clearly detailed and in our revised manuscript we emphasize the details of the analysis conducted to determine the number of clusters used for k-means clustering. To arrive at the k value (the number of clusters) for our k-means clustering analysis of Early and Late target genes we conducted an Elbow Analysis to determine the k value with minimal intra-cluster variation [total within-cluster sum of square (WSS)]. This method calculates the WSS and plots WSS as a function of the number clusters. The WSS decreases as the number of k decreases and the optimal k-value is the point at which the WSS first starts to diminish and the differences as k increases are minimal. When plotted on a line graph, this point is visible as an elbow.

We determined the optimal number of clusters is 4 for both Early and Late active enhancers (Figure 4 —figure supplement 1A). These points are provided in the ‘Co-expressed putative target genes are expressed in spatially distinct areas of the developing cerebellum’ (pp. 15) section in the revised manuscript.

Finally I'm not sure what the GeneRatio measure is supposed to capture. IT would be more straightforward to use an observed vs. expected fold enrichment, as counts can skew depending on number of genes and fold enrichment is standard in GO analysis.

We thank the reviewer for pointing out the pitfalls of the gene ratio measure and we have changed the GO enrichment analysis to an observed vs. expected fold enrichment as suggested by the reviewer. In the revised manuscript, this was conducted for plots displaying the results of our analyses in Figure 4D and 5D. Additionally, we’ve changed the x-axis of Figure 5D to ‘-log10 adjusted p-value’ and terms were sorted based on significant enrichment (ascending adjusted p-value). We think that this improves our findings by emphasizing the most highly enriched GO terms in each set of target genes.

7. Motif analysis: Motifs are often similar for families, so motif to TF matching is imprecise at best. Fine to select best guess TF as authors do, but good to also show motif logo in Figure 3A and include mention that the actual TF could be a different family member unless orthogonal evidence exists.

In agreement with these comments concerning motif families and TF matching, we have revised the manuscript to recognize that the motifs enriched in our cerebellar enhancers may represent the binding of TF families with similar motifs, not just a single best match. To address this point, we now provide a supplementary table of our results and the protein families for the TF with the best matching motif. In addition, we have revised Figure 3A to include the enriched motif logos which provides context for the TF matching. Finally, in the corresponding text of the revised manuscript, we emphasize that the TF bound to an enhancer may potentially be a member of the same protein family, as they typically have similar DNA binding sites. We now provide a Table and text in the Results section ‘Cerebellar enhancers are enriched for neural transcription factor binding sites in an age-dependent manner’ (pp. 12-13) that addresses this point.

8. Intersection with ASD variants is flimsy. The authors should use some presumed negative for non-brain and other brain GWAS to show that the ASD SNP findings are relevant e.g. have stronger signal than presumed negative or less relevant phenotypes.

We thank the reviewer for this recommendation and we believe that the addition of a negative control using non-brain GWAS demonstrates that our ASD SNP enrichment is relevant. As a negative control, we assessed whether our active cerebellar enhancers were enriched for SNPs associated with chronic kidney disease. We indeed find that there is no significant enrichment of chronic kidney disease SNPs, in contrast to ASD SNPs which were enriched in our cerebellar enhancers. This new addition to our analysis can be found in the section ‘Active cerebellar enhancers are enriched for common and de novo genetic variants associated with autism spectrum disorder.’ (pp. 23)

The overall enrichment with de novo mutations in enhancesr does not demonistrate any specific cerebellar relevance unless these enhancers are not active in other brain regions. Other studies have showed general enrichment of de novo mutations in brain enhancers. Not sure whether the rare CNV with the cerebellar gene is useful or convincing of relevance as the N is 1 and the gene expression outside the cerebellum is bot shown or discussed. All this analysis is fine, but shouldn't be treated as a major finding. I think the other parts of this study are much stronger, so de-emphasizing this would not hurt the paper.

While we agree with the reviewer that our enrichment analysis was not conducted exclusively for enhancers with cerebellar-specific activity, we believe that there is value in assessing the enrichment of de novo mutations associated with neurodevelopmental disorders such as ASD. Identifying the collective consequences of mutations in enhancers active in multiple brain regions may lead to a better understanding of the complex presentation of disorders like ASD. In accordance with the reviewer’s suggestion, we emphasize in our revised manuscript that this analysis is meant to be an example of the relevance of our dataset to human development and ASD and demonstrate how it may be utilized to functionally annotate variation associated with neurodevelopmental disorders. We believe that it would be beyond the scope of the current study to delve further into the potential impacts of the CNVs on cerebellar development, however it sets the stage for future studies.

9. The argument that MGI mouse cerebellum phenotype is evidence for validity needs to be compared to some random set in order to evaluate how strong the 98/2261 finding is- e.g. via permutation analysis selecting random sets of 2261 genes to establish background.

We thank the reviewer for this suggestion and we have now conducted a permutation analysis to assess whether the identification of 98 cerebellar genes within our list of 2261 target genes (2261) is evidence of statistical enrichment of cerebellar function. To do so, we generated 10,000 permutations of 2261 genes randomly selected from all genes expressed in the cerebellum and assessing the number of cerebellar genes in each permutation. We found that cerebellar genes were indeed enriched in robust cerebellar enhancer target genes, with an adjusted p-value = 0.0405. This analysis in now included in the section ‘Co-expressed putative target genes are expressed in spatially distinct areas of the developing cerebellum’ (pp. 17) and in Figure 4—figure supplement 2.

Text

1. The authors highlight their study as the first to characterize developmentally active enhancers, both in the abstract and discussion. This is an over-statement at best and the authors would do better to both remove this language and to directly deal with some of the other datasets, as described below.

We recognize that our study is not the first to characterize developmentally active enhancers and in line with the reviewers comment we now frame our work in the context of previous studies and emphasize that our study is best viewed as a valuable contribution to the catalog of active enhancers in the cerebellum. Our study provides a novel and unique dataset previously missing in the literature and can be used as a resource for the discovery of novel regulators of cerebellar development.

2. The introduction should include brief description of the other relevant studies, particularly the mouse ENCODE data and recent snATAC-seq and explain the novelty of approach here. The discussion should devote space to comparing the results from this study and the mouse ENCOE and snATAC-seq specifically.

We thank the reviewer for this suggestion and we now describe the relevant studies that have conducted assessments of enhancer activity in other brain regions (ENCODE, Gorkin et al., 2020) and chromatin accessibility at single cell resolution (snATAC-seq, Sarrapolous et al., 2021). We have included a brief description of these studies in the Introduction (pp. 3-4) of the revised manuscript. This emphasizes the novelty of our approach and results. We also provide a discussion on the comparison of our findings with the ENCODE and snATAC-seq datasets.

Reviewer #3:

Weaknesses:

The authors select three timepoints during cerebellar development for their ChIP-seq analysis. The rational for selecting these specific timepoints is not very clear, especially as they do not cover the entire period of cerebellar development in the mouse.

We address the reviewer’s concern about the rationale for investigating the time periods chosen for ChIP-seq in the revised version of our manuscript by including additional text in the section ‘Enhancer identification during cerebellar development’ (pp. 5) in our results, noting the importance of these periods for the various stages of development across the principle cell types of the cerebellum. Thus, our study is a novel genome-wide assessment of enhancer-associated histone marks throughout key points in the developing cerebellum.

While the resource itself is very valuable, the presentation of example genes that the authors extract from their analyses (Pax3, Bhlhe22, and Cdc42bpb) is rather short and lacking in detail. For example, the authors could have further strengthened their findings by utilising recent single-cell datasets, in particular the single-cell transcriptional atlas of the developing cerebellum published by Carter et al.

This an important point made by the reviewer, particularly in the context of this being a resource paper where the value lies not only in our initial analysis but in it’s utility when combined with other existing datasets. We have strengthened our analyses by comparing our data with previous scRNA-seq datasets quantified in the developing mouse (Carter et al., 2018) and human (Aldinger et al., 2021) cerebellum. Comparisons were conducted specifically for the genes with novel functions in the context of cerebellar development Pax3 and Bhlhe22, where we examine their expression in the cell-types of the cerebellum. We also take the initial steps in validating these findings through immunofluorescent staining at several stages of development and co-staining with cell-type specific markers. We demonstrate the utility of our data as a resource and provide a proper framework for investigating genes with novel functions in the cerebellum. These analyses can be found in sections ‘Early active enhancers are enriched for Pax3 binding sites, a novel marker for GABAergic cells’ (pp. 14-15) and ‘Bhlhe22 is a novel regulator of granule cell development’ pp (21).

The authors use in silico tools to identify transcription factor binding sites enriched in the identified active enhancers. No experimental validation is provided for this analysis. This would have been particularly relevant for transcription factors that have not been previously well-described in cerebellar development including Pax3, which the authors present as a novel marker for GABAergic progenitors.

Based on the recommendations by the Editors of this study and the fact that the corresponding data will now be available as a Resource article, we believe conducting the experimental validation of transcription factor binding is outside the scope of this revision. The goal of our motif enrichment analysis is to predict the regulators of cerebellar enhancer activity. The identification of these TFs provides further insight into the cell-types in which these enhancers may be active and may help predict the cellular processes regulated by cerebellar enhancers. We do demonstrate how our predictions can be validated through a preliminary set of staining experiments for Pax3, a TF enriched in our cerebellar enhancers with unknown function in the cerebellum. This sets the stage for further investigation of Pax3 function, as well as for other novel regulators of cerebellar development identified by our analyses.

Furthermore, it would have been helpful to go into further detail and to provide some predicted candidate genes that might be regulated by Pax3 in GABAergic progenitors.

We thank the reviewer for this suggestion and agree that providing further detail and characterizing the enhancers with Pax3 motifs and their putative target genes may provide a more detailed description of Pax3 function in cerebellar development. To bolster our Pax3 findings, we identified the cerebellar enhancers with a Pax3 motif and their target genes, as per the reviewer’s suggestion. We found that 923 Early active robust cerebellar enhancers with a Pax3 binding motif were highly correlated with a putative target gene. A GO enrichment analysis determined that predicted Pax3 target genes were enriched for biological processes such as ‘regionalization’, ‘pattern specification process’ and ‘neural precursor cell proliferation’ which is in line with the results of our immunofluorescent analysis.

Additionally, we identified several target genes associated with the function of neural progenitors in the cerebellar VZ (ventricular zone) and the specification and differentiation of GABAergic cell types such as Ascl1, Neurog2, Sox11 and Ctnna2 (Figure 4 —figure supplement 2). These findings are detailed in the section ‘Co-expressed putative target genes are expressed in spatially distinct areas of the developing cerebellum’ (pp. 17).

Figure 3: Very limited expression data is shown for Pax3. It would have been helpful to show additional immunostainings at least for P9 (same timepoint as used for the ChIP-seq analysis) if not for further timepoints during cerebellar development. Moreover, it would have been really interesting to compare this to the data available on Pax3 from Carter et al. The authors mention this in their discussion; it would be helpful to include this data in the relevant figure and also to comment on how expression of Pax3 changes during the course of cerebellar development.

In agreement with the reviewer’s suggestions, we have extended the time course of Pax3 immunofluorescent staining to P3 and P9. We found that expression patterns found in P0 extend to P3 but by P9 we cannot detect Pax3 expression. These findings can be found in Figure 3 and Figure 3—figure supplement 1.

I am not sure that I agree with the authors about the fact that 98 out of 2261 identified target genes result in a cerebellar phenotype when knocked out particularly demonstrates the validity of the high throughput approach (page 12). Is this a significant number?

Addressed. See response to Reviewer #2 (Recommendations for the authors) Point 12.

The use of CAGE-sequencing data in this manuscript is very unclear. The authors refer to a previous manuscript (Zhang et al., 2018); however, in this study reports transcriptomic data from only three timepoints (E13, E15, E18) and not 12 timepoints as mentioned in by the authors? The authors should better describe which data they used. Also, given that the CAGE-Seq data is from bulk-sequencing, it would have been much more informative to look at gene expression of their candidate transcription factors in the Carter et al., single-cell dataset. Furthermore, this would have allowed the authors to look beyond the P9 timepoint and thus for a more complete gene expression analysis over cerebellar development. This is also true for Figures 1E,G.

We thank the reviewer for bringing up the CAGE data issue as confusing. This confusion was engendered by the use of an incorrect citation (Zhang et al., 2018). Our apologies, the correct reference is: Ha, T. J., Zhang, P. G. Y., Robert, R., Yeung, J., Swanson, D. J., Mathelier, A.,... Goldowitz, D. (2019). Identification of novel cerebellar developmental transcriptional regulators with motif activity analysis. BMC Genomics, 20(1), 718. doi:10.1186/s12864-019-6063-9. This has been corrected.

The approach of using the enhancer data set to identify novel target genes not previously identified in cerebellar development is unclear to me. Using their described filtering approach, the authors describe genes that differentially regulated by Atoh1. Could these not have been identified from the original dataset by Klisch et al., 2001? What does the initial identification of the enhancers add here?

The Klisch et al., (2001) dataset produced a large number of differentially expressed genes (DEGs) at P7. While ground-breaking in their analysis, it left two big gaps: (1) other time points were not assessed where many developmental events are occurring and (2) this produced about 2000 DEGs making it difficult to identify specific novel gene candidates that are critical for the postnatal differentiation of granule cells. To address these gaps, we (1) analyze enhancer activity and gene expression over embryonic and postnatal development (E12, P0, P9), and (2) due to our analysis of enhancers, which are likely to be regulating genes important for specification and differentiation, this allowed us discover novel TFs not previously identified in cerebellar development. We use our enhancer putative target gene dataset in conjunction with DEGs in the Atoh1 mouse to identify Bhlhe22, a novel TF that we find is important for postnatal granule cell migration. This analysis can be found in the section of the Results entitled ‘Bhlhe22 is a regulator of postnatal granule cell development’ (pp. 18).

Figure 6: As above, the authors should include analysis of Bhlhe22 expression in the Carter et al., single-cell dataset.

Addressed.

The authors normalize gene expression to co-transfected GFP. Can the authors demonstrate that cells that are transfected with GFP were also transfected with the siRNA construct?

We thank the reviewer for the suggestions and the aesthetic adjustments were made by increasing the font size of Figure 7 and bringing greater clarity to our results. We have provided quantitative RT-qPCR data of the reduction of Bhlhe22 in our siRNA treated cultures as evidence of a successful knockdown of Bhlhe22, and believe that immunostaining is not necessary in this regard. Co-transfection with a GFP plasmid and siRNA molecules has been previously been demonstrated as a means for identifying transfected cells in primary neuronal cultures (Tong et al., 2011; Wang et al., 2013; Schmidt et al., 2009; Cosker et al., 2008; Tonges et al., 2011). The co-transfection efficiency of this method has been demonstrated in human embryonic stem cells by Moore et al., 2010 (DOI:10.1186/scrt23), where siRNA molecules targeting GFP were co-transfected with a GFP plasmid by nucleofection, which is the transfection method used in our study. Fluorescence microscopy showed that GFP expression was visibly reduced when GFP-targeting siRNA was co-transfected FACS analysis showed that an siRNA to DNA ratio of 3:1 reduced the percentage of cells expressing GFP by 71.75% when compared to controls indicating high co-transfection efficiency. This has also been demonstrated in nucleofection protocol documentation provided by Lonza, which is the is manufacturer of the Amaxa Nucleofector 2b which is used in our study for transfection (https://bioscience.lonza.com/lonza_bs/US/en/document/21126). In this example, three different cell lines are transfected with GFP plasmid and GFP targeting siRNA. A reduction in GFP expression was observed in all cell lines, indicating high co-transfection efficiency. The co-transfection of a GFP plasmid with siRNA molecules is now common practice in the field, and that cells that are successfully electroporated with GFP typically are also transfected with the siRNA molecules.

Can the authors comment on the evolutionary conservation of their identified enhancer elements? This is particularly relevant for their analysis of variants associated with ASD. Are the loci of interest conserved between mouse and human and thus, do enhancer activity and transcription factor occupancy occur at orthologous regions distal to the same gene promoters?

All of the enhancers analyzed in our enrichment analysis of ASD variation are by definition conserved between mouse and human since the relevant human sequences used for this analysis were converted from the mouse to human genome using the UCSC Genome Browser LiftOver tool. This computational tool is called LiftOver because it converts coordinates in a reference assembly to another genome build. This is commonly used to perform cross species mapping, which takes coordinates of sequences in one species to identify the corresponding coordinates in other species. The majority of robust cerebellar enhancers (6630/7024) were converted from the mouse genome (mm9) to the human genome (hg38). 89.6% of robust cerebellar enhancers with a putative target gene are located at orthologous regions distal to the same gene promoters indicating conservation between mouse and humans and bringing credence to our analysis. The results of this analysis can be found in Supplementary Table 4 and in the revised manuscript in the section ‘Active cerebellar enhancers are enriched for common and de novo genetic variants associated with autism spectrum disorder.’ (pp. 23)

For the identification of CDC42BPB, it would have been nice to link this back to expression data over time and to provide more detail on the potential transcription factors involved.

We thank the reviewer for this suggestion and agree that it is relevant to our study to examine CDC42BPB expression in the developing mouse expression time course. In our revised manuscript, we describe CDC42BPB expression throughout embryonic and postnatal development. We also profile expression in the cell-types of the cerebellum using scRNA-seq in mouse and human samples (Carter et al., 2018; Aldinger et al., 2021). Finally, we confirm that the enhancer predicted to regulate CDC42BPB is associated with the same gene in the mouse and located on the same chromosome and at a similar distal genomic location in the human. (pp. 25)

Finally, given that this is a very important dataset for the scientific community, I wish it would be easier for others to interact with the data rather than just depositing data files into GEO. It would have also been really helpful to provide supplementary tables listing enhancer elements and associated transcription factors as well as nearby genes that the authors have identified. Currently, it is not transparent which lists of genes etc were used for the various analyses shown.

We agree with the reviewer that this is a very important dataset for the scientific community and have made efforts to make our data more accessible. An online resource

(https://goldowitzlab.shinyapps.io/developing_mouse_cerebellum_enhancer_atlas/) was created that can be used to browse, curate and export information regarding cerebellar enhancers and their putative target genes identified in this study. In the revised manuscript, we have also provided additional supplementary tables listing enhancer elements, their associated transcription factors and their associated target genes.

https://doi.org/10.7554/eLife.74207.sa2

Article and author information

Author details

  1. Miguel Ramirez

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Validation, Investigation, Visualization, Methodology, Writing - original draft, Project administration, Writing - review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-6315-4291
  2. Yuliya Badayeva

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Data curation, Software, Formal analysis, Validation, Investigation, Visualization, Methodology, Writing - original draft, Writing - review and editing
    Competing interests
    No competing interests declared
  3. Joanna Yeung

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Formal analysis, Supervision, Validation, Methodology, Project administration
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-0551-5305
  4. Joshua Wu

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Formal analysis, Validation, Investigation, Visualization, Methodology, Writing - review and editing
    Competing interests
    No competing interests declared
  5. Ayasha Abdalla-Wyse

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Data curation, Software, Visualization
    Competing interests
    No competing interests declared
  6. Erin Yang

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Validation, Investigation, Visualization, Methodology
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-5629-2362
  7. FANTOM 5 Consortium

    RIKEN, Wako, Japan
    Contribution
    Data curation, Resources
    Competing interests
    No competing interests declared
  8. Brett Trost

    The Centre for Applied Genomics, The Hospital for Sick Children, Toronto, Canada
    Contribution
    Resources, Data curation, Software, Formal analysis, Validation, Investigation, Visualization, Methodology, Writing - original draft, Writing - review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-4863-7273
  9. Stephen W Scherer

    The Centre for Applied Genomics, The Hospital for Sick Children, Toronto, Canada
    Contribution
    Resources, Data curation, Supervision, Funding acquisition, Investigation, Writing - original draft, Writing - review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-8326-1999
  10. Daniel Goldowitz

    1. Centre for Molecular Medicine and Therapeutics, BC Children’s Hospital Research Institute, Vancouver, Canada
    2. University of British Columbia, Vancouver, Canada
    Contribution
    Conceptualization, Resources, Data curation, Supervision, Funding acquisition, Investigation, Methodology, Writing - original draft, Project administration, Writing - review and editing
    For correspondence
    dang@cmmt.ubc.ca
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-4756-4017

Funding

NSERC Discovery Award

  • Daniel Goldowitz

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Ethics

All studies were conducted according to the protocols approved by the Institutional Animal Care and Use Committee and the Canadian Council on Animal Care at the University of British Columbia.

Senior Editor

  1. Kathryn Song Eng Cheah, University of Hong Kong, Hong Kong

Reviewing Editor

  1. Genevieve Konopka, University of Texas Southwestern Medical Center, United States

Reviewers

  1. Genevieve Konopka, University of Texas Southwestern Medical Center, United States
  2. Alex S Nord, University of California, Davis, United States

Publication history

  1. Received: September 25, 2021
  2. Preprint posted: September 29, 2021 (view preprint)
  3. Accepted: August 5, 2022
  4. Accepted Manuscript published: August 9, 2022 (version 1)
  5. Version of Record published: August 23, 2022 (version 2)

Copyright

© 2022, Ramirez et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 649
    Page views
  • 269
    Downloads
  • 1
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Miguel Ramirez
  2. Yuliya Badayeva
  3. Joanna Yeung
  4. Joshua Wu
  5. Ayasha Abdalla-Wyse
  6. Erin Yang
  7. FANTOM 5 Consortium
  8. Brett Trost
  9. Stephen W Scherer
  10. Daniel Goldowitz
(2022)
Temporal analysis of enhancers during mouse cerebellar development reveals dynamic and novel regulatory functions
eLife 11:e74207.
https://doi.org/10.7554/eLife.74207

Further reading

    1. Developmental Biology
    Javier Solivan-Rivera, Zinger Yang Loureiro ... Silvia Corvera
    Research Article Updated

    Mechanisms that control ‘beige/brite’ thermogenic adipose tissue development may be harnessed to improve human metabolic health. To define these mechanisms, we developed a species-hybrid model in which human mesenchymal progenitor cells were used to develop white or thermogenic/beige adipose tissue in mice. The hybrid adipose tissue developed distinctive features of human adipose tissue, such as larger adipocyte size, despite its neurovascular architecture being entirely of murine origin. Thermogenic adipose tissue recruited a denser, qualitatively distinct vascular network, differing in genes mapping to circadian rhythm pathways, and denser sympathetic innervation. The enhanced thermogenic neurovascular network was associated with human adipocyte expression of THBS4, TNC, NTRK3, and SPARCL1, which enhance neurogenesis, and decreased expression of MAOA and ACHE, which control neurotransmitter tone. Systemic inhibition of MAOA, which is present in human but absent in mouse adipocytes, induced browning of human but not mouse adipose tissue, revealing the physiological relevance of this pathway. Our results reveal species-specific cell type dependencies controlling the development of thermogenic adipose tissue and point to human adipocyte MAOA as a potential target for metabolic disease therapy.

    1. Developmental Biology
    Seok Hee Lee, Xiaowei Liu ... Paolo F Rinaudo
    Research Article Updated

    In vitro fertilization (IVF) has resulted in the birth of over 8 million children. Although most IVF-conceived children are healthy, several studies suggest an increased risk of altered growth rate, cardiovascular dysfunction, and glucose intolerance in this population compared to naturally conceived children. However, a clear understanding of how embryonic metabolism is affected by culture condition and how embryos reprogram their metabolism is unknown. Here, we studied oxidative stress and metabolic alteration in blastocysts conceived by natural mating or by IVF and cultured in physiologic (5%) or atmospheric (20%) oxygen. We found that IVF-generated blastocysts manifest increased reactive oxygen species, oxidative damage to DNA/lipid/proteins, and reduction in glutathione. Metabolic analysis revealed IVF-generated blastocysts display decreased mitochondria respiration and increased glycolytic activity suggestive of enhanced Warburg metabolism. These findings were corroborated by altered intracellular and extracellular pH and increased intracellular lactate levels in IVF-generated embryos. Comprehensive proteomic analysis and targeted immunofluorescence showed reduction of lactate dehydrogenase-B and monocarboxylate transporter 1, enzymes involved in lactate metabolism. Importantly, these enzymes remained downregulated in the tissues of adult IVF-conceived mice, suggesting that metabolic alterations in IVF-generated embryos may result in alteration in lactate metabolism. These findings suggest that alterations in lactate metabolism are a likely mechanism involved in genomic reprogramming and could be involved in the developmental origin of health and disease.