article: Positional information in developing embryos is specified by spatial gradients of transcriptional regulators. One of the classic systems for studying this is the activation of the hunchback (hb) gene in early fruit fly (Drosophila) segmentation by the maternally-derived gradient of the Bicoid (Bcd) protein. Gene regulation is subject to intrinsic noise which can produce variable expression. This variability must be constrained in the highly reproducible and coordinated events of development. We identify means by which noise is controlled during gene expression by characterizing the dependence of hb mRNA and protein output noise on hb promoter structure and transcriptional dynamics. We use a stochastic model of the hb promoter in which the number and strength of Bcd and Hb (self-regulatory) binding sites can be varied. Model parameters are fit to data from WT embryos, the self-regulation mutant hb14F, and lacZ reporter constructs using different portions of the hb promoter. We have corroborated model noise predictions experimentally. The results indicate that WT (self-regulatory) Hb output noise is predominantly dependent on the transcription and translation dynamics of its own expression, rather than on Bcd fluctuations. The constructs and mutant, which lack self-regulation, indicate that the multiple Bcd binding sites in the hb promoter (and their strengths) also play a role in buffering noise. The model is robust to the variation in Bcd binding site number across a number of fly species. This study identifies particular ways in which promoter structure and regulatory dynamics reduce hb output noise. Insofar as many of these are common features of genes (e. g. multiple regulatory sites, cooperativity, self-feedback), the current results contribute to the general understanding of the reproducibility and determinacy of spatial patterning in early development.
  One of the fundamental questions in biology is how embryos develop reproducibly, and it has many aspects. Here, we focus on the reproducibility of the spatial gene expression patterns which determine the body plan. At a broad level, one can ask what the degree of variability is in a population of embryos - the degree to which parameters controlling developmental patterning can vary before major disruptions occur. In recent years, a number of studies have made quantitative comparisons of developmental patterns between embryos in the fruit fly, Drosophila melanogaster, aided by its very well characterized molecular biology. For instance, there has been a recent focus on spatial precision of the maternally-derived Bicoid protein (Bcd; Figure 1, green), which forms an anterior-posterior (AP) concentration gradient [1]–[5]. Bcd is a transcriptional regulator of downstream segmentation genes, and has been studied as a classic example of a positional information gradient, in which alterations in the gradient shift downstream patterns in a concentration-dependent manner [6]–[10]. It has been shown, though, that Bcd has lower spatial precision than its downstream targets: the gap gene hunchback (hb; Figure 1, blue) has a mid-embryo domain boundary at a position some 2 to 7 times less variable than the corresponding Bcd concentration threshold [1], [4], [5], [11], [12]; and the pair-rule gene even-skipped starts out with Bcd-like precision but achieves hb-like precision as its pattern develops [3]. In addition, microfluidic temperature experiments have shown robust downstream patterning following extreme disruption of the Bcd gradient [13], [14]; and even with experimentally flattened Bcd gradients, embryos form gap gene patterns in the correct order [15]. All of this suggests that the initial maternal positional information is modified during development in order for expression patterns to achieve necessary levels of precision. In addition to the precision between embryos, spatial patterns within individual embryos are well-determined, with low cell-to-cell variability (or nucleus-to-nucleus, for the precellular Drosophila blastoderm), despite numerous sources of noise. These include: the state of the DNA; mRNA and protein production; intra- and inter-cellular compartmentalization; and cellular movements and ordering. Error control is likely to occur at each these to limit noise and allow development to proceed. Previous studies have investigated aspects of this, for instance, at the DNA level [16] and overall tissue level [17], [18]). In this paper, we focus on the noise which can arise in mRNA and protein production, due to the inherently random nature of reactions at low copy number. This builds on a now extensive literature of gene expression noise in single celled organisms (e. g. [19]–[24]). But by studying hb patterning - the initial conversion of the smoothly decreasing Bcd gradient into a sharp and precise zygotic pattern segmenting the body (Figure 1) - we focus on spatial noise: what are the dynamics of noise generation in hb patterning (including, in contrast to single cells, randomness in transport between nuclei); and how is noise controlled within the constraints of these dynamics, producing the nucleus-to-nucleus noise levels observed for hb? We show that some degree of Hb between-embryo positional variability can arise purely from randomness in transcription and translation. But the larger issue is that gene expression has a strong potential for amplifying the microscopic randomness of low copy number into indeterminate macroscopic patterns within an embryo (i. e. with indistinct or missing boundaries). In the present work, we investigate what dynamic features and parameter ranges are necessary for hb expression to overcome this, in order to form determinate pattern. A broad distinction can be made between gene expression noise that is external, due to fluctuations in upstream regulator concentrations or global parameters (e. g. rate constants), and internal, due to the random nature of reactions (e. g. how many molecules per unit time are created or destroyed) and transport (how many molecules arrive in or leave a unit volume in a given time). Even in the absence of external sources (i. e. with fixed, non-fluctuating inputs) internal sources will cause fluctuating output. The amplification of external noise can potentially be significant in hierarchical signalling, such as in Drosophila segmentation; but data [4], [25] indicate that Bcd delivers a relatively non-noisy signal to nuclei (discussed further below), which indicates that much of the observed between-nucleus noise in Hb is generated internally, in the process of mRNA and protein production. We directly model the noise production in hb regulation using a chemical master equation approach [26]–[28]. This treats each reaction and transport event with a probability of occurrence per unit time. At the low copy number of many of the species involved in transcriptional regulation, stochastic dynamics predominate, necessitating such a solution method; dynamics generally become more deterministic for copy numbers in the hundreds and above [29]. Several of us were involved in a previous project developing a detailed model of anterior hb expression [30]. This was based on experimental mapping of the hb promoter [31], [32], and simulated regulation in a core region of the proximal promoter responsible for anterior zygotic hb expression (green arrows, Figure 2A). The model included binding/unbinding at 6 Bcd sites (red, Figure 2A) and 2 Hb self-regulatory sites (blue, Figure 2A); Hb production and diffusion; and Bcd translation (at the anterior pole) and diffusion. The model, solved at the deterministic level, successfully predicted Hb boundary position and sharpness for wild-type (WT) and bcd and hb mutants; and showed that sharpness depends on bistable dynamics due to hb self-regulation. Following validation against these macroscopic features, we are now using the model to investigate noise generation in the hb expression dynamics. Figure 2B shows the current version of the model; reactions have been added to explicitly model mRNA synthesis. By simulating regulation at this level of detail we can determine the relative noise contributions of, for instance: binding site number and strength; binding cooperativity; self-regulation; and protein diffusion. Noise is uniquely generated by each of these aspects, and determinate pattern formation depends on their associated parameters (e. g. binding, diffusion, production and decay rates) being within controlled ranges, as well as on the types of dynamics (e. g. binding cooperativity, self-regulation). For the hb promoter, there are a number of experimental tools which allow us to distinguish these contributions, including the hb14F mutant [33], whose Hb protein does not bind DNA; and a series of lacZ reporter constructs driven by fragments of the hb promoter [34]. Data from these embryos, as well as WT, place constraints on the model parameters, allowing us to deduce their relative contributions to the generation and control of hb output noise. Tkacik et al. [35] recently studied the effects of input (Bcd) noise on Hb output. This approach assumed that Hb output exclusively depends on Bcd, in a Hill-type manner. For some cases (hb14F, lacZ constructs), such complete Bcd dependence may apply. For WT, however, Hb self-regulation is significant, and greatly influences the final expression pattern [1], [30]. Very recently, Okabe-Oho et al. [36] published results on a stochastic model of Hb production, following our previous model [30]. They modelled binding of the hb promoter, but only considered 4 possible bound-states (including that the 6 Bcds or the 2 Hbs bind simultaneously as groups). Using all of the bound states (all the b' s and h' s in Figure 2B), has enabled us to build up the binding and production constants from lacZ construct data with from 1 to 6 Bcd binding sites of different strengths, and from 0 to 2 Hb sites (together with WT and hb14F data). This has revealed cooperativity and binding strength effects which could not be addressed in the 4-state system. [35] and [36] focused on input noise due to diffusion of regulators to the promoter. [36] reported that amplification of this type of noise can depend strongly on Hb diffusivity. We see a similar sensitivity to Hb diffusivity, but the reproducibility of Hb boundary sharpness between embryos suggests that diffusivity is highly constrained (i. e. not a variable parameter); and the steep boundary indicates a slow diffusivity, i. e. that it is not optimized for noise reduction. Detection of regulator movement within nuclei is beyond current techniques. But measurements of Bcd concentrations at whole nuclei resolution [4] indicate that nucleus-nucleus Bcd fluctuations should be relatively low. At this spatial resolution (our data is processed into ‘energid’ units, of nucleus plus surrounding cytoplasm), we observe that relative noise is higher in hb mRNA than in protein, and that Hb self-regulation (comparing WT with the hb14F mutant) decreases relative noise in the protein output. Neither of these effects depend on input noise: the former (the noise difference between mRNA and protein) shows the effects of translation; the latter highlights the effects of self-feedback. Our model predicted these effects from dynamic principles. The details of the hb promoter structure matter for determining expression noise. We predict that binding/unbinding noise dominates in the absence of self-regulation, and increased binding site number and strength serves to reduce noise in these cases. In WT, though, self-feedback produces a bistable mechanism: this was previously shown to be critical in boundary sharpness [30]; the present work shows how this mechanism also promotes expression into a more deterministic, low-noise regime. These results suggest how evolution may have incorporated binding sites and self-feedback mechanisms to produce output determinate enough for robust development.
  Published data from WT, the hb14F mutant and the lacZ constructs indicate probable values of the parameters in the model (Figure 2B). We describe below how the data can be used to determine the rates (k-values) sequentially, without the need for a global parameter optimization. This involves deterministic (no noise) solution of the model to match macroscopic features of the data (strength of expression, expression boundary position, boundary sharpness, timing). This parameter set is then used in stochastic solutions of the model to make predictions on noise levels and characteristics, which are corroborated against new experimental data (next section, Stochastic Results). The assumptions involved in the parameter fitting are not expected to affect the noise predictions (see Discussion). The main points on parameter fitting are given here, with further details given in Text S1. With the parameters thus determined, we modelled WT, hb14F and the lacZ constructs to predict the mRNA and protein noise arising from different aspects of transcription and translation.
  The Bcd-Hb system has received a great deal of study over several decades, making it one of the best characterized systems for understanding the mechanistic details of positional specification by gradient reading. We have used binding site information for the hb promoter in conjunction with quantitative imaging to develop and test a stochastic model of expression dynamics. This has allowed us to characterize the noise inherent in gene expression due to the low copy number of DNA (numbers of promoters and binding sites) and mRNA. Identifying how noise is controlled in spatial gene expression patterns is a fundamental problem; the dynamics at the hb promoter provide a model for how this might occur in many patterning events. Modelling the variability in the data, in addition to modelling average features, provides an extra experimental dimension for refining and validating models of gene regulation. Prior work has focused on the effects of Bcd noise on Hb [35]. However, we experimentally observe statistically significant noise reduction in the process of translation and due to Hb self-feedback, which validate model predictions: our study points to the critical aspects of intrinsic noise arising in Hb production, independent of external fluctuations. For the noise levels associated with Bcd' s measured concentration [4], we see little, if any, effect of Bcd concentration fluctuations on hb production. External noise is overshadowed by intrinsic noise arising in binding site occupancy and modulated through transcription and translation. Work in [36] indicated that promoter occupancy noise could be reduced by fast Hb diffusion. However, the steep angle of the Hb boundary suggests a slow diffusivity not optimized for noise reduction. Our work suggests several mechanisms in the hb dynamics which reduce noise. In the absence of Hb self-feedback, as in the hb14F mutant and the lacZ reporter constructs, output is noisier - less shielded from the noise of binding site occupancy. Fitting the binding strengths and production rates associated with 6 individual Bcd sites and 2 Hb sites to the lacZ data has enabled us to predict the degree to which increasing number and strength of binding sites (Bcd and Hb) can buffer promoter occupancy noise. Such basic noise reduction may have evolved fairly independently of other mechanisms; our computations suggest hb noise is similarly controlled in flies with between 4 and 10 Bcd sites. Self-feedback is a major component of WT Hb expression [1], responsible for large changes in protein production and boundary sharpness from hb14F to WT [30]. We observe a significant effect of self-feedback in reducing protein noise levels (Figures 3C vs. 4A; Tables 1 and 2). Our calculations indicate self-feedback is responsible for a change in anterior mRNA levels from roughly 30 (hb14F, Figure 4C) to roughly 200 copies per nucleus (WT, Figure 3B), from a noise-dominated to a more deterministic regime [29], which in turn reduces translational variability in WT. In addition, we predict that self-feedback reduces variability arising from transcription at multiple promoters within nuclei (Table 6). We predict that the noise reduction seen with translation (Figure 3C vs. 3D, Table 1) is similarly due to a concentration difference, since one copy of mRNA makes multiple copies of protein (Figure 3A vs. 3B; Figure 4B vs. 4C; Tables 1 and 2). Translation is expected to produce non-Poisson hb probability distributions, distinct from the Poisson noise expected for Bcd. Such non-Poisson ‘bursting’ noise has been characterized in yeast [23], in which protein variance scales with mean concentration but with a VMR far greater than for a Poisson process, due to amplification of noisy low mRNA copy number. For hb, we expect deviations from Poisson to be stronger in protein than mRNA, and strongest with self-feedback than without. The observation of lower relative noise in WT protein than either WT mRNA or hb14F protein suggests that its higher concentration overcomes the non-Poisson bursting effects. (A purely Poisson-distributed WT Hb would be expected to have even lower noise and show larger differences with both WT mRNA and hb14F protein.) Published data [4], [30], [32], [34], [38], [41], [57] places distinct constraints on any model of the hb promoter, for rates of regulator binding and production. Conclusions can be made from these regarding the relative values of these constants (e. g. relative WT and hb14F Hb levels indicate that WT production rates are 7 times higher; early vs. mid-cycle 14 data indicates Hb production needs to increase concentration by about 50% in 30 minutes), and these conclusions can be followed through for their implications for relative noise differences. If absolute concentrations were higher or lower than estimated, absolute noise levels would be lower or higher, respectively; but the relative results would hold. We have made noise predictions at energid (nucleus plus associated cytoplasm) resolution for mRNA and protein output for: WT; the hb14F mutant; 11 lacZ reporting constructs with various combinations of binding sites; and 4 other species of flies. Enhancing our confidence in the model, we have corroborated model predictions at this resolution for: higher mRNA than protein noise; higher hb14F than WT noise; and the indication of noise decrease due to the 2nd Hb site, comparing the pThb1 and pThb5 constructs. In addition to these, the model illustrates the reduction in noise available from increasing binding site number and strength; the role of moderate cooperativity and slow protein timescales for limiting noise; and the degree to which within-embryo noise can generate between-embryo variability. These results indicate the degree to which hb noise amplitudes are determined by expression dynamics, and how these dynamics produce hb probability distributions distinct from Bcd' s. We have also taken data at intra-nuclear resolution, imaging transcript production from different copies of the promoter within each nucleus. This is at current technical limits of spatial resolution (though recent very high resolution studies with GFP [58], [59] are promising for bcd). Data analysis and modelling at this level shows promise for separating transcriptional noise from other types, such as from inter-nuclear transport. (See [60] for a recent study showing the effect of promoter state on pattern synchrony at this degree of resolution.) Only Bcd and Hb regulation, at the specified binding sites, are considered in the model. The dynamics of this core promoter region reproduce many of the deterministic and stochastic features of Hb activation in the anterior region. Additional Hb sites, such as those in the distal P1 promoter, are known to affect later, posterior expression – any effect of these sites in early expression would be incorporated into the 2 Hb sites of the current model. Additional Bcd sites would not be expected to greatly influence expression, based on Musca and the simulations with up to 10 Bcd sites. The positioning of final expression patterns, especially later in cycle 14, do depend on other factors, such as inhibition by other gap genes (e. g. [61]). (Head gap genes are likely involved in the lower Hb expression at 0–15%EL, e. g. Figure 3C [15].) The region of 40–50%EL, where the Hb boundary interprets a marginal decrease in the Bcd gradient, is promising for quantifying these factors. In preliminary computations, reduced or missing initial Hb gave reduced expression in this parasegment 4 region. Similar phenotypes are observed with some gap mutations [62]. Incorporation of other regulators in the model would permit exploration of their relative contributions to Hb noise in this region. Hb itself forms a morphogenetic gradient, controlling the expression of a number of other gap genes in early segmentation, potentially at higher precision than Bcd [63]. The mechanism of secondary (Hb) gradient formation is in contrast to those for the primary (maternal) gradients, which do not form by spatially-distributed gene expression. The present work, therefore, has focused on how transcription and translation kinetics can be controlled to provide a determinate and precise secondary gradient for specification of the segmentation patterns. While our study has shown the importance of the details of promoter structure and expression dynamics (such as self-feedback) on hb' s expression noise, many of the noise motifs found here will be applicable to other genes (for instance, dependence on binding site number, non-Poisson amplification in transcription and translation). In this way, characterizing hb noise serves as a model for how zygotic gene expression gives rise to the determinate and reliable expression patterns underlying development.
  The model in Figure 2B was computed with the MesoRD software ([27]; http: //mesord. sourceforge. net). This package allows deterministic (used for parameter fitting) or stochastic (used for noise prediction) solution of mechanisms involving reaction and diffusion. A kinetic scheme is entered as elementary reactions (Figure 2B), and rate constants (values in Tables S1, S2, S3) and diffusivities (for Hb and Bcd proteins, see Tables S2, S3 for values) are specified for the model species. Geometry was specified as a one-dimensional series of 100 subvolumes (each a 5µm cube), corresponding to the energids (nucleus plus cytoplasmic neighbourhood) along the AP axis. Computations solve for model species densities in each subvolume, according to the specified reactions and between-subvolume diffusion. In deterministic simulations (for parameter searches), we used a 4th order Runge-Kutta solution method. For matching to data, boundary position and sharpness were determined at half-maximal concentration. For stochastic solutions, MesoRD solves the reaction-diffusion master equation, in which each reaction and diffusion event has a probability (set by the macroscopic rate constants) of occurring in a unit of time. The software implements the next subvolume queuing method [26], [27], [64] to significantly improve memory and processing requirements, making computation possible for the number of species and subvolumes in the hb model. All computations are run in real units (µm, s, etc.). Model parameters were determined by fitting macroscopic features of published data: boundary position and angle, expression levels, and timescales. As described in the Results, building up from the lacZ data to hb14F and WT constrains the values of the binding rates and diffusivities in the model. Further details are given in Text S1. Fly stocks, staining and imaging were as in [30]. Whole mount embryos were imaged by laser confocal scanning microscopy, from WT Oregon-R, hb14F, and the lacZ construct (pThb1,5; [34]) lines. All embryos were heat fixed and immunostained for Hb protein. Fluorescent in situ hybridization (FISH) was used for mRNA determination, for hb and for lacZ, following the method of [65]. Images were collected using an HC PL APO 206 objective and variable digital zoom (1. 2–1. 56). Fluorophores were excited by laser at different wavelengths (488,555, and 647 nm), and detected via a filterless spectral separation system. Channels were scanned sequentially. The microscope was set so that maximum expression was 255 on an 8-bit scale. To reduce photomultiplier noise, each image was scanned sequentially 16 times and the results averaged. Raw images from the confocal microscope, 1024×1024 pixels, were cropped and rotated for standardization. Each energid (each nucleus plus its cytoplasmic neighbourhood) was identified by Voronoi tessellation [66]. See Text S2 for details (Figures S1, S2). Averaged pixel intensities within each energid (Figure S3) were used for comparison to simulation output. Data was used from a 10% DV (dorsoventral) strip, centred on the AP midline, in order to minimize geometric distortion from the embryo periphery. Background fluorescence for these lateral images follows a half ellipsoid (c. f. [67]). We found the parameters of the ellipsoid for each image by a Genetic Algorithms technique (c. f. [68]). An initial visual inspection of the data for each embryo was needed to estimate the approximate height of the background. The quantitative measure of fitting quality for a given set of parameters was based on minimizing the distance between the data points and the ellipsoid surface. Previous approaches used preliminary statistical analysis of the dynamics and positioning of areas of zero specific signal (i. e. areas where all fluorescence was background), requiring large datasets [67]. Our supervised evolutionary search gives a much more efficient method for directly analyzing each image. Singular Spectrum Analysis (SSA) [69], a non-parametric technique with an adaptive filter, has been recently used for separating confocal intensity data into components [70], [71]. Its extension, 2D-SSA [72], was applied to the 2D (AP and DV) intensity surface; the leading components of the decomposition give the pattern' s trend. Noise was then quantified from the difference of each energid' s intensity to the trend value at each position (i. e. local residuals). See Text S3 (Figures S4, S5, S6, S7, S8, S9, S10, S11, and S12) for more detail on the SSA data analysis. Noise measures were calculated from the anterior expressing regions (15–45% EL) as, where data is the average pixel intensity for an energid and trend is the SSA-extracted trend at that position (i. e. the noise measure is a standard deviation for the relative residuals). Noise was calculated similarly for simulation output, from the difference of the stochastic output and the deterministic solution at each position (see Table 1 footnote).
summary: Noise is an intrinsic part of biochemical systems such as gene regulation networks. Noisy gene expression has been well documented in populations of single cells, and is likely a key mechanism in evolutionary change. But in developing embryos, cells within a tissue must overcome such variability in order to provide the uniformity required to coordinate multiple events. Reproducibility and determinacy of the spatial protein patterns preceding tissue differentiation is a critical aspect of development. In this study, we use anterior-posterior (AP) segmentation in the fruit fly (Drosophila) to understand how gene regulation dynamics control noise. One of the earliest AP patterning events is the anterior activation of the hunchback (hb) gene by the maternally-derived Bicoid (Bcd) protein gradient. This interaction has been very well characterized, providing the tools for us to develop a stochastic model of hb gene regulation to make predictions about expression noise, and to corroborate these experimentally. For hb, we show that self-regulation is a critical part of controlling noise, and the multiple Bcd binding sites in the hb promoter also enhance pattern reproducibility. To the degree that such features are shared by other genes, these noise-reducing mechanisms may be common to many pattern forming events.
keywords: computer science/applications
  computational biology/transcriptional regulation
  developmental biology/pattern formation
  biophysics/transcription and translation
  developmental biology/developmental molecular mechanisms
title: Gene Expression Noise in Spatial Patterning: hunchback Promoter Structure Affects Noise Amplitude and Distribution in Drosophila Segmentation