Histone concentration regulates the cell cycle and transcription in early development

Upon fertilization the embryo must produce enough cells to pattern a functioning organism. Many species accomplish this by undergoing rapid cleavage divisions with little transcription. This is followed by zygotic genome activation (ZGA) (Tadros and Lipshitz, 2009; Vastenhouw et al., 2019; Harrison and Eisen, 2015; Blythe and Wieschaus, 2015a). In many organisms, ZGA coincides with the cell cycle slowing that precedes the onset of morphogenesis during a period called the mid-blastula transition (MBT) (Newport and Kirschner, 1982a,b). In Drosophila, the initial divisions are further shortened through omission of cytokinesis, resulting in a syncytium for the first 13 nuclear cycles (NCs), which becomes cellularized when the cycle slows in NC14 (Foe and Alberts, 1983). This switch from transcriptionally silenced, rapid NCs to transcriptionally active, slower NCs is accompanied by changes in the chromatin landscape. In Drosophila, transcription factors bind their consensus sequences (Harrison et al., 2010, 2011; Liang et al., 2008; Driever and Nüsslein-Volhard, 1989), RNA polymerase is gradually recruited to sites of transcription (Gaertner et al., 2012; Chen et al., 2013; Blythe and Wieschaus, 2015b) and nucleosome-free regions open on promoters (Li et al., 2014; Blythe and Wieschaus, 2016). At a larger scale, the genome becomes folded into topologically associated domains, and heterochromatin is differentiated from euchromatin for the first time (Hug et al., 2017; Shermoen et al., 2010; Seller and O'Farrell, 2018).

The early divisions are fueled by maternally supplied RNAs, proteins and metabolites (Newport and Kirschner, 1982a,b; Vastag et al., 2011; Song et al., 2017; Collart et al., 2013; Djabrayan et al., 2019; Liu et al., 2019). Chromatin components, including core histones, are loaded in vast excess of what is required for the first divisions (Horard and Loppin, 2015; Adamson and Woodland, 1974; Woodland and Adamson, 1977; Shindo and Amodeo, 2019). In Drosophila, maternal histone stores are sufficient for progression through the MBT, but zygotic histone production is required for progression past the first post-MBT cell cycle (Günesdogan et al., 2010, 2014; Zhang et al., 2018). It has been suggested that overabundant histones compete with transcription factors for DNA binding to repress transcription in the early cycles (Almouzni et al., 1990, 1991; Prioleau et al., 1994; Almouzni and Wolffe, 1995; Amodeo et al., 2015; Joseph et al., 2017). Indeed, pioneer transcription factors are required to evict nucleosomes at ZGA, and altering the concentration of core histones in the early embryos of both Xenopus and zebrafish can modulate the timing of ZGA (Liang et al., 2008; Harrison et al., 2010, 2011; Lee et al., 2013; Leichsenring et al., 2013; Amodeo et al., 2015; Joseph et al., 2017). Although the loss of embryonic linker histone H1 has been implicated in timing Drosophila ZGA, the role of core histones has not been well characterized (Pérez-Montero et al., 2013).

Here, we examine the effect of histone concentration on the MBT in Drosophila. We manipulate maternally deposited histone concentration and examine the effects on cell cycle and transcription. Histone reduction results in slower and fewer nuclear cycles, whereas histone increase results in additional cycles. We demonstrate that cell cycle changes correspond to large changes in zygotic transcription and maternal mRNA degradation consistent with alterations to ZGA. We identify a subset of transcripts as directly histone sensitive. The promoters of histone-sensitive genes are enriched for the pioneer transcription factor Zelda (ZLD) and the absence of the transcriptional repressor SIN3A and associated histone deacetylases (HDACs). Together, these results demonstrate a direct role for histone concentrations in regulating the MBT, provide a list of transcripts that are both directly and indirectly sensitive to histone concentration, and identify the chromatin features that underlie direct histone sensitivity.

To understand the effects of histone concentration on the MBT we reduced maternally supplied histones by downregulating the gene encoding a crucial histone regulator, Stem-Loop Binding Protein (Slbp) (Sullivan et al., 2001; Dominski et al., 2002; Lanzotti et al., 2002; Iampietro et al., 2014; Lefebvre et al., 2017; He et al., 2018) via maternally driven RNAi (Perkins et al., 2015). Under these conditions, histone H2B was reduced by ∼50% and H3 by ∼60% at the MBT (Fig. S1A). Approximately 50% of embryos laid by Slbp RNAi mothers (henceforth Slbp embryos) that form a successful blastoderm do not undergo the final division and attempt gastrulation in NC13 (Fig. 1A, Movie 1). Another ∼30% exhibit an intermediate phenotype of partial arrest, with only part of the embryo entering NC14 (Movie 2). A minority of Slbp embryos begin gastrulation with all nuclei in NC14. NC12 duration was predictive of NC13 arrest, with NC12 being an average of ∼5 min longer in Slbp embryos that went on to arrest compared with those that did not arrest (P=0.0034035) (Fig. 1C, Table S1, see Materials and Methods).

We first detected cellularization in wild-type (WT) embryos ∼20 min into NC14. Partially arrested Slbp embryos also began cellularization ∼20 min into NC14, with nuclei that arrested in NC13 waiting until the remainder of the embryo had entered NC14 to cellularize. Fully arrested embryos began cellularization ∼20 min into NC13, initiating cellularization one cycle early and ∼20 min earlier in overall developmental time than WT. Despite their reduced cell number, these embryos form mitotic domains and gastrulate without obvious defects, however they die before hatching (Fig. S2A,B).

To examine the effects of increased histone concentration on developmental timing we monitored cell cycle progression in embryos from abnormal oocyte (abo) mutant mothers (henceforth abo embryos). abo is a histone locus-specific transcription factor, the knockdown of which increases the production of replication-coupled histones, particularly H2A and H2B (Berloco et al., 2001). We found that abo increased H2B by ∼90%, whereas total (combined replication-coupled and replication-independent) H3 was not affected in NC14 embryos (Fig. S1B). Approximately 60% of abo embryos displayed fertilization defects or catastrophic early nuclear divisions (Table S1, Tomkiel et al., 1995). Of abo embryos that formed a functioning blastoderm, ∼6% underwent a complete extra nuclear division before gastrulating in NC15, whereas ∼4% underwent a partial extra nuclear division (Fig. 1B, Movies 3,4). Embryos from abo mothers that completed total extra divisions had faster NC14s in which they did not cellularize and spent 40-60 min in NC15 before gastrulating (Fig. 1D). This suggests an alteration of the normal transcription-dependent developmental program. In some cases, the cell cycle program and transcriptional program may be decoupled, evidenced by the fact that some abo embryos attempted to gastrulate while still in the process of division. abo embryos that underwent extra divisions exhibited a range of gastrulation defects including expanded mitotic domains and ectopic furrow formation (Fig. S2A,C, Movies 3,4).

As alterations in histone levels can both decrease and increase the number of divisions before cell cycle slowing, we reasoned that histone levels might affect activation of checkpoint kinase 1 (Chk1, also known as grp), which is required for cell cycle slowing at the MBT (Fogarty et al., 1994, 1997; Sibon et al., 1997, 1999). To test this, we crossed a fluorescent biosensor of Chk1 activity into the Slbp background (Deneke et al., 2016). We found that even in Slbp embryos that did not undergo early gastrulation, Chk1 activity was higher than in WT, consistent with the lengthened cell cycle (Fig. 1E). This result indicates that the observed cell cycle phenotypes in the histone-manipulated embryos are likely mediated through changes in Chk1 activity.

As cellularization and gastrulation require zygotic transcription, we suspected that embryos with altered development likely have altered gene expression. We performed single-embryo RNA-seq on staged Slbp embryos that remained in NC13 for more than 30 min (Fig. 2A). We compared these with either NC-matched (NC13) or time-matched (NC14) WT embryos. To control for maternal effects of Slbp RNAi, we collected pre-blastoderm stage WT and Slbp embryos. We found that the Slbp embryos undergo ZGA one NC earlier than WT. We identified ∼5000 genes that are differentially expressed (see Materials and Methods for details) between Slbp and WT NC13, with ∼60% being upregulated (Fig. 2B, Tables S2-S5,S19,S21). The upregulated genes have largely previously been identified as new zygotic transcripts, including cell cycle regulators such as fruhstart (frs, also known as Z600) and signaling molecules such as four-jointed (fj), whereas the downregulated genes are enriched for maternally degraded transcripts (Fig. 2C,D, Table S6,S16,S22) (Lott et al., 2011; De Renzis et al., 2007). We believe this represents a coherent change in ZGA timing instead of global transcription dysregulation, as 98% of the genes that are overexpressed in Slbp are expressed before the end of NC14 in our control or previously published datasets (Table S20) (Lott et al., 2011). Indeed, the transcriptomes of histone-depleted embryos that stopped in NC13 are more similar to WT NC14 than WT NC13, which suggests a role for cell cycle elongation in ZGA (Fig. S3A, Yuan et al., 2016). Nonetheless, ∼1500 genes are differentially expressed between Slbp NC13 and WT NC14 without accounting for differences in ploidy. Of these, the majority of the ∼1000 overexpressed genes are again associated with zygotic transcription, and downregulated genes associated with maternal products (Fig. S3E, Tables S2-S6,S14,S15). Thus, ZGA is even further accelerated in the histone knockdown than can be explained by purely time alone.

As ZGA is accelerated by histone depletion, we asked whether ZGA would be delayed in the histone overexpression mutant. We performed RNA-seq on pools of abo and WT embryos collected 15-30 min into NC14 (Fig. 3A). We identified >1000 genes that are differentially expressed between abo and WT, with approximately equal numbers of genes up- and downregulated (Fig. 3B, Fig. S4, Tables S2,S7-S9,S16,S19,S21). As expected, the downregulated genes in abo are enriched for previously identified zygotically expressed transcripts (Lott et al., 2011; De Renzis et al., 2007) (Fig. 3D, Tables S6,S14,S16,S22), and upregulated transcripts are enriched for maternally deposited genes. Thus, histone overexpression delays the onset of ZGA.

Zygotic genes, the transcription of which is upregulated by histone depletion and downregulated by histone overexpression, contain many important developmental and cell cycle regulators including: frs, hairy (h), fushi tarazu (ftz) and odd-skipped (odd) (Figs 2C,3C, Tables S4,S7,S8). Conversely, the maternally degraded transcripts that are destabilized by histone depletion and stabilized by histone overexpression include several cell cycle regulators such as Cyclin B (CycB), string (Stg, also known as Cdc25^string) and Myt1 (Tables S4,S7,S8). Therefore, histone concentration can modulate the expression and stability of specific cell cycle regulators, which may contribute to the onset of MBT.

As histone concentration has previously been implicated in sensing the nuclear-cytoplasmic (N/C) ratio (Amodeo et al., 2015), we compared the genes that are changed in both the histone under- and overexpression embryos with those that had previously been found to be dependent on either the N/C ratio or developmental time (Lu et al., 2009). Both previously identified N/C ratio-dependent and time-dependent genes (Lu et al., 2009) followed the same general trends as the total zygotic genesets (De Renzis et al., 2007), indicating that histone availability cannot explain these previous classifications (Fig. S6, Table S6,S22).

Next, we sought to disentangle the effects of cell cycle length from transcription in the histone overexpression mutant. We performed single-embryo time-course RNA-seq on abo and WT embryos collected 3 min before mitosis of NC10-NC13 and 3 min into NC14 (Fig. 3E). In addition, we collected unfertilized embryos (henceforth NC0) of both genotypes to control for differences in maternal contribution. Even with a stringent selection process that accounted for cell cycle time and embryo health (see Materials and Methods, Fig. S5A,B, Tables S2,S10-S13), we identified a small set of robustly upregulated (179) and downregulated (260) genes across NC10-NC14 (Tables S16,S19,S21). Of the newly transcribed genes, we detected 111 genes with delayed transcription, including frs (Fig. 3F) and only 37 that are upregulated (Table S14). We confirmed these results using qPCR (Fig. 3G). When compared with previous datasets, zygotic genes tend to be underexpressed, as was the case for the pooled abo dataset; however, the majority of these enrichments are not statistically significant (Fig. S6D). Nonetheless the majority of these underexpressed genes are expressed during NC14 in WT (Lott et al., 2011; Fig. S7). This geneset, in combination with the time-matched Slbp comparison, enables further examination of the chromatin features that underlie histone sensitivity for transcription independent of cell cycle changes.

To identify chromatin features associated with histone sensitivity, we compared the presence of 143 modENCODE chromatin signals near the transcriptional start site (TSS±500 bp) of genes whose expression is altered by changes in histone concentration independent of cell cycle time (Fig. S3, Figs 3E,4A). We found a clear pattern of unique chromatin features for the histone-sensitive genes, compared with all newly transcribed genes, that is highly similar between the histone over- and underexpression experiments (Fig. 4B, Tables S17,S18). The pioneer transcription factor ZLD, known to be important for nucleosome eviction during ZGA, is enriched in the promoters of histone-sensitive genes. Insulator proteins such as BEAF-32 and CP190 are depleted in histone-sensitive genes (Chen et al., 2013; Li et al., 2014). Promoters of histone-sensitive genes also show a strong reduction for SIN3A, a transcriptional repressor associated with cell cycle regulation (Pile et al., 2002; David et al., 2008). SIN3A is known to recruit HDACs to TSSs, and almost all HDACs also show significant de-enrichment at the TSSs of histone-sensitive genes (Kadosh and Struhl, 1998; Silverstein and Ekwall, 2005). Taken together, these marks make up a unique chromatin signature that may sensitize a locus to changes in histone concentration, as is likely for pioneer factors such as ZLD. Other aspects of this signature may indicate that these genes are subsequently subject to later developmental regulation, as indicated by H3K4me3 and H3K27me3 (Li et al., 2014; Chen et al., 2013).

In this study, we have demonstrated that histone concentration regulates the timing of the MBT in Drosophila, resulting in both early gastrulation and extra pre-MBT divisions from histone reduction and increase, respectively. Histone concentration also regulates ZGA. Thousands of genes are prematurely transcribed in histone-depleted embryos and hundreds of genes are delayed in histone-overexpressing embryos. The majority of these genes appear to be downstream of changes in cell cycle duration, suggesting a model in which histones directly regulate cell cycle progression. In other cell types, histone loss halts the cell cycle via accumulation of DNA damage and stalled replication forks (Ye et al., 2003; Prado and Aguilera, 2005; Groth et al., 2007; Gunjan and Verreault, 2003). In the early embryo, changes in histone availability may similarly create replication stress to directly or indirectly activate Chk1 as we have shown. In turn, Chk1 would inhibit Stg and/or Twine to slow the cell cycle (Di Talia et al., 2013; Farrell and O'Farrell, 2013; Deneke et al., 2016; Royou et al., 2008; Fasulo et al., 2012; Price et al., 2000; Ji et al., 2004; Stumpff et al., 2004; Shimuta et al., 2002). This mechanism is supported by previous observations that loss of zygotic histones causes the downregulation of Stg in the first post-MBT cell cycle (Günesdogan et al., 2014). In this case, the observed transcriptional changes would be independent or downstream of the altered cell cycle.

Alternatively, direct changes in transcription downstream of histone availability may feed into the cell cycle. In bulk, histone-sensitive transcripts might underlie the replication stress that has been previously proposed to slow the cell cycle at the MBT (Blythe and Wieschaus, 2015b). Consistent with this, the cell cycle lengthening and partial arrest phenotypes observed in mutant RNA Pol II embryos occur at a similar frequency to those we observe as the result of histone depletion (Sung et al., 2013). Another possibility is that specific histone-sensitive transcripts are responsible for cell cycle elongation. One promising candidate for a histone-sensitive cell cycle regulator is the N/C ratio-sensitive CDK inhibitor frs, as zygotic transcription of frs plays a crucial role in stopping the cell cycle at the MBT (Großhans et al., 2003). In contrast, tribbles, an N/C ratio-dependent inhibitor of Twine that has also been implicated in cell cycle slowing, does not show a consistent response between histone perturbations (Farrell and O'Farrell, 2013). In this previously proposed model, maternal histone stores may compete with pioneer transcription factors to set the timing of transcription initiation (Amodeo et al., 2015; Joseph et al., 2017). Indeed, the central Drosophila pioneer transcription factor ZLD is enriched at the promoters of histone-sensitive genes. Moreover, we have identified a broader set of chromatin features that may sensitize individual loci to changes in histone concentrations. These include less obvious candidates for global early transcriptional regulators, such as SIN3A, HDACs and class I insulator proteins that may protect transcripts from changes in histone concentrations. Our work highlights the importance of histone concentration in regulating the timing of MBT and provides evidence that promoters of histone-sensitive genes possess a unique chromatin signature. However, future studies will be required to isolate the specific downstream effectors that respond to changes in histone concentrations in the early embryo.

Slbp RNAi (BSC: 51171), abo1 (BSC: 2525), OvoD/bTub85 (BSC: 2149) and nls-RFP (BSC: 31418) lines were obtained from the Bloomington Drosophila Stock Center. Oregon-R (Ore-R) and P(mat-tub-Gal4)mat67; P(mat-tub-Gal4)mat15 were a gift from Eric Wieschaus (Princeton University, USA). The Chk1 localization sensor was a gift from Stefano Di Talia (Duke University, USA). All fly stocks were maintained through standard methods at 25°C unless otherwise specified and grown on a standard cornmeal media. Slbp embryos were produced by crossing nls-RFP; P(mat-tub-Gal4)mat67; P(mat-tub-Gal4)mat15 driver virgin females to UAS-Slbp RNAi males at 18°C. The resulting female progeny were placed into cups with Ore-R males at 18°C and their embryos used for experiments. WT for Slbp experiments was nls-RFP; P(mat-tub-Gal4)mat67; P(mat-tub-Gal4)mat15 driver virgin females crossed to Ore-R males at 18°C. The resulting female progeny were placed into cups with Ore-R males at 18°C, and their embryos used as WT controls for Slbp experiments. abo embryos were produced by collecting nls-RFP; abo¹ homozygous females from an nls-RFP; abo¹/SM1 stock line. Females were placed into cups with Ore-R males at 25°C, and their embryos were used for experiments. For abo control experiments nls-RFP females were collected and placed into cups with Ore-R males at 25°C, and their embryos were used as WT. We found no significant deviations in cell cycle duration between WT embryos laid at 18°C or 25°C once imaging began at 22°C, indicating that any temperature-dependent effects on the cell cycle were mitigated by imaging at a constant temperature. Unfertilized embryos were collected from crossing nls-RFP; +; and nls-RFP; abo1; homozygous virgin females to st1 βTub85DD ss1 es/TM3, Sb1 males derived from the OvoD/bTub85 (BSC: 2149). Slbp Chk1 sensor embryos were produced by crossing nls-RFP;Pmat-tub-Gal4)mat67; P(mat-tub-Gal4)mat15 driver virgin females to y,w; UAS-Slbp RNAi; Cdc25C[183-251]-EGFP males at 18°C. The resulting female progeny were placed into cups with Ore-R males at 18°C and their embryos used for experiments. WT for Chk1 sensor experiments was nls-RFP; P(mat-tub-Gal4)mat67; P(mat-tub-Gal4)mat15 driver virgin females crossed to y,w ; +; Cdc25C[183-251]-EGFP males at 18°C. The resulting female progeny were placed into cups with y,w males at 18°C, and their embryos used as WT controls for Slbp Chk1 experiments.

Embryos were dechorionated using 4% sodium hypochlorite and mounted on a 35 mm coverslip dish (MatTek) and covered with water. Cell cycle observations and RNA collections for Slbp and abo embryos were taken at 22°C using a Nikon Ti-E spinning disk confocal microscope with a 20×1.3 NA air objective at 60 s/frame acquisition (1022×1500 pixel area). WT and Slbp embryos expressing a Chk1 activity sensor were acquired at 24°C using a Nikon A1R-Si HD confocal microscope with a 60×1.4 NA oil objective at 20 s/frame (500×248 pixel area). Images were processed with Nikon NIS-Elements and FIJI (Schindelin et al., 2012).

The duration of an NC was calculated from the number of frames between nuclear envelope breakdown in at least 50% of the nuclei in the embryo to 50% nuclear envelope breakdown in the next NC. Gastrulation was determined to be when multiple cells began movement away from the single tissue sheet.

Statistical significance for nuclear cycle duration between WT and Slbp embryos for each NC was determined using two-way ANOVA performed using R (3.4). Statistical parameters such as sample numbers, mean and adjusted P-value for multiple comparisons are included as follows: For NC12 – Slbp early (n=7, mean=23.9 min) versus WT (n=6, mean=16.7 min), P=0.0000583; Slbp early (n=7, mean=23.9 min) versus Slbp normal (n=4, mean=19.5 min), P=0.0034035; Slbp normal (n=4, mean=19.5 min) versus WT (n=6, mean=16.7 min), P=0.1308121. For NC13 – Slbp early (n=4, mean=86.8 min) versus WT (n=7, mean=25.6 min), P<0.0000001; Slbp normal (n=4, mean=33.3 min) versus WT (n=7, mean=25.6 min), P=0.0056589; Slbp early (n=4, mean=86.8 min) versus Slbp normal (n=4, mean=33.3 min), P<0.0000001.

WT and Slbp embryos were imaged from NC11 to NC13 mitosis with a Chk1 activity sensor (Deneke et al., 2016). Chk1 activity was determined from the cytoplasmic-to-nuclear intensity ratio of the Chk1 localization sensor in a 150×75 pixel area. The nuclear signal was taken from segmented nuclei eroded by 2 pixels to ensure nuclear signal was analyzed. The cytoplasmic signal was taken from inverted nuclear masks dilated by 2 pixels. Four embryos per genotype were analyzed.

cDNA from RNA isolated from single embryos collected 3 min into NC14 was made with random primer mix using ProtoScript First Strand cDNA Synthesis Kit following the manufacturer's protocol (New England Biolabs, E6560L). qPCR was conducted on an ABI ViiA7 using the FG Taqman GEX master mix (Thermo Fisher Scientific, 4369016) and the following gene expression assays: frs (DM01822845; VIC) normalized to RpL32 (DM02151887; FAM).

For western blotting, protein extracts were collected from WT, Slbp and abo embryos 55 min after pole bud formation, corresponding to early NC14 in the WT and individually stage-confirmed by halocarbon oil (Sigma-Aldrich, H8773). Embryos were washed with deionized (DI) water, dechorionated with 4% sodium hypochlorite, washed again with DI water, then lysed in ice-cold embryo lysis buffer (as per Günesdogan et al., 2014). Laemmli buffer (Bio-Rad Laboratories, 1610737) was added to each sample, and they were incubated at 95°C for 10 min. Five embryos were collected per sample. Proteins were separated on a TGX 12% acrylamide gel (Bio-Rad Laboratories), stain-free dye was crosslinked under UV for 1 min, and transferred to a low fluorescence PVDF membrane. Membranes were incubated overnight in rabbit anti-H3 antibody (1:2000; Abcam, ab1791) and mouse anti-H2B antibody (1:2000, Abcam, ab52484), washed, and then incubated for a minimum of 1 h in Alexa Fluor 647-conjugated goat anti-rabbit IgG antibody (1:2500; Invitrogen, A-21244) and/or Alexa Fluor 488-conjugated goat anti-mouse IgG antibody (1:2500; Invitrogen, A-11001). Fluorescence was detected using a gel imager (Bio-Rad ChemiDoc MP) and quantified in Image Lab. H2B and H3 signal were normalized to total protein using the bright ∼45 kDa band in the stain-free channel which corresponds to vitellogenin. For the Slbp and WT comparison, the normalized H2B and H3 protein concentrations were averaged for each genotype and the average WT concentration for each protein was scaled to a value of 1. Error was calculated using a linear model in R (lm in base R) to account for gel differences and the mean genotype effects and the associated standard error were extracted. (Fig. S1A,B).

Input RNA for RNA-seq and qPCR were collected from individual and pooled embryos laid by tightly staged WT (as defined above), abo homozygous and Slbp RNAi mothers gathered from apple juice agar plates with yeast. Individual embryos were placed into an RNAse-free tube, were lysed and 100 μl lysis buffer (Applied Biosystems, KIT0214) was added. Embryos were then flash-frozen in liquid nitrogen and stored at −80°C. Pooled embryos were placed in 100 μl RNAlater (Invitrogen, AM7020) and stored at 4°C. When enough embryos were collected, RNAlater was removed and embryos were processed as above. See supplementary Materials and Methods for further details on RNA collection.

The integrity of total RNA samples was assessed on a Bioanalyzer 2100 using RNA 6000 Pico chips (Agilent Technologies).

For single embryo WT and Slbp RNA samples, additional ribosomal RNA depletion was applied before the Smart-seq2 library preparation using the Ribo-Zero rRNA Removal (Human, Mouse, Rat) Kit (Illumina, MRZH11124). The cDNA samples, RNA-seq and libraries were examined on the Bioanalyzer (Agilent Technologies), using DNA HS chips for size distribution, and quantified by Qubit fluorometer (Invitrogen). Different DNA barcodes were added to each sequencing library, and the libraries to be sequenced together were pooled at equal molar amount. The RNA-seq libraries were sequenced on Illumina HiSeq 2500 Rapid flow-cells as single-end 75nt reads, following the standard Illumina protocol. Raw sequencing reads were filtered using Illumina HiSeq Control Software and only the Pass-Filter reads were used for further analysis.

For RNA samples from single time-course WT and abo embryos, the poly-A-containing RNA transcripts were converted to cDNA using an oligo-dT adaptor and further amplified by PCR following the Smart-seq2 method (Picelli et al., 2014). Illumina sequencing libraries were made from the amplified cDNA samples using the Nextera DNA library prep kit (Illumina, FC-121-1031).

For RNA samples from pooled WT and abo embryos, poly-A-containing RNA transcripts were enriched using oligo-dT bead and further converted to cDNA and Illumina sequencing library using PrepX RNA-seq library kit on the automated Apollo 324 NGS Library Prep System (Wafergen Biosystems) following the manufacturer's protocol.

Analyses on RNA-seq and modENCODE data were performed on the high-performance computing cluster (64-bit Springdale Linux) at the Lewis-Sigler Institute for Integrative Genomics, Princeton University, using the appropriate packages within the conda environment and package management system. All statistical analyses were performed using R (3.4), Bioconductor packages (3.8; Huber et al., 2015), and the conda package management system (4.5.11). In addition, we used custom Unix, Perl, Awk and Sed scripts as necessary. See supplementary Materials and Methods for greater detail on these analyses.

We are grateful to Gary Laevsky and the Molecular Biology Confocal Imaging Facility; Wei Wang, Lance Parsons, Robert Leach and the Lewis-Sigler Institute Genomics Core Facility; Gordon Gray and the Drosophila Media Core Facility at Princeton University; and Michael Denieu and Christopher H. Chandler for technical support. We thank Nareg Djabrayan, Eric Wieschaus, Martin Wühr, Stas Shvartsman and Jan Skotheim for discussion and critical reading of the manuscript. We are grateful to Shelby Blythe for helpful discussion and technical assistance. Stocks obtained from the Bloomington Drosophila Stock Center (National Institutes of Health, P40OD018537) were used in this study.