<title>Introduction</title> <p>How an embryo achieves pattern and form from an initially undifferentiated state has fascinated people at least since the time of Aristotle. Scientific advances on this began over a century ago, with, for example, the experiments of Hans Driesch on sea urchin embryos ##UREF##0##[1]##, from which he proposed that the embryo has a coordinate system specifying cellular position; and from the experiments of Ethel Browne ##UREF##1##[2]##, who showed that a piece of hydra mount induced a secondary axis when grafted into the body of another hydra. These and other subsequent results were synthesized by Lewis Wolpert in 1969 ##REF##4390734##[3]## into a definition of positional information. According to this concept, the spatial asymmetries of concentration gradients of chemical signals (morphogens) provide positional information during cellular differentiation; each cell (or nucleus) reads its position from the local morphogen concentration and differentiates accordingly. Wolpert's concept of morphogen gradients has become a central tenet of developmental biology ##REF##17254964##[4]##–##REF##16352710##[6]##. Modern molecular techniques have demonstrated numerous cases of protein concentration patterns in embryogenesis, and many have been shown to act as morphogens. In the late 1980's, the Bicoid (Bcd) protein gradient was characterized and its concentration-dependent effect on downstream target genes in <italic>Drosophila</italic> was demonstrated ##REF##3383245##[7]##–##REF##2567637##[9]##. This has since become one of the most studied examples of morphogen gradient signaling in developmental biology ##UREF##3##[10]##,##REF##14744427##[11]##.</p> <p>Reaction-network models have been successfully applied to describe a great variety of systems in physics, chemistry, and biology ##UREF##4##[12]##–##REF##9888852##[14]##. Along with this, many mathematical tools have been developed to support such applications. With these tools, one can show that certain reaction networks may exhibit multiple stationary states, for particular ranges of their rate constants. Bistability is a special case, in which the system can evolve to either of two asymptotically stable steady states (concentration levels). Under certain conditions, spatial patterning or oscillations can arise ##REF##7805628##[15]##–##UREF##6##[17]##. In biology, bistability has long been established in control of the cell cycle and other oscillations ##REF##9486845##[18]##,##REF##17667953##[19]##, and also recently reported in an artificial gene regulation network ##REF##17170763##[20]##. In <italic>Drosophila</italic>, spatial bistability has been proposed for dorso-ventral patterning ##REF##16864795##[21]##,##REF##15759004##[22]##.</p> <p>In early embryogenesis, the diffusion of Bcd protein, translated from mRNA localized at the anterior end of the egg, forms an exponential concentration gradient, establishing the anterior–posterior (AP) axis (##FIG##0##Figure 1A and 1C##) ##REF##3383244##[8]##,##REF##3686007##[23]##,##REF##2901954##[24]##. Bcd is a transcriptional regulator, and through its asymmetric distribution acts as a morphogen, governing the positions at which the downstream gap genes will be activated. In combination with cross-regulation between these genes, the initial Bcd asymmetry is propagated and refined, establishing the first stage of embryo segmentation ##REF##2567637##[9]##, ##UREF##7##[25]##–##REF##10910359##[32]##. It is still not well characterized, however, what mechanisms interpret the smooth Bcd positional information into sharp and precisely positioned downstream target gene expression.</p> <p>
<italic>hunchback</italic> (<italic>hb</italic>) is one of the first gap genes activated by Bcd, with strong anterior expression and a sharp on–off boundary at mid-embryo (##FIG##0##Figure 1B and 1C##) ##REF##2567637##[9]##, ##UREF##9##[33]##–##REF##2797150##[35]##. Anterior <italic>hb</italic> activation depends on Bcd, as shown by Struhl et al ##REF##2567637##[9]## and Driever et al ##REF##2911348##[34]##, and on its own self-regulation, as already reported by Treisman et al ##REF##2797150##[35]## and Margolis et al ##REF##7555732##[36]##; many Bcd and Hb binding sites have been identified in the <italic>hb</italic> promoter region, as reported by Treisman et al., among others ##REF##2797150##[35]##–##REF##8199047##[37]##. Hb has maternal (<italic>hb<sub>mat</sub></italic>) and zygotic contributions, and provides positional information for other gap genes, such as <italic>Krüppel</italic> (<italic>Kr</italic>), <italic>knirps</italic> (<italic>kni</italic>), and <italic>giant</italic> (<italic>gt</italic>), and for the homeotic gene <italic>Ultrabithorax</italic> (<italic>Ubx</italic>) ##REF##11518507##[38]##–##REF##10741965##[41]##. Removal of both maternal and zygotic <italic>hb</italic> expression results in severe deletions and polarity reversals of the most anterior segments ##REF##3803711##[42]##. In normal development, Hb expression drops from highest to lowest over about 10% egg length (##FIG##0##Figure 1B and 1C##); Considerable attention has been focused on what molecular mechanism generates this Hb sharpness. Early on, it was shown that a <italic>hb</italic> enhancer element of 300 base pairs (bp), containing 6 Bcd binding sites, is sufficient to reproduce the regulatory activity of Bcd on <italic>hb</italic>
##REF##2567637##[9]##,##REF##2911348##[34]##. It was shown that Bcd binds to these sites cooperatively and it was hypothesized that, due to this cooperativity, a small change in Bcd concentration across some threshold could produce a large change in <italic>hb</italic> promoter occupancy, generating the on–off expression pattern ##REF##2567637##[9]##, ##REF##2911348##[34]##, ##REF##2502714##[43]##–##REF##16150708##[46]##. However, these studies did not establish that cooperativity is sufficient to generate Hb border sharpness.</p> <p>To quantify the degree of Bcd's cooperativity, Ma et al. ##REF##8620846##[44]## used a six-Bcd site fragment of the <italic>hb</italic> promoter in a DNase I footprint assay, and found a Hill coefficient of about 3.6; Burz et al. ##REF##9774343##[47]##, using a gel-shift assay with a 230 bp <italic>hb</italic> enhancer, found a Hill coefficient of 3.0. From quantified <italic>in vivo</italic> patterns of Bcd and Hb proteins, Gregor et al. ##REF##17632062##[48]##, estimated a higher value for this coefficient, of around 5 (though the effects of <italic>hb</italic> self-regulation were neglected, addressed further in the <xref ref-type="sec" rid="s3">Discussion</xref>); and suggested that it could support the proposal of Crauk and Dostatni ##REF##16271865##[49]## that Hb expression is entirely determined by Bcd cooperative binding. However, systems with such high Hill coefficients would be expected to show temperature sensitivity. Houchmandzadeh et al. ##REF##11845210##[50]## showed that the Bcd gradient is strongly affected by temperature changes of 20°C, but that the Hb pattern is largely unaffected. Dependence on Bcd with Hill coefficients between 3 and 5 would be expected to show far greater effects on Hb than are observed, indicating that Bcd cannot be the only factor regulating the Hb border. The insufficiency of Bcd cooperativity to produce Hb sharpness is also supported by the findings of Simpson-Brose et al. ##REF##8087852##[51]##, who showed that synergy between Hb and Bcd is necessary to establish the expression patterns of the gap genes, including <italic>hb</italic> itself.</p> <p>To address these issues, we have taken a combined experimental and theoretical approach to understand how the <italic>hb</italic> gene converts the positional information of the smooth Bcd gradient into a sharp expression pattern. We used wild-type (WT) embryos to experimentally determine how Hb position and sharpness change in time; and we measured how these quantities are affected in embryos mutant for Bcd cooperative binding and for <italic>hb</italic> self-regulation, and by use of an artificial promoter with 6 Bcd and 2 Hb binding sites. We also developed a predictive reaction–diffusion model of <italic>hb</italic> transcription, incorporating both Bcd cooperative binding and <italic>hb</italic> self-regulation. By fitting this model to wild-type Bcd and Hb patterns, we determined kinetic parameters of the model, such as binding constants. With these parameters, our model successfully reproduces the dynamic development of the Hb pattern. By reducing Bcd binding constants or the number of Bcd binding sites, our model reproduces the same mutant phenotypes as our experiments, and predicts a loss of sharpness for a <italic>hb</italic> self-regulation mutant, which we experimentally verified. By applying dynamical systems theory to the model, we show that Hb sharpness is due to spatial bistability stemming from <italic>hb</italic> self-regulation. This means that Hb does not have a single steady-state concentration continuously dependent on Bcd, but that at a threshold Bcd concentration, two stable steady states become available to Hb, and a small change in Bcd concentration can create a dramatic shift in Hb concentration. This shift between steady states is responsible for the sharpness of the Hb boundary. The position of the Bcd threshold is controlled by Bcd cooperative binding, but this mechanism itself is not sufficient to generate <italic>hb</italic>'s expression sharpness.</p>
<title>Materials and Methods</title> <title>Dataset</title> <p>We stained for Bcd and Hb proteins in WT Oregon-R embryos, as well as in the <italic>hb</italic> mutant <italic>hb<sup>14F</sup></italic>
##REF##8001780##[52]##, and two <italic>bcd</italic> mutants ##REF##16150708##[46]## (<italic>bcd<sup>K57R</sup></italic>, <italic>bcd<sup>E1</sup>/+</italic>). lacZ expression for the pThb5 construct (driven by a fragment of the <italic>hb</italic> promoter; ##REF##2502714##[43]##) was visualized by two methods for staining β-galactosidase mRNA (##TAB##1##Table 2##). The simultaneous Hb protein and RNA visualization was also done in WT Oregon-R embryos.</p> <title>Preparation and Staining</title> <p>As outlined in ##TAB##1##Table 2##, three different staining procedures were used for obtaining expression patterns. For all procedures, embryos were dechorionated; heat fixed in NaCl 0.7%+Triton-X100 0.05% for 3 seconds and immediately chilled in ice; and devitellinized with methanol.</p> <p>For protein staining ##REF##15834586##[76]##, embryos were incubated with guinea pig and rat primary antibodies to Hb and Bcd, respectively, followed by secondary antibodies labeled with Alexa Fluor 647 (to Hb) and 488 (to Bcd; Molecular Probes). All antibody incubations and washes were done in PBS+0.1% Tween-20. Blocking was done in Western Blocking Reagent (Roche), diluted 5 times. All secondary antibodies were preabsorbed by incubating them with 0- to 12-h-old WT embryos for at least 2 h at 4 C.</p> <p>For the lacZ embryos, we used simultaneous immunostaining to Hb and Bcd and in situ hybridization. With FISH, we followed the method of Janssens et al. ##REF##15834586##[76]##: a lacZ riboprobe was prepared with a 2.5-kb PvuII lacZ fragment blunt-cloned into the EcoRV site of pBluescriptIIKS+ (gift from S. Small), labeled with fluorescein by transcription using T3 polymerase. After hybridization, lacZ mRNA was visualized by sequential incubation with rabbit antibody to fluorescein (Molecular Probes), followed by antibody to rabbit labeled with Alexa Fluor 488 (Molecular Probes). The embryos were simultaneously stained for Hb and Bcd proteins, as in the previous paragraph, using secondary antibody labeled with Alexa Fluor 555 to detected Bcd. Alternately (##SUPPL##10##Figure S6##), some lacZ embryos were measured via enzymatic staining (whole mount in situ hybridization): β-galactosidase mRNA was hybridized in situ with a digoxygenin-labeled DNA probe, following standard protocols ##REF##2476281##[67]##. The hybridization products were detected with a phosphatase-coupled antibody against digoxygenin. For simultaneous determination of Hb protein and RNA, we used the same FISH procedure as above, sequentially using guinea pig and rabbit antibodies to Hb and fluorescein, respectively, and secondary antibodies to guinea pig and rabbit labeled with Alexa 647 and 488, respectively.</p> <p>Following fixation and staining, embryos were mounted in 40 ml mounting medium (Prolong Antifade by Invitrogen) and covered with a 22×30 mm cover glass (No. 1½).</p> <title>Confocal Microscopy</title> <p>Following the methods of Janssens et al. ##REF##15834586##[76]##, whole-embryo images were taken using a laser confocal scanning microscope (Leica TCS SP2). Images were collected using an HC PL APO 20× objective and variable digital zoom (1.2–1.5×). 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. To reduce image noise from the photomultiplier tubes, each embryo was scanned sequentially 16 times and the results averaged.</p> <p>The settings of the microscope were adjusted for each gene product such that pixels expressed at maximum intensity were 255 on the 8-bit scale. Initial image size before processing was 1024×1024 pixels. Raw images were averaged, cropped and rotated. This standardization allowed us to compare levels of gene expression at different times, or in different experiments performed on different days ##REF##15059825##[77]##.</p> <title>Processing of Images</title> <p>For embryos triply-stained for segmentation proteins, the extraction of AP intensity profiles is well established ##REF##15834586##[76]##,##REF##15059825##[77]##. With such data, a nuclear mask can be created, and intensity data mapped to nuclei (next section). Co-staining for Bcd and Hb proteins and β-galactosidase mRNA presents greater challenges: signal strength and quality are very different for proteins and RNA; and the anterior localization of Bcd and Hb make identification of posterior nuclei very difficult. We developed a non-mask method for profile extraction for these experiments (section after next).</p> <title>Image Processing via Nuclear Masks</title> <p>For embryos stained for three segmentation proteins, the three images are used to generate a ‘pixel maximum’ image, of the brightest pixels among the images. On this image, pixels are then classified as belonging to a nucleus or not, by edge-detection of bright nuclei against dark background. An error-correction step repairs any ‘fused’ nuclei. With the resulting nuclear mask, dorso-ventral, AP coordinates, and average fluorescence level of the three gene products can be mapped to individual nuclei. Intensity profiles are extracted from a central 10% strip of nuclei along the AP axis ##REF##15834586##[76]##,##REF##15059825##[77]##.</p> <title>Direct Image Processing (No Nuclear Mask)</title> <p>For lacZ embryos, co-stained for Hb and Bcd proteins and β-galactosidase mRNA, nuclei cannot be reliably identified, especially in the posterior (preliminary nuclear staining in a fourth channel shows much crosstalk). For these experiments, we directly extract the pixel intensities in a 10% strip (corresponding area to above). For high-intensity protein staining the signal is strong, but for low-intensity RNA staining we must recover expression from a noisy signal (next section). A one dimensional (1D) AP profile was created from the strips, by averaging intensities in each DV pixel column from the central 10% strip along the AP axis. In addition to some between-pixel noise, the resulting profiles show noise in nuclear order and in the distribution of stained material between nuclei and cytoplasm. Minimization of these two sources of noise is described in the next section. To test the quality of our direct method, we manually made nuclear masks for several co-stained lacZ embryos using the multiple ROI feature in <italic>ImageJ</italic> software ##UREF##20##[78]##. ROIs are circles with radii comparable to the nuclear radii in a given image. Each ROI was positioned manually to outline a given nucleus. Nuclear-resolution AP profiles from this method are of comparable quality to pixel-resolution profiles from our direct extraction method.</p> <title>Signal–Noise Decomposition</title> <p>Noise in intensity profiles can influence model-fitting and statistical analysis of expression patterns ##REF##15059825##[77]##,##UREF##20##[78]##. To obtain clear expression patterns, we used singular spectrum analysis (SSA ##UREF##21##[79]##), a non-parametric technique with an adaptive filter. This allowed us to remove experimental (e.g., photomultiplier tube) noise and noise due to variability in nuclear order and in nuclear-cytoplasmic distribution of gene products. We used the methods of Golyandina et al ##UREF##22##[80]##, and software developed by Nina Golyandina and Theodore Alexandrov ##UREF##23##[81]##.</p> <title>Background Removal</title> <p>Non-specific binding of antibodies to biological material results in background fluorescence in our images. For triple-stained protein images it has been shown ##REF##15711806##[82]## that this background is a paraboloid. For every image we calculated the parameters of this paraboloid from regions of the embryo in which a particular gene is not expressed, then transformed original fluorescence at or below this background to zero. For lacZ embryos simultaneously stained for protein and mRNA it is unclear whether background has a comparable shape; in these cases, we use a simple flat background, subtracting the minimum raw intensity off all values.</p> <title>Quantifying Border Position and Sharpness</title> <p>An advantage of direct image processing is the large number of data points (around 1000) and smoothness of each profile. This makes it possible to apply standard calculus techniques to characterize the profiles: we define the Hb domain border as the inflection point, and sharpness as the first derivative at that position. With normalized intensity data (0–100% scale), this slope can be expressed as an angle of inclination (as in ##FIG##1##Figure 2##). These techniques can be applied to data, as well as HSR model results.</p> <title>Temporal Classification</title> <p>In addition to confocal scanning, all embryos were observed along the dorsal edge with Differential Interference Contrast (DIC) optics. Distances were measured from the egg surface to the invaginating membrane, and from the surface to the cortex. The ratio of membrane depth to cortex depth was used to estimate embryo age in minutes, using a published standard curve ##REF##3151484##[63]##.</p> <title>The Zero Eigenvalue Analysis</title> <p>The Zero Eigenvalue Analysis ##UREF##11##[56]##–##UREF##13##[58]## is a very efficient method, because the search for bistability is reduced to the solution of a system of equalities and inequalities (see Equation S2.22 in ##SUPPL##3##Text S2##) that are easier to find than a direct solution of the polynomial equation describing the stationary states (see Equations S1.3′–9′, S1.15′, and S1.16′). This technique readily allows one to find the set of kinetic parameters that produce bistability, and gives two steady state solutions, which can be used to easily make the bifurcation diagram, like that shown in ##FIG##3##Figure 4##. Finding bistability with direct solution of a polynomial requires solutions that are different, real and positive. This is frequently not convenient for degree higher than 2 ##UREF##24##[83]##, and not analytically solvable for degree higher than 4. Zero Eigenvalue Analysis can be applied to such higher degree systems. For example, Li ##UREF##13##[58]## has used this method to determine multiplicity of stationary states in the famous Goldbeter and Lefever allosteric model ##REF##4263005##[84]##, consisting of 14 species, 32 reactions, and 27 kinetic constants.</p>
<title>Results</title> <title>Mutant Genotypes That Affect Hb Sharpening and Positioning</title> <p>In order to investigate the relative contributions of self-regulation and Bcd cooperativity to Hb sharpness and position, we performed immunohistochemistry assays in wild-type and mutant embryos in nuclear cleavage cycle 14A (##FIG##0##Figure 1## and ##SUPPL##5##Figure S1##; see also <xref ref-type="sec" rid="s4">Materials and Methods</xref>). We used systematic image processing approaches to extract gene expression patterns from confocal microscopy images, to determine embryo ages and to quantify pattern position and sharpness.</p> <p>The effect of loss of <italic>hb</italic> self-regulation can be directly studied with the loss of function <italic>hb<sup>14F</sup></italic> allele ##REF##8001780##[52]##, which codes for a mutant protein having no DNA-binding capacity (i.e., no self-regulation). We scanned 39 embryos expressing this allele and found that embryos homozygous for the <italic>hb<sup>14F</sup></italic> allele (##FIG##0##Figure 1D##) have a strong decrease in sharpness (21.8% reduction, from 80.2° in WT to 62.8°), and small shift in position (9.3% decrease, from 47.0% EL in WT to 42.6% EL; ##TAB##0##Table 1## summarizes sharpness and position for all experiments). The <italic>hb<sup>14F</sup></italic> homozygotes were easily identified by low signal intensities ##REF##2846287##[53]## (see ##SUPPL##6##Figure S2## for non-normalized profiles). Heterozygote (<italic>hb</italic>
<sup>14F</sup>/+) and WT embryos were not easily distinguished, and are shown together in ##FIG##0##Figure 1E##. The means for this group show little or no change from WT: sharpness shows a 1.3% change (from 80.2° to 81.3°), and position shows a small shift, from 47.0% EL to 44.3% EL (5.7% change; ##FIG##0##Figure 1E##). In order to compare Hb patterns between different assays, in which absolute intensity varies, all experimental and theoretical profiles in ##FIG##0##Figure 1## are normalized. With this, our measure of sharpness is determined by the AP projection of the profile (see ##SUPPL##7##Figure S3##); this is covered in more detail in the <xref ref-type="sec" rid="s3">Discussion</xref>. Note that all <italic>hb<sup>14F</sup></italic> homozygote embryos display lower sharpness than any WT embryo older than 8 minutes in cycle 14 (the timing of sharpness maturation is discussed further in ##FIG##1##Figure 2## and in the <xref ref-type="sec" rid="s3">Discussion</xref>).</p> <p>Despite Bcd's role in AP positioning during early <italic>Drosophila</italic> development, <italic>bcd</italic> heterozygotes (<italic>bcd<sup>E1</sup>/+</italic>) are highly viable embryos, in which, for example, among 593 embryos only 4% were unhatched, and no larval head defects were found ##REF##16150708##[46]##. In these mutants, we found that Hb was shifted anteriorly (16.6% change, from 47.0% EL in WT to 39.2% EL), as previously reported ##REF##3383244##[8]##,##REF##11845210##[50]##, but sharpness was not affected (0.5% change; from 80.2° in WT to 80.6°; ##FIG##0##Figure 1F##).</p> <p>Mutant <italic>bcd</italic> genes encoding proteins specifically defective in cooperative DNA binding have been isolated by Hanes and collaborators ##REF##16150708##[46]##,##REF##11124901##[54]##, using a genetic screen in yeast. It was shown that these mutations do not disrupt the DNA recognition or transcriptional activity of Bcd. We used one of these mutants, <italic>bcd<sup>K57R</sup></italic>, to analyze the effect of Bcd cooperative binding on Hb sharpness and position. In embryos with one dosage of the <italic>bcd<sup>K57R</sup></italic> allele (in a null Bcd background), incomplete penetrance gives two outcomes: embryos with a weak response to the Bcd defect (##FIG##0##Figure 1G##; showing a small anterior shift compared to <italic>bcd<sup>E1</sup>/+</italic>, from 39.2% EL to 37.6% EL); and embryos with a strong response, which have a large anterior shift compared to <italic>bcd<sup>E1</sup>/+</italic> (##FIG##0##Figure 1H##; from 39.2% EL to 24.8% EL). In both cases, sharpness is not reduced (0.5% and 2.5% change from WT, for weak and strong mutants, respectively).</p> <p>Driever et al. ##REF##2502714##[43]## used a series of lacZ constructs to describe the effect of high and low affinity binding sites for the establishment of the localized zygotic expression domains. These constructs have a lacZ coding sequence attached to different fragments of the <italic>hb</italic> promoter sequence. To make flies carrying these constructs, they are introduced into the fly genome, but the WT background is maintained. In one of these constructs, pThb5, the promoter is a portion of the native sequence, having 6 Bcd and 2 Hb binding sites (see ##SUPPL##8##Figure S4##). These constructs do not show self-regulation, since the protein encoded by them (β-galactosidase) has no transcriptional activity. pThb5 expression has significantly reduced sharpness (from 80.2° in WT to 72.7°, ##TAB##0##Table 1##; see ##SUPPL##9##Figure S5## and ##SUPPL##10##Figure S6## for lacZ expression patterns) and a slight anterior shift in position (from 47.0% in WT to 45.7% EL). In addition, we found that pThb5 expression is sharper than <italic>hb</italic>
<sup>14F</sup> (72.7° vs. 63.9°), indicating that Hb protein, normally expressed in these flies, can increase the sharpness of pThb5 expression (This effect of the Hb WT expression was predicted by our model, as shown below.)</p> <title>A Reaction–Diffusion Network Model Reproduces the Dynamics of Hb Pattern Formation</title> <p>To investigate what causes the changes in positioning and sharpness shown in ##FIG##0##Figure 1## and ##TAB##0##Table 1##, we developed a predictive reaction–diffusion model. This <italic>hunchback</italic> self-regulatory (HSR) model captures both Bcd cooperative binding and <italic>hb</italic> self-regulation, with six Bcd sites and two Hb sites driving <italic>hb</italic> expression. This model readily reproduced the phenotypes of the WT and mutant experiments, allowing us to make predictions and gain new understanding of the molecular mechanisms producing the measured macroscopic patterns. The model reactions are summarized in ##FIG##2##Figure 3##. <italic>hb</italic> expression requires two steps: Bcd and Hb protein binding to <italic>hb</italic> promoter (reversible reactions (2<italic>n</italic>, <italic>n</italic> = 1, …, 6) and (14, 16), respectively); and Hb protein synthesis (irreversible reactions (1+2<italic>n</italic>, <italic>n</italic> = 1, …, 8)). Reaction (1) represents Bcd production; reactions (18) and (19) represent Bcd and Hb degradation, respectively.</p> <p>We used Fick's Law to describe Bcd and Hb diffusion, and the Law of Mass Action ##UREF##10##[55]## for the reactions in ##FIG##2##Figure 3##, to derive a system of coupled partial differential equations (PDE's; see ##SUPPL##2##Text S1##) for the species <italic>B</italic>, <italic>H</italic>, <italic>b<sub>0</sub>–b<sub>6</sub></italic>, <italic>h<sub>0</sub>–h<sub>2</sub></italic> (only species B and H are allowed to diffuse). The dynamics of the system are completely described by initial conditions, diffusion coefficients and the kinetic constants, found by fitting model output to experimentally measured expression patterns.</p> <p>We used a finite difference method to solve the model PDE's (see ##SUPPL##2##Text S1##) and a steepest descent method to determine the <italic>k</italic> parameters, by fitting the total Bcd ([<italic>B</italic>]<sub>T</sub>) and Hb ([<italic>H</italic>]<sub>T</sub>) concentrations (##SUPPL##2##Text S1##, Equations S1 and S2) to the respective patterns of an embryo in mid-nuclear cleavage cycle 14 (about 36.4 minutes into the cycle, ##FIG##1##Figure 2A##; see <xref ref-type="sec" rid="s4">Materials and Methods</xref> for age determination method). The Bcd gradient was fit first, by using a zero initial concentration and fitting the model to Bcd experimental data (##FIG##1##Figure 2A##) without the Hb reactions (##FIG##2##Figure 3##, reactions 3, 5, 7, 9, 11, and 13–18), and with Bcd production (##FIG##2##Figure 3##, reaction 1) only at 9% EL (the position of maximum Bcd level in the data). The Hb reactions were then taken into account, and the model fit to the Hb experimental data (##FIG##1##Figure 2A##). With the <italic>k</italic>'s determined (##SUPPL##0##Table S1##), we simulated the dynamics of the HSR model (##FIG##1##Figure 2B##), using zero initial Hb concentration (embryos lacking maternal Hb develop normally ##REF##3803711##[42]##; also, final concentrations are largely unaffected by low initial [H<sub>T</sub>]—see <xref ref-type="sec" rid="s3">Discussion</xref>). The HSR model qualitatively reproduced the time development of the Hb pattern (##FIG##1##Figure 2B##, 59 embryos of different ages), even though it was fit to only one pair of Bcd-Hb patterns (##FIG##1##Figure 2A##; see ##SUPPL##11##Figure S7## for a direct comparison between data and model time evolution profiles). The match was best with respect to sharpness increase (##FIG##1##Figure 2D##); for position, the model shifted more than observed (##FIG##1##Figure 2C##; see also ##SUPPL##11##Figure S7A and S7B##), possibly reflecting the simplified aspect of the model, like the number of Hb or Bcd sites. ##FIG##1##Figure 2## indicates that developmental age, after the transient behavior in the first 10 minutes of cycle 14, is not a significant factor in sharpness or position.</p> <title>Spatial Bistability, Produced by <italic>hb</italic> Self-Regulation, is Responsible for the Sharp Hb Border</title> <p>Once we established that the HSR model accounts for wild-type expression, we could analyze its dynamics to determine what is responsible for converting the smooth Bcd spatial distribution into the sharp Hb pattern. We chose Zero Eigenvalue Analysis ##UREF##11##[56]##–##UREF##14##[59]## as a technique for searching parameter values that produce bistability in our model. This method establishes a set of restrictions (the sign compatible relations, see ##SUPPL##3##Text S2##) which must be met in order for bistability to occur. Using this technique on a simplified HSR model ##UREF##15##[60]##, with reduced number of Bcd binding sites (reactions 4 to 13 removed, ##FIG##2##Figure 3##) and normal Hb binding, we demonstrated that the model does have bistable solutions (see ##SUPPL##3##Text S2##). Although this analysis was performed for a well-mixed system, where concentrations are assumed to be uniform, the bistable behavior is also found in numerical solutions of the full model, where spatial distributions of concentrations are considered. The bifurcation diagram (##FIG##3##Figure 4A##), plotting [<italic>H</italic>]<sub>T</sub> steady-state solutions for various k<sub>0,B</sub> and [<italic>B</italic>]<sub>T</sub>, shows that for Bcd concentrations within the bistable region (green line), Hb concentration has two stable solutions (blue lines), separated by an unstable solution (red line). This bistability causes a very large change in Hb concentration (ΔHb, from one stable branch to the other) over a small change in Bcd concentration (ΔBcd) as it passes through a threshold (at the anterior boundary of the bistable region). In the <italic>Drosophila</italic> embryo, the Bcd gradient provides different concentrations along the anterior–posterior axis, which work like different initial conditions in the well-mixed system. It creates a sharp Hb border (##FIG##0##Figure 1C##) at the position where Bcd crosses this threshold. Comparison between <italic>hb</italic> RNA and protein profiles shows that RNA patterns are sharper than protein ones (84.2° and 82.5° mean sharpness, respectively; ##SUPPL##12##Figure S8## and ##TAB##0##Table 1##). These findings combine with our results above, to indicate that Hb sharpness results from spatially bistable dynamics, which is due to <italic>hb</italic> self-regulation. In the HSR model, if the <italic>hb</italic> self-regulation reactions (##FIG##2##Figure 3##, reactions 14–17) are removed, the network loses bistability and has only a single steady state (##FIG##3##Figure 4B##), in which Hb varies smoothly with Bcd.</p> <p>With Hb bistability characterized in the model, we can proceed to simulating the macroscopic behaviors (i.e., expression phenotypes) shown in ##FIG##0##Figure 1##, by altering binding strengths and site numbers, to reproduce the corresponding mutant genotype. For example, to model <italic>hb</italic> expression in the <italic>bcd</italic>
<sup>E1</sup>/+ heterozygotes we reduced the Bcd source term (k<sub>0,B</sub>) to 66.6% of original (see ##SUPPL##1##Table S2##), finding anteriorly shifted Hb pattern (16.5% change from WT simulation) without disrupting sharpness (1.0% change from WT), in agreement with data (##FIG##0##Figure 1F##, heavy blue curve; ##TAB##0##Table 1##). Hb sharpness was maintained because reduction of Bcd production by this amount does not change the model's bistable phase portrait, ##FIG##3##Figure 4A##, but does shift the position of the Bcd threshold anteriorly. To simulate loss of self-regulation, we removed reactions 14–17 (##FIG##2##Figure 3##) from the full model (6B2H, for 6 Bcd and 2 Hb sites), to give 6B0H sites. The model predicts a loss of sharpness (19.0% change from WT; from 84.3° to 68.3°; ##FIG##0##Figure 1D##, heavy blue curve, ##TAB##0##Table 1##) in qualitative agreement with <italic>hb<sup>14F</sup></italic> homozygote experimental data (which showed a 20.3% change from WT). The bifurcation diagram for 6B0H sites (##FIG##3##Figure 4B##) shows that this loss of Hb sharpness is due to loss of bistability, since Hb concentration becomes a smoothly decreasing function of Bcd concentration.</p> <p>The pThb5 construct contains an estimated six active Bcd sites and two Hb sites (##SUPPL##8##Figure S4##), but it does not exhibit self-regulation because the protein coded by it has no transcriptional activity. To reproduce the lacZ expression of pThb5 we derived an extra set of reactions by replacing H with L (LacZ) in reactions (1+2<italic>n</italic>, <italic>n</italic> = 1, …, 8; and 18), and added these reactions to the full model; with this model (6B2H_lacZ) we have no self-activation for lacZ, but we still take into account both Hb sites in the lacZ promoter. We found no shift in position, but did find an 11.1% loss of sharpness from the WT simulation (##TAB##0##Table 1##, ##SUPPL##7##Figure S3##), in agreement with the experimental loss of sharpness from WT to lacZ (9.3%). This indicates that the loss of sharpness for pThb5 is caused by the lack of self-regulation. The lacZ experimental patterns are sharper than <italic>hb</italic>
<sup>14F</sup> (72.7° and 63.9°, respectively), probably due to the Hb sites in the construct promoter region. Our model predicts this effect, with sharpness for 6B2H_lacZ (74.9°) higher than sharpness for 6B0H (68.3°).</p> <title>Bcd Controls Hb Position, not Sharpness</title> <p>The loss of sharpness for <italic>hb<sup>14F</sup></italic> homozygotes, ##FIG##0##Figure 1D## and ##TAB##0##Table 1##, demonstrates that cooperative Bcd binding is not sufficient to generate the sharp Hb border, since Bcd cooperativity is not affected in this <italic>hb</italic> allele. Bcd cooperative binding does play an important role in Hb pattern positioning, however, as demonstrated by the <italic>bcd<sup>K57R</sup></italic> results (##FIG##0##Figure 1G and 1H##; ##TAB##0##Table 1##; see also ##REF##2502714##[43]##), in which pattern position is altered without affecting sharpness. We simulated the reduced cooperativity in <italic>bcd<sup>K57R</sup></italic> embryos ##REF##11124901##[54]## by reducing the Bcd binding constants in reactions (2<italic>n</italic>, <italic>n</italic> = 1, …, 6) by dropping <italic>factor</italic> to 67% and 95% of original, to simulate the strong and weak mutants, respectively (##FIG##2##Figure 3## caption, ##SUPPL##1##Table S2##). ##FIG##3##Figure 4C## shows that this reduction in cooperativity shifts the bistable region towards regions of high Bcd concentration (##FIG##3##Figure 4C##, red and black lines), anteriorly shifting the pattern without disrupting its sharpness, in agreement with the data (##FIG##0##Figure 1G##, heavy blue line; ##FIG##0##Figure 1H##, heavy blue line; ##TAB##0##Table 1##). These results show that Bcd cooperative binding controls the position at which the Hb bistable switch occurs.</p> <p>Though not sufficient for sharpness, Bcd cooperativity is necessary for Hb bistability to produce sharpness. We can demonstrate this by simulating a strong decrease in cooperativity <italic>in silico</italic> by a decrease in the number of Bcd binding sites (removing reactions 10–13, giving 4B2H sites) without affecting self-regulation. These simulations show a strong reduction in <italic>hb</italic> activation, giving both a strong anterior shift and a drop in sharpness (##FIG##3##Figure 4D##). The small box in ##FIG##3##Figure 4D## shows that bistability was not disrupted, since self-regulation was not affected, but the Bcd threshold was shifted to a very high concentration not reached by the Bcd gradient (indicated by red line). This result indicates that Bcd cooperative binding is necessary for <italic>hb</italic> activation to reach its bistable threshold, which in turn is necessary for sharpness to occur.</p> <p>The above results show that small disruptions of Bcd cooperative binding result in positional shifts, without loss of sharpness, while large enough disruptions of cooperative binding also disrupt Hb sharpness, since the bistable switch is not reached. However, the bistable switch itself can only be produced by <italic>hb</italic> self-regulation.</p>
<title>Discussion</title> <p>
<italic>hb</italic> is one of the most studied developmental genes in <italic>Drosophila melanogaster</italic>. Since its discovery, many aspects of its transcriptional regulation and roles in establishing the segmented body plan of the embryo have been studied. We have shown, using a predictive reaction–diffusion model as well as a series of experimental data, that <italic>hb</italic> reads out the positional information of the morphogenetic Bcd gradient with bistable kinetics, resulting from <italic>hb</italic> self-regulation. The sharp Hb border is generated by the transition from one stable Hb concentration to the other, and the AP position where this switch occurs is determined by a threshold in the Bcd concentration, which establishes the <italic>hb</italic> activation level through cooperative binding. Our model reproduces the time development of the sharp border formation and predicts that if <italic>hb</italic> self-regulation is removed, bistability will be lost and sharpness disrupted. The expression pattern of the <italic>hb</italic> self-regulation mutant, <italic>hb<sup>14F</sup></italic>, confirms this prediction. Our analysis confirms earlier observations of the role of Bcd cooperative binding in positioning, but establishes that this mechanism is not sufficient for sharpness. Here, we discuss these conclusions in more detail in the context of the literature.</p> <title>Bcd Cooperativity</title> <p>The role of Bcd cooperative binding on Hb positioning has been demonstrated since 1989 ##REF##2567637##[9]##,##REF##2911348##[34]##,##REF##2502714##[43]##. For example, Driever et al. ##REF##2502714##[43]## used a selection of artificial lacZ constructs, each containing some portion of the native <italic>hb</italic> regulatory sequence, with varying numbers of high and low affinity Bcd binding sites. They showed that reducing the number or strength of Bcd binding sites shifted the lacZ expression anteriorly, demonstrating the role of cooperativity for pattern positioning. Increasing the binding strength or number of sites gave posterior shifts and sharper borders, suggesting that cooperativity could also be responsible for Hb sharpness. However, even for the construct with the highest number of sites (12; 6 strong and 6 weak), which showed the strongest expression level, pattern was not as sharp as wild-type Hb. While a role for cooperative binding in sharpening was suggested by these results, the authors noted this could not be firmly concluded from their data ##REF##2502714##[43]##. Struhl at al. ##REF##2567637##[9]## also observed shallower than endogenous Hb borders with a series of similar lacZ constructs. Based on their results, Driever et al. ##REF##2502714##[43]## proposed a ‘gradient-affinity model’, wherein target genes with high affinity binding sites, like <italic>hb</italic>, would be efficiently expressed even at low Bcd concentrations, and target genes containing low affinity binding sites would be positioned in more anterior positions.</p> <p>Our model reproduces these positioning effects of Bcd cooperative binding, as shown in the simulations in which binding site strength was varied (e.g., ##FIG##0##Figure 1F–H##). It also explains that the border of lacZ expression patterns is not as sharp as wild-type Hb because of the loss of self-regulation in such constructs (##SUPPL##9##Figure S5## and ##SUPPL##10##Figure S6##). These results indicate that bistability can play a role in the gradient-affinity model, since they show that changing the cooperativity level shifts the Hb pattern (##FIG##0##Figure 1G and 1H##) but does not change its sharpness, allowing an on–off expression boundary to be placed at many positions in the embryo.</p> <p>Instead of changing cooperativity by changing the binding sites in artificial constructs, Burz and Hanes ##REF##11124901##[54]## generated several Bcd cooperativity mutants, such as the <italic>bcd</italic>
<sup>K57R</sup> used in this study. They showed that this mutant is stable <italic>in vivo</italic> (in yeast cells) and is not affected in its DNA recognition, nuclear entry, or transcriptional activity characteristics ##REF##16150708##[46]##,##REF##11124901##[54]##, and <italic>in situ</italic> hybridization showed that localization and expression of <italic>bcd<sup>K57R</sup></italic> mRNA is normal ##REF##16150708##[46]##. Through analyzing the expression pattern of the gap genes <italic>hb</italic>, <italic>gt</italic> and <italic>Kr</italic> in this mutant, Lebrecht et al. ##REF##16150708##[46]## showed that cooperative Bcd binding is critical for embryonic patterning. That study also reported a reduction in Hb sharpness, contrary to what we report here (see ##FIG##0##Figure 1G and 1H## and ##TAB##0##Table 1##). In ##REF##16150708##[46]##, sharpness (slope) was calculated on non-normalized data, which makes their results susceptible to artificial variations in gene expression levels that can occur at many steps in a staining assay, such as embryo fixation and permeabilization. Using non-normalized data and measuring the slope by the quotient between profile intensity, Δ<italic>y</italic>, by the distance along AP axis, Δ<italic>x</italic>, introduces variability in slope due to variability in Δ<italic>y</italic>, see ##SUPPL##7##Figure S3##. Here, we present new results for Hb sharpness, computing slope with normalized data (Δ<italic>y</italic> = 1) to reduce the contribution of these experimental errors and to compare mutants with different levels of expression. Our method, similar to that used by Crauk and Dostatni ##REF##16271865##[49]##, depends only on how far it takes Hb to drop from ‘on’ to ‘off’ expression (Δ<italic>x</italic>). Recently, Gregor et al. ##REF##17632061##[61]## have shown that immunofluorescent signals are proportional to protein concentration plus a nonspecific background. This indicates that normalizing immunofluorescent signals provides an equivalent result to normalizing real concentrations. With this approach, we do not find a significant difference in sharpness between <italic>bcd</italic>
<sup>K57R</sup>, <italic>bcd</italic>
<sup>WT</sup>, or <italic>bcd</italic>
<sup>E1</sup>/+.</p> <p>Gregor et al. ##REF##17632062##[48]## recently presented quantitative data comparing Bcd and Hb intensities from whole embryos, and analyzed the precision of this input/output relation. They fit the Hill equation (Equation 1) to the Bcd/Hb input/output relation and estimated that Bcd binds to the <italic>hb</italic> promoter with a Hill coefficient of 5; somewhat higher than the in vitro values ##REF##8620846##[44]##,##REF##9774343##[47]##, but within the known number of Bcd binding sites ##REF##2567637##[9]##,##REF##2911348##[34]##. However, they neglected the contribution of <italic>hb</italic> self-regulation in establishing the levels of Hb protein, so the value reported for Bcd<sub>1/2</sub> would not produce half-maximal Hb synthesis in the absence of self-regulation (see ##SUPPL##4##Text S3## for more details). In other words, to reach maximum Hb production without self-regulation would require a higher Bcd concentration, and a higher value for Bcd<sub>1/2</sub>. Calling this corrected value Bcd<sub>1/2</sub>
<sup>Cor</sup>, Bcd<sub>1/2</sub>
<sup>Cor</sup>&gt;Bcd<sub>1/2</sub>, which means that ln(Bcd/Bcd<sub>1/2</sub>
<sup>Cor</sup>)&lt;ln(Bcd/Bcd<sub>1/2</sub>), and using Equations 3 and 4 shows that <italic>n</italic>
<sup>Cor</sup>&gt;<italic>n</italic> (see ##SUPPL##4##Text S3## for derivation of these equations). This indicates that a corrected Hill coefficient (<italic>n</italic>
<sup>Cor</sup>) should be higher than that reported by Gregor at al. and likely higher than the six Bcd binding sites known for <italic>hb</italic> regulation ##REF##2567637##[9]##,##REF##2911348##[34]##, making a claim that Bcd cooperativity determines Hb sharpness unlikely. In ##REF##17632062##[48]## it was argued that the effect of additional factors, uncorrelated with Bcd, would require <italic>hb</italic> readout of Bcd concentration to be even more reliable than reported. However, this argument does not justify neglecting <italic>hb</italic> self-regulation for establishing levels of Hb protein, because Hb precision (measured by standard deviation) is not correlated with Hb protein levels, as shown in ##REF##17632062##[48]##. In other words, <italic>hb</italic> self-regulation may not influence the precision of the readout process, but it does determine protein levels, and this role cannot be neglected in calculating the Hill coefficient.



</p> <title>The Requirement of Bistability for Sharpness</title> <p>Our HSR model, though fit to WT data, predicts the loss of sharpness we found experimentally in the self-regulatory mutant, <italic>hb</italic>
<sup>14F</sup> (##FIG##0##Figure 1D## and ##SUPPL##5##Figure S1##). The experimental results directly support the need for <italic>hb</italic> self-regulation for sharp pattern development. <italic>hb</italic>
<sup>14F</sup> is a lack of function mutant, generated by Lehmann et al. ##REF##3803711##[42]##, which forms normal mRNA ##REF##8087852##[51]## but has a truncated protein with no DNA binding capacity. The protein is stable, persisting into central nervous system development ##REF##11525736##[62]##, and has been visualized with Hb antibody ##REF##2846287##[53]## staining at lower intensity than WT (##SUPPL##6##Figure S2##; we measured intensity at 10–20% WT, comparable to the 10% reported in ##REF##2846287##[53]##). Since the Bcd protein or its binding are not affected in these mutants, these embryos clearly show that Bcd cooperative binding is not sufficient for producing Hb sharpness. Our results with <italic>hb</italic>
<sup>14F</sup> agree with the observations of Houchmandzadeh et al. ##REF##11845210##[50]## that expression in the <italic>hb<sup>6N</sup></italic> allele (also with non-functional protein) suggested a role for self-regulation in sharpening. Our model shows that loss of self-regulation disrupts the bistable behavior in <italic>hb</italic>
<sup>14F</sup> expression, resulting in the loss of sharpness. Similarly, expression of the pThb5 lacZ construct shows reduced sharpness in comparison to WT (##SUPPL##9##Figure S5##, ##SUPPL##10##Figure S6##, and ##TAB##0##Table 1##), since the protein coded by lacZ (β-galactosidase) is not self-regulatory. Construct sharpness is greater than <italic>hb</italic>
<sup>14F</sup> sharpness, however. The model predicts this, by taking into account that native, patterned Hb protein can bind in the construct promoter (compare ##FIG##0##Figure 1D## to ##SUPPL##9##Figure S5## and ##SUPPL##10##Figure S6##; see also ##TAB##0##Table 1## and ##SUPPL##8##Figure S4##).</p> <p>Crauk and Dostatni ##REF##16271865##[49]## recently reported sharp expression for a lacZ construct containing only three (strong) Bcd binding sites. We found lowered sharpness for the pThb5 construct (##FIG##0##Figure 1H##), which contains six Bcd binding sites (three strong, three weak; as well as two Hb and Kr sites). Since this is opposite to what any cooperative effect should be for increasing sites, the differences are likely to be methodological. We used both whole mount fluorescent <italic>in situ</italic> hybridization (FISH; ##SUPPL##9##Figure S5##) and traditional in situ hybridization (##SUPPL##10##Figure S6##; as used in ##REF##16271865##[49]##) to visualize lacZ expression. Both methods gave similar measures of sharpness (##TAB##0##Table 1##), but the enzymatic staining is more susceptible to signal saturation and tends not to be proportional to RNA concentration. Crauk and Dostatni ##REF##16271865##[49]## also reported reduced sharpness in embryos with truncated Bcd proteins, <italic>Bcd-ΔC</italic> and <italic>Bcd-ΔQC</italic> (with specific defects in protein activity). In light of the <italic>hb</italic> transcriptional dynamics found in our analysis, we believe such alterations to Bcd could cause transcription to remain sub-threshold for bistable activation, similar to ##FIG##3##Figure 4D##.</p> <title>The Method Used for Embryo Staging</title> <p>All embryos used in ##FIG##0##Figure 1## and ##TAB##0##Table 1## were in nuclear cleavage cycle 14A (precellular blastoderm), within time classes T5 to T6 (26 to 39 min into this cycle), during which Hb levels are at their highest (##FIG##0##Figure 1D## shows two T7 embryos, to show the normal posterior patterning in these mutants). We staged each embryo by established methods ##REF##3151484##[63]##, following dorsal membrane invagination measured from images obtained by Differential Interference Contrast (DIC) optics. For comparing Hb dynamics and model simulation (##FIG##1##Figure 2##), we used embryos in the first 36.4 minutes of cycle 14A. While Hb expression dynamically amplifies over this period (##FIG##1##Figure 2B##), the mature sharpness is reached within 5–10 minutes, after which it is stable (##FIG##1##Figure 2D##). Note that ##FIG##1##Figure 2D## shows that Hb pattern in all WT embryos older than 8 min in cycle 14 are sharper than any embryo in ##FIG##0##Figure 1D##.</p> <title>Maternal Hb</title> <p>The earliest Hb protein pattern in the embryo is of maternal origin. Before nuclear cleavage cycle 8, maternal <italic>hb</italic> mRNA is distributed uniformly throughout the egg, but its translation is repressed by the posterior Nanos (Nos) protein gradient, resulting in a smooth anterior gradient of maternal Hb (Hb<sub>mat</sub>) protein. This is gradually substituted by the zygotically expressed Hb, starting in cycle 11 ##REF##2716822##[27]##, ##UREF##9##[33]##, ##REF##2450283##[64]##–##REF##1720354##[69]##. To see the effect of these early Hb distributions on cycle 14 dynamics, we ran HSR model simulations with an initial Hb pattern taken from cycle 13 data (##SUPPL##13##Figure S9##; parameters and data from same assay as in ##FIG##0##Figures 1C## and ##FIG##1##2A##). Simulation results were the same as in ##FIG##1##Figure 2##, indicating no effect from the Hb initial condition. This behavior is related to the bistable behavior of Hb: the diagram in ##SUPPL##14##Figure S10## shows that, inside the bistable region, only relatively high initial conditions (above the unstable branch) can produce high Hb concentrations. It indicates that the initial concentrations of Hb, determined by zygotic production, are low and not sufficient to carry the system through the transition from lower to upper stable branch. This agrees with previous results showing that embryos lacking maternal Hb develop normally ##REF##3803711##[42]##,##REF##2704418##[65]##,##REF##2704419##[66]##. Like the <italic>Drosophila</italic> embryos, the HSR model is robust to variability in Hb initial concentration.</p> <title>The Method of Constructing Reaction-Network Models</title> <p>The method we have used to construct the HSR network, describing <italic>hb</italic> regulation by Bcd cooperative binding and <italic>hb</italic> self-regulation, can be readily applied to other genetic regulatory systems in <italic>Drosophila</italic> or other organisms, since the regulatory interactions are general. We avoided using a Hill kinetics approach to model cooperativity because this would require some questionable assumptions, such as all six sites being equal, which is counter to published Bcd binding data ##REF##2502714##[43]##,##REF##9774343##[47]##, and bound simultaneously, which is highly improbable. One advantage of using the Hill equation could be its few number of parameters; however, using relations k<sub>b(n-1),bn</sub> = <italic>factor</italic>
<sup>n</sup>.k<sub>b0,b1</sub> and k<sub>bn,b(n-1)</sub> = k<sub>b1,b0</sub> for <italic>n</italic> = 2, …, 6 in reactions (2+2<italic>n</italic>, <italic>n</italic> = 1, …, 5) allowed us to describe cooperativity with just three parameters, k<sub>b0,b1</sub>, k<sub>b1,b0</sub> and <italic>factor</italic>. The effects of more or less binding sites and more or less transcriptional regulators can easily be built into our kind of model. Our method allows for a direct link between macroscopic pattern formation and its molecular basis. As well, such a model is amenable to mathematical analysis with modern nonlinear techniques, which have developed rapidly in recent years ##UREF##6##[17]##, ##UREF##10##[55]##, ##UREF##16##[70]##–##UREF##19##[73]##. In the present example, using such techniques to search for multiple steady states allowed us to identify the bistability inherent in the self-regulatory reactions, and determine the model parameters necessary for triggering this.</p> <title>The Bistable Behavior</title> <p>In reaction networks, bistability is frequently verified by changing the initial concentration of one species, in the well-mixed system, and checking the concentrations of all other species when the system reaches the stationary state. In a monostable regime, small variations in the initial concentration generally produce small variations in the stationary state. However, if the concentration is in the vicinity of a threshold, where the transition from the monostable to the bistable regimes occurs, small changes in the initial condition can produce large variation in the stationary state, because the concentrations of the species can follow a completely different trajectory (i.e., sequence of intermediate concentrations), ending up in a very different stationary state (##FIG##3##Figure 4A##). In our spatially-patterned case, the anterior–posterior Bcd gradient provides many different concentrations that work like different initial conditions in the well-mixed system; in such a way that at the position where Bcd crosses the threshold (##FIG##3##Figure 4A##), the Hb stationary concentration changes abruptly, producing the sharp Hb border (##FIG##1##Figure 2A##).</p> <p>The origin of bistability in the Bcd-Hb system is a consequence of the positive feedback of <italic>hb</italic> self-regulation. If Hb production is not high enough, self-regulation can only produce a small change in Hb production, and the consequent increase in degradation counteracts almost all increase in Hb production; this regime occurs in the posterior half of the embryo. If Hb production is more effectively increased, by increasing Bcd concentration, the positive feedback can produce a certain additional amount of Hb protein, which can be sufficient to start increasing Hb production more efficiently. If this occurs, the additional amount of Hb will increase the feedback even more strongly, ending up in a completely different regime, having higher Hb concentrations; this regime occurs in the anterior half of the embryo.</p> <title>Concluding Remarks</title> <p>Our data and model show that positioning and sharpness of the Hb pattern are separable processes. With the <italic>hb<sup>14F</sup></italic> allele and the pThb5 construct, we show that sharpness can be disrupted with self-regulation defects; and our theoretical analysis suggests this is due to loss of bistability. Earlier work has suggested many of the shifting and sharpening effects we find here. However, there has been debate about the relative roles of the transcriptional regulators: some studies have suggested a role for <italic>hb</italic> self-regulation in sharpening ##REF##11845210##[50]##, while others indicate that it could be completely controlled by Bcd ##REF##17632062##[48]##,##REF##16271865##[49]##. It has also been known that the number of Bcd binding sites in the <italic>hb</italic> promoter affects pattern position ##REF##2567637##[9]##,##REF##2911348##[34]##,##REF##2502714##[43]##. Our data and model offer a synthesis: positioning is largely dependent on the Bcd occupation states of the <italic>hb</italic> promoter, but sharpening is a result of bistability in the <italic>hb</italic> activation dynamics, caused by <italic>hb</italic> self-regulation. Bcd cooperativity, through affecting <italic>hb</italic> transcription, determines the threshold at which bistability occurs, but is not itself sufficient for sharpening.</p> <p>In 1977, Lewis at al. ##REF##859349##[74]## used theoretical arguments to suggest that bistable control can account for the interpretation of gradients in positional information. More recently, bistability has been found in many complex biological processes ##REF##9888852##[14]##–##REF##10672420##[16]##, ##REF##9486845##[18]##–##REF##16864795##[21]##,##REF##15208707##[75]## and spatial bistability has been proposed in dorso-ventral patterning in <italic>Drosophila</italic>
##REF##16864795##[21]##,##REF##15759004##[22]##. Here, we have combined experiments, modeling and analysis to suggest that this dynamic feature may also be central to AP patterning, and that for <italic>hb</italic> transcription bistability arises from the convergence of two regulatory mechanisms (Bcd cooperative binding and <italic>hb</italic> self-regulation). This provides a specific mechanism to the earlier indication that Bcd and Hb synergy is required for <italic>Drosophila</italic> gap patterning ##REF##8087852##[51]##. Moreover, in agreement with Lewis et al. ##REF##859349##[74]##, our findings indicate that bistability may be central to threshold-dependent reading mechanisms of the positional information established by smooth maternal signals.</p> <p>Our approach, of developing a kinetic transcriptional model from molecular data such as binding sites and regulatory interactions (repression or activation), using dynamical systems theory to determine the model dynamics, and confirming the model predictions against quantitative experiments, could be used for uncovering regulatory mechanisms in many other pattern formation systems, in fruit flies and in other organisms.</p>
<p>Conceived and designed the experiments: FJPL. Performed the experiments: FJPL. Analyzed the data: FJPL FMCV DMH PMB AVS. Contributed reagents/materials/analysis tools: FJPL FMCV PMB AVS. Wrote the paper: FJPL DMH AVS. Built the model: FJPL FMCV PMB.</p> <p>During embryonic development, the positional information provided by concentration gradients of maternal factors directs pattern formation by providing spatially dependent cues for gene expression. In the fruit fly, <italic>Drosophila melanogaster</italic>, a classic example of this is the sharp on–off activation of the <italic>hunchback</italic> (<italic>hb</italic>) gene at midembryo, in response to local concentrations of the smooth anterior–posterior Bicoid (Bcd) gradient. The regulatory region for <italic>hb</italic> contains multiple binding sites for the Bcd protein as well as multiple binding sites for the Hb protein. Some previous studies have suggested that Bcd is sufficient for properly sharpened Hb expression, yet other evidence suggests a need for additional regulation. We experimentally quantified the dynamics of <italic>hb</italic> gene expression in flies that were wild-type, were mutant for <italic>hb</italic> self-regulation or Bcd binding, or contained an artificial promoter construct consisting of six Bcd and two Hb sites. In addition to these experiments, we developed a reaction–diffusion model of <italic>hb</italic> transcription, with Bcd cooperative binding and <italic>hb</italic> self-regulation, and used Zero Eigenvalue Analysis to look for multiple stationary states in the reaction network. Our model reproduces the <italic>hb</italic> developmental dynamics and correctly predicts the mutant patterns. Analysis of our model indicates that the Hb sharpness can be produced by spatial bistability, in which <italic>hb</italic> self-regulation produces two stable levels of expression. In the absence of self-regulation, the bistable behavior vanishes and Hb sharpness is disrupted. Bcd cooperative binding affects the position where bistability occurs but is not itself sufficient for a sharp Hb pattern. Our results show that the control of Hb sharpness and positioning, by <italic>hb</italic> self-regulation and Bcd cooperativity, respectively, are separate processes that can be altered independently. Our model, which matches the changes in Hb position and sharpness observed in different experiments, provides a theoretical framework for understanding the data and in particular indicates that spatial bistability can play a central role in threshold-dependent reading mechanisms of positional information.</p> <title>Author Summary</title> <p>Pattern formation during embryonic development, or morphogenesis, is one of the most intriguing problems in biology, entailing the sequence of processes by which a relatively simple system, the fertilized egg, becomes a mature organism. In these processes, the genetic information, stored at the molecular scale in the DNA, is translated into the macroscopic spatial expression patterns that precede the tissue–organ scale of body organization. It can also be understood as a flux of information from the genetic to the organ–system level. In the fruit fly <italic>Drosophila melanogaster</italic>, one of the early processes during its embryonic development is the formation of the sharp Hunchback protein pattern. To generate this pattern, the <italic>hunchback</italic> gene interprets the position-dependent information in the shallow maternal Bicoid gradient and converts it into the sharp Hunchback protein pattern. We propose that bistability in the dynamics of <italic>hunchback</italic> gene regulation can account for this information reading process, and we show that this bistable mechanism can be produced by the ability of this gene to regulate its own expression. The solution of this problem offers new approaches to understand the phenomenon of morphogenesis.</p>
<title>Supporting Information</title>
<p>We thank John Reinitz for fruitful discussions on the model and sharing of lab facilities; Bloomington Stock Center for <italic>hb<sup>14F</sup></italic> flies; Steve Hanes for the <italic>bcd</italic>
<sup>K57R</sup> flies; Steve Small for <italic>hb</italic> RNA probes; Ah-Ram Kim for lacZ probes; Gerald Thomsen for revising the manuscript; Carlos Alonso, Ah-Ram Kim, and Lisa Prazak for instructions on experiments; Svetlana Surkova for helping embryo temporal classification; and Eliana Abdelhay, José Pontes, Hector Barrabin, and Nathalie Dostatni for discussions on results.</p>
<fig id="pcbi-1000184-g001" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1000184.g001</object-id><label>Figure 1</label><caption><title>Control of Hb sharpness and position.</title><p>(A,B) An embryo in mid-nuclear cleavage cycle 14A, with immunostaining to Bicoid (Bcd) (A) and Hunchback (Hb) (B) proteins. Anterior pole on left, dorsal side on top. (C) Fluorescence profiles versus anteroposterior position for A (green) and B (blue). Hb position, 47.1% EL; sharpness, 82.7°. (D) Hb profiles for homozygotes of <italic>hb<sup>14F</sup></italic>, an allele coding a non-DNA binding Hb protein ##REF##8001780##[52]##, showing dramatically reduced sharpness (63.9°). Heavy blue line: <italic>hunchback</italic> self-regulatory (HSR) model (##FIG##2##Figure 3##) prediction for absence of self-regulation. (E) <italic>hb<sup>14F</sup></italic> heterozygotes and wild-type together, showing similar position (44.3% EL) and sharpness (81.3°). See ##SUPPL##6##Figure S2## for non-normalized data. (F) Hb profiles from <italic>bcd<sup>E1</sup></italic>/+ embryos (Bcd mRNA half-dosage). Heavy blue line: simulation for this background, by reducing Bcd synthesis in the HSR model. (G,H) Hb profiles from embryos expressing one copy of <italic>bcd<sup>K57R</sup></italic>, an allele affecting Bcd cooperativity ##REF##11124901##[54]##, gives two outcomes: a small anterior shift ((G); 3.0% change from <italic>bcd<sup>E1</sup></italic>/+); and a strong anterior shift ((H); 36.7% change from <italic>bcd<sup>E1</sup></italic>/+). Heavy blue lines: HSR model simulations for weakly and strongly reduced Bcd cooperative binding ((G) and (H), respectively). Maximum intensities are normalized to one, to allow comparison of profile sharpness from different experiments. All embryos are between 26 and 39 minutes into nuclear cleavage cycle 14A, as determined by membrane invagination and the relative position from surface to cortex (see <xref ref-type="sec" rid="s4">Materials and Methods</xref>). In (D) there are two T7 embryos, showing normal posterior pattern (not used for sharpness or position measurements). Individual embryo images are shown in ##SUPPL##5##Figure S1##.</p></caption></fig> <fig id="pcbi-1000184-g002" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1000184.g002</object-id><label>Figure 2</label><caption><title>Dynamics of Hunchback (Hb) pattern formation.</title><p>(A) Fitting the <italic>hunchback</italic> self-regulatory (HSR) model (##FIG##2##Figure 3##) to the same profiles as in ##FIG##0##Figure 1C##. Red and black lines: experimental Bicoid (Bcd) and Hb, respectively. Blue and green lines: fitting for total Bcd and Hb concentrations (Equations S1.1–2, in ##SUPPL##2##Text S1##), respectively. (B) 47 Hb wild-type profiles at different ages (red lines). Black line is the oldest one (same as in (A)). Blue lines: Dynamic simulation of HSR model, using kinetic constants determined in (A). Green line is the same as in (A). (C) pattern positions change slightly over time. Red dots and blue line: Hb boundary positions measured from embryos and from simulation in (B), respectively. (D) The pattern quickly achieves mature sharpness, with little change after 10 minutes. Red dots: border inclination (from embryos in (B)) versus embryo age. Blue line: Model prediction for the progressive increase of border inclination, for the simulations in (B). Black arrows indicate the embryo used to fit the model (shown in (A)). See ##SUPPL##11##Figure S7## for a direct comparison between data and model time evolution profiles.</p></caption></fig> <fig id="pcbi-1000184-g003" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1000184.g003</object-id><label>Figure 3</label><caption><title>The <italic>hunchback</italic> (<italic>hb</italic>) self-regulatory (HSR) model to simulate <italic>hb</italic> transcriptional activation by Bicoid (Bcd) and self-regulation.</title><p>(1): Bcd (B) synthesis from a source term; (2<italic>n</italic>, <italic>n</italic> = 1, …, 6): Bcd binding to <italic>hb</italic> promoter; (14, 16): Hb (H) binding to <italic>hb</italic> promoter; (1+2<italic>n</italic>, <italic>n</italic> = 1, …, 6) and (15, 17): Hb synthesis; (18 and 19): Hb and Bcd decay, respectively. b<sub>n</sub> and h<sub>n</sub> represent the fragments of <italic>hb</italic> promoters containing 6 Bcd and 2 Hb sites, respectively; subscripts <italic>n</italic> indicate how many Bcd or Hb molecules are bound. 0 denotes either inert or constant concentration species (e.g. mRNA). k<sub>b0,b1</sub> indicates the transition from b<sub>0</sub> to b<sub>1</sub> states; k<sub>b1,H</sub> indicates production of H from b<sub>1</sub>, and so on. We introduced cooperativity by taking k<sub>b(n-1),bn</sub> = factor<sup>n</sup>.k<sub>b0,b1</sub> and k<sub>bn,b(n-1)</sub> = k<sub>b1,b0</sub> for <italic>n</italic> = 2, …, 6. In addition, we set k<sub>bn,H</sub> = (1+Synt_Factor).k<sub>b(n-1),H</sub> for <italic>n</italic> = 2, …, 6 to account for the effect of multiple protein binding to the gene promoter. These relations strongly reduced the number of parameters, and the model fitting was not sensitive to changes in them.</p></caption></fig> <fig id="pcbi-1000184-g004" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1000184.g004</object-id><label>Figure 4</label><caption><title>Bifurcation diagrams and simulations for the <italic>hunchback</italic> (<italic>hb</italic>) self-regulatory (HSR) model.</title><p>(A) Bistability: solid blue curves are the stable steady states for Hb total concentration [H]<sub>T</sub> (y-axis) for a given Bicoid (Bcd) total concentration [B]<sub>T</sub> (bottom axis), or as Bcd production rate (k<sub>0,B</sub>, top axis). Red curve is unstable steady states ([H]<sub>T</sub> will evolve away from such concentrations). Hb concentration will evolve to the low steady states (short dashed lines) at low Bcd concentration, but abruptly steps up to the high steady states (long dashed lines) as Bcd concentration, [B]<sub>T</sub>, moves out of the bistable region (31.2&lt;[B]<sub>T</sub>&lt;40.5; green line on bottom axis). (B) Loss of bistability in the absence of self-regulation: bifurcation diagram for HSR model without self-regulation reactions (14–17) shows no bistability. (C) Reduction in Bcd cooperative binding shifts bistable region toward high [B]<sub>T</sub> regions. Black and Red: bifurcation diagrams for heavy blue lines in ##FIG##0##Figure 1G and 1H##, respectively (simulations for reduced Bcd cooperative binding (##FIG##2##Figure 3##; ##SUPPL##1##Table S2##); blue line, wild-type diagram (same as in (A)); (D) Simulation for reduced number of Bcd binding sites (removing reactions 10–13) shows disrupted sharpness since Bcd concentrations (red line in the bifurcation diagram, small box) do not reach the bistable threshold (green dot).</p></caption></fig>
<table-wrap id="pcbi-1000184-t001" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1000184.t001</object-id><label>Table 1</label><caption><title>Mean positions and sharpness (inclination) for all embryos and simulations presented.</title></caption><alternatives><table frame="hsides" rules="groups"><colgroup span="1"><col align="left" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/></colgroup><thead><tr><td align="left" rowspan="1" colspan="1">Embryos Background/Simulation</td><td align="left" rowspan="1" colspan="1">Mean Position (% EL)</td><td align="left" rowspan="1" colspan="1">Position Std. Dev.</td><td align="left" rowspan="1" colspan="1">Mean Inclination (degrees)</td><td align="left" rowspan="1" colspan="1">Inclination Std. Dev.</td></tr></thead><tbody><tr><td align="left" rowspan="1" colspan="1">WT Hb pattern used for fitting</td><td align="left" rowspan="1" colspan="1">47.1</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">82.7</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">HSR Model (6B2H) Fit to WT</td><td align="left" rowspan="1" colspan="1">48.0</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">84.3</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">WT <italic>hb</italic>
</td><td align="left" rowspan="1" colspan="1">47.0</td><td align="left" rowspan="1" colspan="1">1.5</td><td align="left" rowspan="1" colspan="1">80.2</td><td align="left" rowspan="1" colspan="1">1.7</td></tr><tr><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>14F</sup> heterozygotes and WT</td><td align="left" rowspan="1" colspan="1">44.3</td><td align="left" rowspan="1" colspan="1">1.8</td><td align="left" rowspan="1" colspan="1">81.3</td><td align="left" rowspan="1" colspan="1">1.3</td></tr><tr><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>14F</sup> homozygotes</td><td align="left" rowspan="1" colspan="1">41.9</td><td align="left" rowspan="1" colspan="1">3.3</td><td align="left" rowspan="1" colspan="1">63.9</td><td align="left" rowspan="1" colspan="1">3.9</td></tr><tr><td align="left" rowspan="1" colspan="1">HSR model (6B0H) for <italic>hb</italic>
<sup>14F</sup> homozygotes</td><td align="left" rowspan="1" colspan="1">49.0</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">68.3</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">
<italic>bcd</italic> heterozygous, <italic>bcd</italic>
<sup>E1</sup>/+</td><td align="left" rowspan="1" colspan="1">39.2</td><td align="left" rowspan="1" colspan="1">0.2</td><td align="left" rowspan="1" colspan="1">80.6</td><td align="left" rowspan="1" colspan="1">1.3</td></tr><tr><td align="left" rowspan="1" colspan="1">HSR (6B2H) model for <italic>bcd</italic> heterozygotes.</td><td align="left" rowspan="1" colspan="1">40.1</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">83.4</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">
<italic>bcd</italic>
<sup>K57R</sup> weak</td><td align="left" rowspan="1" colspan="1">37.6</td><td align="left" rowspan="1" colspan="1">2.3</td><td align="left" rowspan="1" colspan="1">80.6</td><td align="left" rowspan="1" colspan="1">1.3</td></tr><tr><td align="left" rowspan="1" colspan="1">HSR model (6B2H) for <italic>bcd</italic>
<sup>K57R</sup> weak</td><td align="left" rowspan="1" colspan="1">38.0</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">83.8</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">
<italic>bcd</italic>
<sup>K57R</sup> strong</td><td align="left" rowspan="1" colspan="1">24.8</td><td align="left" rowspan="1" colspan="1">2.9</td><td align="left" rowspan="1" colspan="1">82.2</td><td align="left" rowspan="1" colspan="1">1.2</td></tr><tr><td align="left" rowspan="1" colspan="1">HSR model (6B2H) for <italic>bcd</italic>
<sup>K57R</sup> strong</td><td align="left" rowspan="1" colspan="1">24.4</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">84.3</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">pThb5 lacZ (FISH)</td><td align="left" rowspan="1" colspan="1">45.7</td><td align="left" rowspan="1" colspan="1">1.8</td><td align="left" rowspan="1" colspan="1">72.7</td><td align="left" rowspan="1" colspan="1">0.9</td></tr><tr><td align="left" rowspan="1" colspan="1">pThb5 lacZ (<italic>in situ hybridization</italic>)</td><td align="left" rowspan="1" colspan="1">47.3</td><td align="left" rowspan="1" colspan="1">2.8</td><td align="left" rowspan="1" colspan="1">74.9</td><td align="left" rowspan="1" colspan="1">2.4</td></tr><tr><td align="left" rowspan="1" colspan="1">HSR model for pThb5 (6B2H_lacZ)</td><td align="left" rowspan="1" colspan="1">48.0</td><td align="left" rowspan="1" colspan="1">NA</td><td align="left" rowspan="1" colspan="1">74.9</td><td align="left" rowspan="1" colspan="1">NA</td></tr><tr><td align="left" rowspan="1" colspan="1">WT <italic>hb</italic> RNA (simultaneous FISH Immunofluorescence)</td><td align="left" rowspan="1" colspan="1">45.0</td><td align="left" rowspan="1" colspan="1">1.0</td><td align="left" rowspan="1" colspan="1">84.2</td><td align="left" rowspan="1" colspan="1">1.4</td></tr><tr><td align="left" rowspan="1" colspan="1">WT Hb protein (simultaneous FISH Immunofluorescence)</td><td align="left" rowspan="1" colspan="1">45.4</td><td align="left" rowspan="1" colspan="1">1.0</td><td align="left" rowspan="1" colspan="1">82.5</td><td align="left" rowspan="1" colspan="1">1.1</td></tr></tbody></table></alternatives></table-wrap> <table-wrap id="pcbi-1000184-t002" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1000184.t002</object-id><label>Table 2</label><caption><title>Methods for obtaining expression patterns, for the specified numbers of WT and mutant embryos.</title></caption><alternatives><table frame="hsides" rules="groups"><colgroup span="1"><col align="left" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/><col align="center" span="1"/></colgroup><thead><tr><td align="left" rowspan="1" colspan="1">Method</td><td colspan="6" align="left" rowspan="1">Genes</td></tr><tr><td align="left" rowspan="1" colspan="1"/><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>WT</sup>
<italic>bcd</italic>
<sup>WT</sup>
</td><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>14F</sup>
<italic>bcd</italic>
<sup>WT</sup>
</td><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>WT</sup>
<italic>bcd</italic>
<sup>K57R</sup>
</td><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>WT</sup>
<italic>bcd</italic>
<sup>E1</sup>/+</td><td align="left" rowspan="1" colspan="1">
<italic>lacZ</italic>
</td><td align="left" rowspan="1" colspan="1">
<italic>hb</italic>
<sup>WT</sup>
<italic>RNA hb</italic>
<sup>WT</sup>
<italic>prot.</italic>
</td></tr></thead><tbody><tr><td align="left" rowspan="1" colspan="1">FISH<xref ref-type="table-fn" rid="nt102">a</xref>
</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">9</td><td align="left" rowspan="1" colspan="1">13</td></tr><tr><td align="left" rowspan="1" colspan="1">Immunofluorescence<xref ref-type="table-fn" rid="nt103">b</xref>
</td><td align="left" rowspan="1" colspan="1">50</td><td align="left" rowspan="1" colspan="1">39</td><td align="left" rowspan="1" colspan="1">23</td><td align="left" rowspan="1" colspan="1">7</td><td align="left" rowspan="1" colspan="1">
<xref ref-type="table-fn" rid="nt105">d</xref>
</td><td align="left" rowspan="1" colspan="1">13</td></tr><tr><td align="left" rowspan="1" colspan="1">Enzymatic staining<xref ref-type="table-fn" rid="nt104">c</xref>
</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">–</td><td align="left" rowspan="1" colspan="1">21</td><td align="left" rowspan="1" colspan="1">
<xref ref-type="table-fn" rid="nt105">d</xref>
</td></tr></tbody></table></alternatives></table-wrap>
<disp-formula><label>(1)</label></disp-formula> <disp-formula><label>(2)</label></disp-formula> <disp-formula><label>(3)</label></disp-formula> <disp-formula><label>(4)</label></disp-formula>
<supplementary-material content-type="local-data" id="pcbi.1000184.s001"><label>Table S1</label><caption><p>Kinetic constants for model fitting and simulations.</p><p>(0.02 MB PDF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s002"><label>Table S2</label><caption><p>Kinetic constants to simulate a weak and strong response to half-dosage of <italic>bcd</italic>
<sup>K57R</sup>.</p><p>(0.01 MB PDF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s003"><label>Text S1</label><caption><p>Model equations and fitting procedure.</p><p>(0.03 MB PDF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s004"><label>Text S2</label><caption><p>Analysis for multiple stable stationary states.</p><p>(0.14 MB PDF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s005"><label>Text S3</label><caption><p>Obtaining the logarithmic Hill equation.</p><p>(0.02 MB PDF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s006"><label>Figure S1</label><caption><p>Some individual embryo images for the profiles in ##FIG##0##Figure 1##. Each of the overlays in ##FIG##0##Figure 1D–H## is composed of intensity profiles, along the AP axis, from individual embryo images. Here, we show some examples of these individuals, with expression patterns for <italic>hb</italic>; individual embryos on top and profiles on bottom. (A) An embryo homozygous for <italic>hb<sup>14F</sup></italic>, one of the profiles used in ##FIG##0##Figure 1D##. (B) An embryo from the overlay in ##FIG##0##Figure 1E##. (C) An embryo with a half dosage of Bicoid (Bcd) mRNA (<italic>bcd<sup>E1</sup></italic>/+; see ##FIG##0##Figure 1F##); (D) An embryo with a weak <italic>bcd<sup>K57R</sup></italic> phenotype (see ##FIG##0##Figure 1G##); (E) An embryo with a strong <italic>bcd<sup>K57R</sup></italic> phenotype (see ##FIG##0##Figure 1H##).</p><p>(1.59 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s007"><label>Figure S2</label><caption><p>Non-normalized profiles for 39 embryos expressing the <italic>hb<sup>14F</sup></italic> allele. The <italic>hb<sup>14F</sup></italic> homozygotes can be easily identified by low signal intensities, as described in the literature ##REF##2846287##[53]##. WT and heterozygotes could not be easily distinguished in this way.</p><p>(0.55 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s008"><label>Figure S3</label><caption><p>The influence of experimental error on sharpness measurement. Δ<italic>y</italic>
<sub>ex</sub> is the variability due to staining procedure, such as embryo fixation and permeabilization, and Δ<italic>x</italic> is the projection on the AP axis. Even though P1 and P2 have the same AP projection (Δ<italic>x</italic>) they will have different sharpness measurements (β&gt;α), but this difference is caused by the experimental error Δ<italic>y</italic>
<sub>ex</sub>. With normalized patterns (Δ<italic>y</italic>
<sub>ex</sub> = 0 and Δ<italic>y</italic> = 1), sharpness measurement will be determined only by the AP projection.</p><p>(0.18 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s009"><label>Figure S4</label><caption><p>Hb promoter and lacZ artificial construct. A fragment of 4776 bp from the <italic>hb</italic> gene. Hb (blue) and Bcd (red) sites were identified by DNAse footprinting in ##REF##2797150##[35]## and ##REF##2502714##[43]##, respectively. Black arrows indicate the transcription initiation sites. Green arrows indicate the position of the fragment used in the pThb5 lacZ construct, which includes the six Bcd and both Hb sites. To avoid an unnecessarily extensive model, we took into account only two Hb sites.</p><p>(0.26 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s010"><label>Figure S5</label><caption><p>Loss of sharpness shown by FISH images for the mRNA expression pattern of the pThb5 lacZ construct, driven by a six Bcd and two Hb site promoter. (A) An individual embryo on top and its corresponding profile on bottom. (B) Overlay of 9 profiles. The pattern shows a slight anterior shift (45.7% EL from 47% EL in WT) and 9.3% decrease in sharpness (72.7° from 80.2° in WT) Heavy line: Model simulation taking into account both Hb sites on this construct (6B2H_lacZ); pattern positioned at 48.0% EL with 74.9° inclination.</p><p>(0.50 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s011"><label>Figure S6</label><caption><p>Demonstrating reduced sharpness in the pThb5 lacZ construct with the alternate method of enzymatic in-situ hybridization (compare to FISH, used in ##SUPPL##9##Fig. S5##). (A) Embryo image taken by Differential Interference Contrast (DIC) optics on top and corresponding intensity profile on bottom. (B) 21 profiles taken by the same procedure. Pattern position gives 47.3% EL and 74.9° sharpness. Heavy blue line is the simulation of HSR model (6B2H_lacZ) model for this construct (48.0% EL and 74.9° inclination). The observed reduction in sharpness is independent of the method used for visualizing mRNA.</p><p>(0.63 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s012"><label>Figure S7</label><caption><p>Temporal evolution of the Hb pattern. Multicolored and dark blue lines indicate experimental data and HSR model, respectively. Embryos are the same as in ##FIG##2##Figure 3B–D##. The embryo ages are indicated in the legend. The oldest embryo, 36.4 min, is the same used to fit the model (##FIG##2##Fig. 3A##). Each plot covers a temporal window of 9.1 min: (A) Embryos having 9.1 minutes or less in cycle 14; (B), (C) and (D) show embryos having ages from 9.2 to 18.2, 18.3 to 27.3 and 27.4 to 36.4 min, respectively. The earliest computed pattern corresponds to age 2.3 min. There is an interval of roughly 2.3 min between each computed pattern; four of them are shown in each plot. See <xref ref-type="sec" rid="s4">Materials and Methods</xref> for embryo temporal classification.</p><p>(1.16 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s013"><label>Figure S8</label><caption><p>
<italic>hb</italic> RNA and protein profiles. Simultaneous fluorescence immunostaining and in situ hybridization showing <italic>hb</italic> RNA and protein. (A) RNA (above) and protein (below) signals from two different embryos (right and left); (B) RNA (solid lines) and protein (dashed lines) from left (blue) and right (red) embryos in (A), respectively. (C) Comparison between 13 RNA and protein patterns, left and middle chart, respectively. The difference in sharpness is indicated, on right chart, by two straight lines having inclinations of 84.2° (red) and 82.5° (black), which correspond to RNA and protein mean inclinations, respectively (see ##TAB##0##Table 1##). The high sharpness for RNA data supports our result that Hb sharpness results from spatially bistable dynamics.</p><p>(1.11 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s014"><label>Figure S9</label><caption><p>The HSR model is robust for different Hb initial concentrations. Running the HSR model with an initial condition of zero Hb concentration (yellow line) or a Hb pattern taken from cycle 13 (gray line) gives the same result (solid blue and dotted red lines, respectively), indicating no effect from different initial conditions. This is in agreement with previous results reporting that embryos lacking maternal Hb develop normally ##REF##3803711##[42]##, ##REF##2704418##[65]##–##REF##2704419##[66]##.</p><p>(0.14 MB TIF)</p></caption></supplementary-material> <supplementary-material content-type="local-data" id="pcbi.1000184.s015"><label>Figure S10</label><caption><p>Bistability diagram showing the origin of Hb robustness for different initial conditions. Green line on the Bcd concentration axis indicates the bistable region, where for each Bcd concentration there are three stationary states available for Hb concentration. Blue and red lines are the stable and unstable branches of stationary states, respectively. Black circles indicate initial conditions outside the bistable region (green line); they evolve to the upper stable branch. Red and green circles indicate initial conditions inside the bistable region: All initial Hb conditions below the unstable branch (red circles) evolve to the lower stable branch, having low Hb concentrations; and initial Hb concentrations above the unstable branch (green circles) evolve to the upper stable branch, having high Hb concentrations. The arrows indicate the temporal evolution.</p><p>(0.22 MB TIF)</p></caption></supplementary-material>
<table-wrap-foot><fn id="nt101"><p>All embryos are between 26 and 39 minutes into nuclear cleavage cycle 14A. NA means not applicable.</p></fn></table-wrap-foot> <table-wrap-foot><fn id="nt102"><label>a</label><p>Whole mount fluorescent in situ hybridization.</p></fn><fn id="nt103"><label>b</label><p>Fluorescence immunostaining.</p></fn><fn id="nt104"><label>c</label><p>Whole mount in situ hybridization with digoxygenin-labeled RNA.</p></fn><fn id="nt105"><label>d</label><p>All lacZ embryos were co-stained for Bcd and Hb proteins.</p></fn></table-wrap-foot> <fn-group><fn fn-type="COI-statement"><p>The authors have declared that no competing interests exist.</p></fn><fn fn-type="financial-disclosure"><p>FJPL was supported by the Brazilian National Research Council (CNPq)/Ministry of Science and Technology (MCT, Brazil), the joint National Science Foundation (NSF)/National Institute of General Medical Sciences (NIGMS) BioMath program, grant 1-R01-GM072022, and National Institutes of Health grant 2-ROI-RR07801. FMCV and PMB were supported by the CNPq/MCT, Brazil. DMH and AVS were supported by the joint NSF/NIGMS BioMath program, grant 1-R01-GM072022.</p></fn></fn-group>
<graphic xlink:href="pcbi.1000184.g001"/> <graphic xlink:href="pcbi.1000184.g002"/> <graphic id="pcbi-1000184-t001-1" xlink:href="pcbi.1000184.t001"/> <graphic xlink:href="pcbi.1000184.g003"/> <graphic xlink:href="pcbi.1000184.g004"/> <graphic xlink:href="pcbi.1000184.e001.jpg" mimetype="image" position="float"/> <graphic xlink:href="pcbi.1000184.e002.jpg" mimetype="image" position="float"/> <graphic xlink:href="pcbi.1000184.e003.jpg" mimetype="image" position="float"/> <graphic xlink:href="pcbi.1000184.e004.jpg" mimetype="image" position="float"/> <graphic id="pcbi-1000184-t002-2" xlink:href="pcbi.1000184.t002"/>
<media xlink:href="pcbi.1000184.s001.pdf"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s002.pdf"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s003.pdf"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s004.pdf"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s005.pdf"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s006.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s007.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s008.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s009.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s010.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s011.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s012.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s013.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s014.tif"><caption><p>Click here for additional data file.</p></caption></media> <media xlink:href="pcbi.1000184.s015.tif"><caption><p>Click here for additional data file.</p></caption></media>
[{"label": ["1"], "element-citation": ["\n"], "surname": ["Driesch"], "given-names": ["H"], "year": ["1893"], "article-title": ["Entwicklungsmechanische studien VI. \u00dcber einige fragen der theoretischen morphologie."], "source": ["Z Wiss Zool"], "volume": ["55"], "fpage": ["34"], "lpage": ["62"]}, {"label": ["2"], "element-citation": ["\n"], "surname": ["Browne"], "given-names": ["EN"], "year": ["1909"], "article-title": ["The production of new hydranths in hydra by the insertion of small grafts."], "source": ["J Exp Zool"], "volume": ["7"], "fpage": ["1"], "lpage": ["37"]}, {"label": ["5"], "element-citation": ["\n"], "surname": ["Lander", "Nie", "Wan"], "given-names": ["AD", "Q", "FYM"], "year": ["2002"], "article-title": ["Do morphogen gradients arise by diffusion?"], "source": ["Dev Biol"], "volume": ["247"], "fpage": ["471"], "lpage": ["471"]}, {"label": ["10"], "element-citation": ["\n"], "surname": ["Wolpert"], "given-names": ["L"], "year": ["2002"], "source": ["Principles of Development. 2nd ed"], "publisher-loc": ["New York"], "publisher-name": ["Oxford University Press"]}, {"label": ["12"], "element-citation": ["\n"], "surname": ["Bonn", "Meunier", "Greffier", "Al-Kahwaji", "Kellay"], "given-names": ["D", "J", "O", "A", "H"], "year": ["1998"], "article-title": ["Bistability in non-Newtonian flow: rheology of lyotropic liquid crystals."], "source": ["Phys Rev E"], "volume": ["58"], "fpage": ["2115"], "lpage": ["2118"]}, {"label": ["13"], "element-citation": ["\n"], "surname": ["Sessoli", "Gatteschi", "Caneschi", "Novak"], "given-names": ["R", "D", "A", "MA"], "year": ["1993"], "article-title": ["Magnetic bistability in a metal-ion cluster."], "source": ["Nature"], "volume": ["365"], "fpage": ["141"], "lpage": ["143"]}, {"label": ["17"], "element-citation": ["\n"], "surname": ["Feinberg", "Horn"], "given-names": ["M", "FJM"], "year": ["1974"], "article-title": ["Dynamics of open chemical systems and algebraic structure of underlying reaction network."], "source": ["Chem Eng Sci"], "volume": ["29"], "fpage": ["775"], "lpage": ["787"]}, {"label": ["25"], "element-citation": ["\n"], "surname": ["J\u00e4ckle", "Tautz", "Schuh", "Seifert", "Lehmann"], "given-names": ["H", "D", "R", "E", "R"], "year": ["1986"], "article-title": ["Cross-regulatory interactions among the gap genes of "], "italic": ["Drosophila"], "source": ["Nature"], "volume": ["324"], "fpage": ["668"], "lpage": ["670"]}, {"label": ["26"], "element-citation": ["\n"], "surname": ["Frohnhofer", "Nussleinvolhard"], "given-names": ["HG", "C"], "year": ["1986"], "article-title": ["Organization of anterior pattern in the "], "italic": ["Drosophila", "bicoid"], "source": ["Nature"], "volume": ["324"], "fpage": ["120"], "lpage": ["125"]}, {"label": ["33"], "element-citation": ["\n"], "surname": ["Tautz", "Lehmann", "Schnurch", "Schuh", "Seifert"], "given-names": ["D", "R", "H", "R", "E"], "year": ["1987"], "article-title": ["Finger protein of novel structure encoded by Hunchback, a 2nd member of the gap class of "], "italic": ["Drosophila"], "source": ["Nature"], "volume": ["327"], "fpage": ["383"], "lpage": ["389"]}, {"label": ["55"], "element-citation": ["\n"], "surname": ["Horn", "Jackson"], "given-names": ["F", "R"], "year": ["1972"], "article-title": ["General mass action kinetics."], "source": ["Arch Ration Mech Anal"], "volume": ["47"], "fpage": ["81"], "lpage": ["116"]}, {"label": ["56"], "element-citation": ["\n"], "surname": ["Li"], "given-names": ["HY"], "year": ["1998"], "article-title": ["Zero eigenvalue analysis for the determination of multiple steady states in reaction networks."], "source": ["Z Naturforsch [A]"], "volume": ["53"], "fpage": ["171"], "lpage": ["177"]}, {"label": ["57"], "element-citation": ["\n"], "surname": ["Li"], "given-names": ["HY"], "year": ["1999"], "article-title": ["Determination of multiple steady states in an active membrane transport model."], "source": ["Z Naturforsch [A]"], "volume": ["54"], "fpage": ["245"], "lpage": ["250"]}, {"label": ["58"], "element-citation": ["\n"], "surname": ["Li"], "given-names": ["HY"], "year": ["1998"], "article-title": ["Determination of multiple steady states in a family of allosteric models for glycolysis."], "source": ["J Chem Phys"], "volume": ["109"], "fpage": ["8485"], "lpage": ["8493"]}, {"label": ["59"], "element-citation": ["\n"], "surname": ["Chuang", "Chao", "Ho", "Li"], "given-names": ["GS", "AC", "PY", "HY"], "year": ["2004"], "article-title": ["Computational multiple steady states in enzymatically catalyzed oxidation of monophenols by tyrosinase in an isothermal CFSTR."], "source": ["Biochem Eng J"], "volume": ["21"], "fpage": ["161"], "lpage": ["170"]}, {"label": ["60"], "element-citation": ["\n"], "surname": ["Lopes"], "given-names": ["FJP"], "year": ["2002"], "source": ["Analise da Formacao de Padroes no Desenvolvimento da "], "italic": ["Drosophila melanogaster"], "publisher-loc": ["Rio de Janeiro"], "publisher-name": ["Universidade Federal do Rio de Janeiro"]}, {"label": ["70"], "element-citation": ["\n"], "surname": ["Horn"], "given-names": ["F"], "year": ["1972"], "article-title": ["Necessary and sufficient conditions for complex balancing in chemical-kinetics."], "source": ["Arch Ration Mech Anal"], "volume": ["49"], "fpage": ["172"], "lpage": ["186"]}, {"label": ["71"], "element-citation": ["\n"], "surname": ["Feinberg"], "given-names": ["M"], "year": ["1987"], "article-title": ["Chemical-reaction network structure and the stability of complex isothermal reactors. 1. The deficiency-zero and deficiency-one theorems."], "source": ["Chem Eng Sci"], "volume": ["42"], "fpage": ["2229"], "lpage": ["2268"]}, {"label": ["72"], "element-citation": ["\n"], "surname": ["Craciun", "Feinberg"], "given-names": ["G", "M"], "year": ["2005"], "article-title": ["Multiple equilibria in complex chemical reaction networks: I. The injectivity property."], "source": ["SIAM J Appl Math"], "volume": ["65"], "fpage": ["1526"], "lpage": ["1546"]}, {"label": ["73"], "element-citation": ["\n"], "surname": ["Li", "Ho"], "given-names": ["HY", "PY"], "year": ["2000"], "article-title": ["Subnetwork analysis for the determination of multiple steady states in complex reaction networks."], "source": ["Ind Eng Chem Res"], "volume": ["39"], "fpage": ["3291"], "lpage": ["3297"]}, {"label": ["78"], "element-citation": ["\n"], "surname": ["Moodley", "Murrell"], "given-names": ["K", "H"], "year": ["2004"], "article-title": ["A colour-map plugin for the open source, Java based, image processing package, ImageJ."], "source": ["Comput Geosci"], "volume": ["30"], "fpage": ["609"], "lpage": ["618"]}, {"label": ["79"], "element-citation": ["\n"], "surname": ["Tsonis", "Kumar", "Elsner", "Tsonis"], "given-names": ["AA", "P", "JB", "PA"], "year": ["1996"], "article-title": ["Wavelet analysis of DNA sequences."], "source": ["Phys Rev E"], "volume": ["53"], "fpage": ["1828"], "lpage": ["1834"]}, {"label": ["80"], "element-citation": ["\n"], "surname": ["Golyandina", "Nekrutkin", "Zhigliavskii\u0131 \u02d8"], "given-names": ["N", "VV", "A"], "year": ["2001"], "source": ["Analysis of Time Series Structure: SSA and Related Techniques"], "publisher-loc": ["Boca Raton (Florida)"], "publisher-name": ["Chapman & Hall/CRC"]}, {"label": ["81"], "element-citation": ["\n"], "surname": ["Alexandrov", "Golyandina", "Spirov"], "given-names": ["T", "N", "A"], "year": ["2008"], "article-title": ["Singular spectrum analysis of gene expression profiles of early Drosophila embryo: exponential-in-distance patterns."], "source": ["Res Lett Signal Process"], "volume": ["2008"], "fpage": ["5"]}, {"label": ["83"], "element-citation": ["\n"], "surname": ["Fuks", "Tabachnikov"], "given-names": ["DB", "S"], "year": ["2007"], "source": ["Mathematical Omnibus: Thirty Lectures on Classic Mathematics"], "publisher-loc": ["Providence (Rhode Island)"], "publisher-name": ["American Mathematical Society"]}]
PLoS Comput Biol. 2008 Sep 26; 4(9):e1000184