Author response: A gene regulatory network for neural induction
Year: 2022
Type: peer-review
Abstract: Article Figures and data Abstract Editor's evaluation Introduction Results Discussion Materials and methods Data availability References Decision letter Author response Article and author information Metrics Abstract During early vertebrate development, signals from a special region of the embryo, the organizer, can redirect the fate of non-neural ectoderm cells to form a complete, patterned nervous system. This is called neural induction and has generally been imagined as a single signalling event, causing a switch of fate. Here, we undertake a comprehensive analysis, in very fine time course, of the events following exposure of competent ectoderm of the chick to the organizer (the tip of the primitive streak, Hensen’s node). Using transcriptomics and epigenomics we generate a gene regulatory network comprising 175 transcriptional regulators and 5614 predicted interactions between them, with fine temporal dynamics from initial exposure to the signals to expression of mature neural plate markers. Using in situ hybridization, single-cell RNA-sequencing, and reporter assays, we show that the gene regulatory hierarchy of responses to a grafted organizer closely resembles the events of normal neural plate development. The study is accompanied by an extensive resource, including information about conservation of the predicted enhancers in other vertebrates. Editor's evaluation In this manuscript, Trevers and colleagues undergo a detailed genome-wide exploration of the mechanisms of neural induction in chick embryos. They describe the gene regulations governing the patterning of extra-embryonic ectoderm into neural ectoderm upon the graft of an early Hensen's node ectopically, an assay for neural induction and neural commitment. The data are assembled into a Gene Regulatory Network of 175 transcription factors and their projected interactions, based on a fine-scale temporal analysis. This study will be an important resource for the field of neural induction. https://doi.org/10.7554/eLife.73189.sa0 Decision letter eLife's review process Introduction One of the most influential studies in developmental biology was the discovery, 100 years ago, that a small region of the vertebrate embryo, named the ‘organizer’, can induce ectodermal cells that do not normally contribute to the neural plate to form a complete, patterned nervous system (Spemann, 1921; Spemann and Mangold, 1924). In amphibians, where these experiments were initially conducted, the ‘organizer’ resides in the dorsal lip of the blastopore. A few years later, Waddington demonstrated that an equivalent region exists in birds (duck and chick) and mammals (rabbit) (Waddington, 1933; Waddington and Schmidt, 1933; Waddington, 1934; Waddington, 1936; Waddington, 1937): the tip of the primitive streak, a structure known as Hensen’s node (Hensen, 1876). This interaction between the organizer and the responding ectoderm, which causes the latter to acquire neural plate identity, has been termed ‘neural induction’ (Spemann and Mangold, 1924; Nieuwkoop, 1952; Gallera, 1971b; Saxén, 1980; Gurdon, 1987; Storey et al., 1992; Stern, 2005). Neural induction has often been imagined as a single event, ‘switching’ the fate of the responding tissue from non-neural to neural. But it is clear from timed grafting and subsequent removal of organizer transplants that the responding ectoderm requires about 12 hr of exposure to the organizer to acquire neural identity in a stable way (‘commitment’) (Gallera and Ivanov, 1964; Gallera, 1971b; Streit et al., 1998). Recent work has also revealed that the expression of many genes changes over time after grafting an organizer, suggesting that the process has considerable complexity (Streit et al., 1997; Streit et al., 1998; Streit and Stern, 1999; Streit et al., 2000; Sheng et al., 2003; Stern, 2005; Albazerchi and Stern, 2007; Papanayotou et al., 2008; Gibson et al., 2011; Pinho et al., 2011; Papanayotou et al., 2013; Trevers et al., 2018). What happens during this 12 hr period? Is it possible to define distinct steps, and perhaps identify the molecular events that represent ‘induction’ and ‘commitment’? Surprisingly, given the long time since the discovery of neural induction, these questions have hardly been addressed. To begin to answer them requires a precise approach to identify and model the key interactions between genes and transcriptional regulators inside the responding cells that accompany their responses to signals from the organizer over time. Grafting an organizer to a region of ectoderm that does not normally contribute to neural tissue (but is competent to do so) provides the opportunity to study the progress of neural induction relative to ‘time-zero’, the moment when the organizer is first presented to the tissue. It also allows for the separation of neural inductive events from other processes that occur adjacent to the normal neural plate (e.g. mesendoderm formation and patterning), because sites that are remote enough from the normal neural plate and are competent to respond to an organizer graft can generate a complete neural tube without induction of mesoderm, and without recruiting cells from the host’s neural plate (Hornbruch et al., 1979; Dias and Schoenwolf, 1990; Storey et al., 1992). Here, we have taken advantage of these properties, together with recent major technological advances in transcriptomics and epigenetic analysis, to dissect the molecular events that accompany neural induction in the chick embryo in fine time course and to generate the first detailed gene regulatory network (GRN) for this process. The GRN comprises 175 transcriptional regulators and 5614 predicted interactions between them over the course of neural induction. We then compare the spatial and temporal properties of these changes with development of the normal neural plate of the embryo using in situ hybridization, single-cell RNA-sequencing (scRNAseq) and reporter assays. The latter also allow us to test the activity of some of the key gene regulatory elements. We present a comprehensive resource allowing visualization and querying of all of these interactions and regulatory elements on a genome-wide level. Together, our study offers a global view of the genetic hierarchy of transcriptional regulatory interactions during neural induction and normal neural plate development over time. Results Transcriptional profiling identifies responses to neural induction in time course In the chick embryo, a graft of Hensen’s node to the inner area opaca, an extraembryonic region of competent non-neural ectoderm (Gallera and Ivanov, 1964; Dias and Schoenwolf, 1990; Storey et al., 1992; Streit et al., 1998), induces the formation of a mature, patterned neural tube after about 15 hr of culture from the host ectoderm (Figure 1A). To characterize the events over this period, the induction of several neural markers was assessed. OTX2 is first expressed in the pre-primitive-streak stage epiblast at EGKXII-XIII (Albazerchi and Stern, 2007; Pinho et al., 2011) and is induced by 7/7 grafts after 3 hr (Figure 1B–F). SOX2 is first expressed at HH4+/5 (Rex et al., 1997; Streit et al., 1997; Streit et al., 1998; Sheng et al., 2003) and requires 9 hr for induction (8/8) (Figure 1G–K). SOX1 starts to be expressed weakly in the neural plate around HH7-8 (Stavridis et al., 2010; Uchikawa et al., 2011) and is induced by 50% of grafts after 12 hr (Figure 1L–P). This time course confirms that a sequence of events occurs in response to a grafted node over 0–12 hr (Figure 1—figure supplement 1A), culminating in the acquisition of neural plate/tube identity (Streit et al., 1997; Streit et al., 1998; Streit and Stern, 1999; Streit et al., 2000; Sheng et al., 2003; Stern, 2005; Albazerchi and Stern, 2007; Papanayotou et al., 2008; Gibson et al., 2011; Pinho et al., 2011; Papanayotou et al., 2013; Trevers et al., 2018). Figure 1 with 1 supplement see all Download asset Open asset Transcriptional profiling identifies responses to neural induction in time course. (A) Hensen’s node (HN) was grafted from HH3-4 donors to the inner area opaca (yellow) of hosts at the same stage. An ectopic neural tube expressing SOX2, OTX2, and HOXB1 is induced after 15 hr of culture. (B–F) Expression of neural markers compared to their time course of induction by a node graft: OTX2; first expressed in pre-streak epiblast (EGKXII-XIII) and induced by grafts after 3 hr. (G–K) SOX2; first expressed in the neural plate at HH5-6 and induced after 9 hr. (L–P) SOX1; first expressed in the forming neural tube at HH7-8 and induced after 12 hr. (Q) Identifying transcriptional responses to a node graft. HN was grafted from HH4- donors to HH4- hosts. The HN graft was removed and underlying ‘induced’ (IN) and contralateral ‘uninduced’ (UN) ectoderm isolated after 5, 9, or 12 hr. Uninduced ‘0 hr’ ectoderm from HH4- embryos was also dissected. Samples were analysed by RNA-sequencing (RNAseq). (R–T) Induced and corresponding uninduced tissues were compared at each time point to identify differentially expressed transcripts. Volcano plots show upregulated (green) and downregulated (blue) transcripts. (U–W) Venn diagrams of 482 genes encoding transcriptional regulators (U) that are upregulated (V) or downregulated (W) at different time points. Scale bars: B: 1mm; C: 250 μm; D: 250 μm; E: 250 μm. These scale bars apply to all other figures with embryos at equivalent stages throughout the paper. Figure 1—source data 1 RNAseq galGal3 analysis key and differentially expressed transcripts. https://cdn.elifesciences.org/articles/73189/elife-73189-fig1-data1-v2.xlsx Download elife-73189-fig1-data1-v2.xlsx Figure 1—source data 2 RNAseq galGal4 analysis key and differentially expressed transcripts. https://cdn.elifesciences.org/articles/73189/elife-73189-fig1-data2-v2.xlsx Download elife-73189-fig1-data2-v2.xlsx Figure 1—source data 3 Nanostring probe codeset. https://cdn.elifesciences.org/articles/73189/elife-73189-fig1-data3-v2.xlsx Download elife-73189-fig1-data3-v2.xlsx To uncover genes whose expression changes over this period, we performed RNAseq of the responding tissue at three time points following a graft of Hensen’s node: 5 hr (to identify early ‘pre-neural’ responses), 9 hr (when SOX2 expression defines neural specification), and 12 hr, as the host tissue starts to express SOX1. Hensen’s node was grafted from HH4- chick donors to a region of competent epiblast (inner area opaca) of chick hosts at the same stage. Embryos were cultured for the desired period of time, after which the graft was removed and the underlying ‘induced’ tissue was collected. Control (‘uninduced’) ectoderm from the corresponding region on the contralateral side of the same embryos (Figure 1Q) was also collected, as well as competent area opaca at HH4-, representing a ‘0 hr’ starting point for the time course. When ‘induced’ and ‘uninduced’ reads at each time point are compared by DESeq analyses, 8673 differentially expressed transcripts were identified (see volcano plots, Figure 1R–T). Of these, 4130 were upregulated (enriched in ‘induced’ tissue) and 4543 were downregulated (depleted in ‘induced’ tissue) relative to the ‘uninduced’ counterpart. To construct a GRN we focused on 482 transcription factors and chromatin modifiers that are differentially expressed within these data. Of these ‘transcriptional regulators’, 255 are upregulated, 203 are downregulated, and 24 have more complex expression dynamics (Figure 1U and Figure 1—figure supplement 1B). Grouping transcripts based on the timing of their response (Figure 1V–W) reveals that 120 are differentially expressed throughout (53 upregulated and 67 downregulated), while others are associated with particular time points. Therefore, transcriptional responses to organizer grafts involve changes in the expression of many genes over a relatively short period. Due to this complexity, we chose to increase the time resolution of the analysis. NanoString nCounter was used to quantify the gene expression changes of transcriptional regulators at six time points (1, 3, 5, 7, 9, and 12 hr after the node graft) in ‘induced’ tissue compared to ‘uninduced’ ectoderm. By consolidating the data from RNAseq and NanoString, a set of refined expression profiles was established for 213 transcriptional regulators – 156 that are enriched and 57 that are depleted in induced tissues (Figure 3—source data 1). These represent the core components of our GRN. Epigenetic changes identify chromatin elements associated with neural induction We next sought to identify the regions of non-coding chromatin that govern these transcriptional responses to signals from the organizer. Histone modifications can regulate gene expression by altering chromatin structure. For example, H3K27ac is associated with actively transcribed genes and their enhancers, whereas transcriptionally inactive regions are often marked by H3K27me3 (Tiwari et al., 2008; Heintzman et al., 2009; Creyghton et al., 2010; Kharchenko et al., 2011; Rada-Iglesias et al., 2011; Tolhuis et al., 2011; Zentner et al., 2011; Bonn et al., 2012). Therefore, ChIPseq and ATACseq (Buenrostro et al., 2013; Buenrostro et al., 2015; Corces et al., 2017) were performed on induced and uninduced ectoderm following 5, 9, or 12 hr of a node graft, to detect histone and chromatin conformation changes during neural induction. H3K27ac and H3K27me3 enriched chromatin (ChIPseq) were identified genome-wide by comparison to matched genomic input samples. Chromatin sites were then categorized according to their histone signatures across induced and corresponding uninduced tissues at each time point (Figure 2—figure supplement 1A). Sites that become acetylated and/or demethylated in induced tissues, compared to uninduced, were considered to undergo ‘activation’ (Indices 1–3) while those that become deacetylated and/or methylated undergo ‘repression’ (Indices 4–6). Chromatin marked by both H3K27ac (activation) and H3K27me3 (repression) marks in either tissue were described as ‘poised’ (Indices 7–12). The remaining indices (13–16) define sites that do not change marks between the two conditions. Chromatin sites undergoing activation are enriched for H3K27ac marks in induced tissues at each time point (Figure 2B–D). Sites of repression are more varied: those belonging to this category at 5 hr lose acetylation, whereas several sites become methylated in induced tissues at 9 and 12 hr. Very few sites are ‘poised’ at any time. Overall, as neural induction progresses, the number of sites undergoing activation increases, along with a reduction in those that do not change state. In contrast, ATACseq suggests that there is little difference between induced and uninduced samples in terms of chromatin accessibility (Figure 2—figure supplement 2) at 5 hr, but the differences become more marked at later time points. At 12 hr, the ATACseq profile resembles that of the endogenous neural plate of the embryo (Figure 2—figure supplement 2). Figure 2 with 3 supplements see all Download asset Open asset Epigenetic changes identify chromatin elements associated with neural induction. (A) Hensen’s node (HN at HH4-) was grafted to hosts of the same stage. Embryos were cultured for 5, 9, or 12 hr before HN was removed and induced (IN) and contralateral uninduced (UN) ectoderm collected. ChIPseq was performed for H3K27ac and H3K27me3. (B–D) Putative regulatory sites were predicted according to the H3K27 profiles of IN and UN tissues at each time point (see Figure 2—figure supplement 1). They include sites undergoing ‘activation’ or ‘repression’, being ‘poised’, or showing no change. Heat maps illustrate the enrichment of H3K27ac (blue) and H3K27me3 (purple) in IN and UN tissues within ±2.0 kb from the peak centre (PC) for each annotated group. Graphs plot the distributions of H3K27ac and H3K27me3 enrichment for each group. Figure 2—source data 1 ChIPseq tissue dissociation details. https://cdn.elifesciences.org/articles/73189/elife-73189-fig2-data1-v2.xlsx Download elife-73189-fig2-data1-v2.xlsx Figure 2—source data 2 ATACseq indices. https://cdn.elifesciences.org/articles/73189/elife-73189-fig2-data2-v2.xlsx Download elife-73189-fig2-data2-v2.xlsx Next, we focused on epigenetic changes associated with the 213 genes encoding transcriptional regulators for which refined expression profiles were established by NanoString analysis (see above). ‘Constitutive’ CTCF-bound sites (putative insulators) flanking these genes were obtained from chicken CTCF-ChIPseq data (Khan et al., 2013; Kadota et al., 2017). They were used to predict the boundaries within which we considered H3K27 marks as being associated with the gene of interest; these regions (‘loci’) are up to 500 kb in length and may contain several genes (Figure 2—figure supplement 1B). Within each locus, the H3K27ac or H3K27me3 enriched peaks were categorized as before (Figure 2—figure supplement 1A). The 213 selected loci contain a total of 6971 sites that change in response to a node graft (Indices 1–10; Figure 2—figure supplement 1A). We noticed that the flanking regions (proximal promoters) of transcriptional regulators whose expression increases are enriched for H3K27ac marks in induced tissue (Figure 2—figure supplement 3). In contrast, there is no difference between H3K27ac and H3K27me3 marks around the flanking regions of repressed transcriptional regulators at 5 and 9 hr – they only become enriched for methylation at 12 hr. Our analysis therefore identifies putative regulatory elements that are controlled epigenetically during the process. A GRN for neural induction Having identified transcriptional responses to signals from the organizer and their accompanying chromatin changes, we combined them to construct a GRN to describe the time course of these events and to illustrate predicted interactions between transcriptional regulators. The open source platform BioTapestry (Longabaugh et al., 2005; Longabaugh et al., 2009; Paquette et al., 2016) (http://www.biotapestry.org) is used extensively to provide a visually intuitive representation of developmental networks (Davidson, 1990; Arnone and Davidson, 1997; Davidson et al., 1998; Betancur et al., 2010; Simões-Costa and Bronner, 2015; Thiery et al., 2020). Importantly, this software allows dynamic changes in the interactions between regulators and their target genes to be represented. A BioTapestry model for regulatory gene interactions To integrate these complex time course data into a GRN that models the interactions during neural induction, we developed a custom computational pipeline (Figure 3). The 213 transcriptional regulators previously identified were selected as candidate components of the network because of their responses to a node graft in both RNAseq and NanoString data, and their associated chromatin changes in ChIPseq. Interactions were predicted by screening putative regulatory sites that belong to Indices 1–3 and 7 (activation) (Figure 2—figure supplement 1C–E) for the presence of putative binding sites for transcription factors that are differentially expressed in the dataset, and expressed at the appropriate time point to fit with a putative regulatory role (i.e. simultaneously or just before) for the target locus. Figure 3 Download asset Open asset Constructing a BioTapestry model for regulatory gene interactions: computational workflow. Computational pipeline to combine transcriptomic and epigenetic time course data, to build a gene regulatory network (GRN) for neural induction. The output data are available to view in the UCSC browser. Figure 3—source data 1 Integrated GRN gene expression profile. https://cdn.elifesciences.org/articles/73189/elife-73189-fig3-data1-v2.xlsx Download elife-73189-fig3-data1-v2.xlsx We then used the regulatory rules shown in Figure 4A–B to predict positive or negative regulatory events between regulators and their putative targets, overlaying the predicted transcription factor binding profiles with the time course of expression from RNAseq and NanoString. At each time point, changes in expression of a putative regulator that could lead to the same change in a candidate target are represented as positive interactions. Negative, or inhibitory, interactions are predicted when changes in expression of the regulator could cause the opposite change in a downstream target. Genes lacking both regulatory inputs and outputs with other genes in the network were excluded. Figure 4 Download asset Open asset Constructing a BioTapestry model for regulatory gene interactions: the regulatory logic. (A) Network interactions were modelled using chromatin sites that belong to Indices 1–7 and the defined rules between regulators and their potential targets. At sites that become active in induced tissues (Indices 1, 3, and positive interactions are when a putative input and target are and negative interactions are modelled when an input and target have expression These rules are at sites that undergo repression in induced tissues (Indices 5, and that interactions are not predicted when the potential regulators of these repressed sites are not expression and chromatin profiles were combined to predict interactions between gene regulatory network (GRN) RNA-sequencing and NanoString expression data provide a of 213 transcriptional regulators that are upregulated and or downregulated and after 1, 3, 5, 7, 9, or 12 hr of a node graft. For these candidate GRN CTCF-ChIPseq data was used to predict CTCF-bound of up to 500 Within these, sites that are enriched for H3K27ac or H3K27me3 were identified by ChIPseq performed on induced and uninduced tissues at 5, 9, or 12 hr. Chromatin sites were categorized to Figure 2—figure supplement and those belonging to Indices 1, 3, and 7 were selected to build a These were for transcription factor binding that to other GRN This identifies genes A and as potential regulators of target and genes and as potential regulators of target At each time point, the of interactions was predicted according to the regulatory rules in interactions are when a putative input and target are interactions are modelled when an input and target have expression interactions are modelled in BioTapestry using for positive interactions and for interactions. are shown in at time points when they are expressed or when they are A GRN comprising 5614 predicted interactions between 175 components is using is by the time point when it first starts to express at hr 1 hr 3 hr 5 hr 7 hr 9 hr and 12 hr Figure data 1 of GRN transcription factor binding sites within regulatory sites at and following a node graft. Download The GRN (Figure interactions between 175 was represented using BioTapestry 1 and 2). are predicted to positive and/or negative which are and can time points. For example, the non-neural gene is predicted to be repressed by after just 1 hr of a node graft, by repression and other factors over the remaining time course. the other is predicted to be induced by after 1 hr, and 2). In 5614 interactions are predicted to occur between these the and dynamic sequence of regulatory events that are by signals from the organizer. regulatory sites and their dynamics into the network BioTapestry networks represent each gene with regulatory and each regulator as a single input into the target it is known that gene expression is controlled by each with spatial and temporal In the elements were identified within non-coding regions that are across by their to appropriate expression in This approach was for example, to identify many enhancers SOX2 expression et al., 2003; et al., – the of expression of most genes is to be controlled by elements. Our epigenetic analysis that in all there are putative regulatory sites associated with genes that change which can be in different chromatin (Figure To illustrate the of different chromatin we a for which is in a relatively gene with only a few putative regulatory elements associated with The changes that occur as sites undergo epigenetic ‘activation’ or ‘repression’ Indices Figure were modelled across the 12 hr time course. BioTapestry a for regulatory elements of the candidate regions of in Figure at 5, 9, and 12 hr represented as (Figure This which elements contain binding sites for each transcription One in Figure is initially it does not to inputs from other GRN This is then methylated at 9 and 12 hr as a number of other elements become to and In situ (Figure confirms that is not expressed in the pre-primitive-streak stage epiblast or in the neural plate at expression enriched later in the neural plate and neural as predicted by the network regulatory This the complexity of regulatory dynamics that gene expression at a relatively locus. it is not to represent the network including all the elements in their different in this way using all identified elements and their predicted activity are presented in the associated UCSC together with the predicted binding sites to the network within each a is also given in Figure data Figure 5 Download asset Open asset A regulatory sites and their (A) UCSC view of RNA-sequencing and ChIPseq associated with is upregulated after 9 and 12 hr of neural induction The regulatory is defined by flanking CTCF-bound putative regulatory sites within this were predicted based on the ChIPseq H3K27 This sites that undergo activation (Indices 7, in repression (Indices in or show no change in factor binding sites by network components are for components that are expressed at the same time point and for those that are A BioTapestry was from these predicted binding sites and expression 3, active at 12 hr, is not shown in the as there is no that is expressed at 12 hr and predicted to to during neural is initially not 1 is active but is not predicted to be by other gene regulatory network (GRN) GRN is according to the time point when it first starts to express at hr 1 hr 3 hr 5 hr 7 hr 9 hr and 12 hr is upregulated after 9 sites 2 and undergo activation and could be by transcription factors that are also 1 expression is at 12 regulators to sites 5, and is no predicted to be In situ transcripts in the neural plate at and neural tube at but not at and of core transcriptional regulators of the network Genes with and are often considered to be in a The of a network A is the number of interactions from A to other components in the network A whereas the of a network A is a that the of A based on the at which interactions between of other genes do only by A and To predict
Source: eLife
Language: English
Subfield: Molecular Biology
Domain: Life Sciences
Open Access status: gold
Estimated APC paid: $2,500
OpenAlex ID: W4360923093