Genetic studies have revealed that segment determination in Drosophila melanogaster is based on hierarchical regulatory interactions among maternal coordinate and zygotic segmentation genes. The gap gene system constitutes the most upstream zygotic layer of this regulatory hierarchy, responsible for the initial interpretation of positional information encoded by maternal gradients. We present a detailed analysis of regulatory interactions involved in gap gene regulation based on gap gene circuits, which are mathematical gene network models used to infer regulatory interactions from quantitative gene expression data. Our models reproduce gap gene expression at high accuracy and temporal resolution. Regulatory interactions found in gap gene circuits provide consistent and sufficient mechanisms for gap gene expression, which largely agree with mechanisms previously inferred from qualitative studies of mutant gene expression patterns. Our models predict activation of Kr by Cad and clarify several other regulatory interactions. Our analysis suggests a central role for repressive feedback loops between complementary gap genes. We observe that repressive interactions among overlapping gap genes show anteroposterior asymmetry with posterior dominance. Finally, our models suggest a correlation between timing of gap domain boundary formation and regulatory contributions from the terminal maternal system.
THE segmented body plan of Drosophila melanogaster becomes determined during the first 3 hr of embryogenesis (Simcox and Sang 1983). The genetics of segment determination in the Drosophila blastoderm is very well understood (see Akam 1987; Ingham 1988, for review). Saturation mutagenesis screens have enabled the isolation of a complete or almost complete set of segmentation genes (Nüsslein-Volhard and Wieschaus 1980; Nüsslein-Volhard et al. 1987). The zygotic segmentation gene network is a hierarchical dynamical system whose regulatory layers consist of gap, pair-rule, and segment-polarity genes (Nüsslein-Volhard and Wieschaus 1980). Initial conditions for zygotic segmentation gene expression are given by gradients of the maternal proteins Bicoid (Bcd; Figure 1, A and D), Hunchback (Hb; Figure 1, B and E), and Caudal (Cad; Figure 1, C and F; see St. Johnston and Nüsslein-Volhard 1992, for review). Further maternal input is provided by the terminal maternal system, which regulates segmentation gene expression in the pole regions of the embryo through the zygotic terminal gap genes tailless (tll) and huckebein (hkb; Weigel et al. 1990). In this study, we focus on the regulation of the gap genes hunchback (hb), Krüppel (Kr), knirps (kni), and giant (gt), which are expressed in broad domains during the late blastoderm stage (Figure 1, G–L).
Detailed genetic and molecular studies have yielded considerable information on the regulatory interactions underlying gap gene expression. Still, our current knowledge of gap gene regulation is incomplete. This is partly due to the ambiguity or absence of experimental data on particular regulatory interactions. However, it is also due to methodological issues concerning the inference of regulatory interactions based on the qualitative study of mutant gene expression in multicellular organisms (cf. Reinitz and Sharp 1995). These issues are rooted in the complexity and the essentially quantitative nature of the dynamical mechanisms of spatial pattern formation. Each blastoderm nucleus has different initial concentrations of maternal gene products and hence different initial conditions for zygotic gene expression. This leads to widely and qualitatively different dynamics of zygotic gene expression in different nuclei despite the fact that the underlying regulatory network is the same in each nucleus. A change in the initial conditions in maternal mutants, or in the regulatory circuitry in zygotic mutants, can have unexpected and counterintuitive effects, making interpretation of mutant gene expression patterns a highly nontrivial task in all but the simplest cases.
We illustrate the difficulties in interpreting mutant expression patterns with the example of the regulatory effect of Hb on Kr. The anterior boundary of the central Kr domain is shifted anteriorly in hb mutants (Jäckle et al. 1986), while Kr expression is generally weaker than that in wild-type embryos (Pankratz et al. 1989). Moreover, embryos overexpressing hb show posterior expansion of the central Kr domain (Hülskamp et al. 1990). Finally, Kr expression is absent in embryos lacking both Bcd and Hb, but is restored in a concentration-dependent manner by reintroducing increasing dosages of Hb (Struhl et al. 1992; Schulz and Tautz 1994). It has been proposed that these effects are due to a dual regulatory role of Hb with activation of Kr at low and repression of Kr at high concentrations of Hb (Hülskamp et al. 1990; Struhl et al. 1992; Schulz and Tautz 1994).
However, the above observations can be explained equally well by indirect activation of Kr through Kni. The expression domain of kni, which encodes a repressor of Kr (Jäckle et al. 1986; Hoch et al. 1992), expands anteriorly in hb mutants (Hülskamp et al. 1990), explaining reduced levels of Kr. The slightly altered posterior gt domain in hb mutants (Eldon and Pirrotta 1991) further complicates interpretation, since Gt is a repressor of both Kr (Kraut and Levine 1991b) and kni (Eldon and Pirrotta 1991; Capovilla et al. 1992). Expression of Kr is restored to high levels in hb kni double mutants (Harding and Levine 1988), further supporting an indirect mechanism. Moreover, embryos overexpressing hb lack kni expression altogether (Kraut and Levine 1991b), and posterior extension of the Kr domain in these embryos resembles Kr expression in kni mutants (Jäckle et al. 1986). Finally, kni is widely expressed in embryos lacking Bcd and Hb, but is repressed in a concentration-dependent manner when Hb is reintroduced (Struhl et al. 1992), which suggests that Kr derepression in these embryos is due to increasing repression of kni.
The above example reveals three main problems for inferring regulatory mechanisms from qualitative mutant expression data. These are the problems of consistency, uniqueness, and completeness.
Consistency of a proposed regulatory mechanism can be established only by keeping track of all regulatory inputs to a specific gene (cf. Reinitz and Sharp 1995). In the case of Kr, this involves at least five different regulatory contributions (by Bcd, Cad, Hb, Kni, and Gt). Current experimental approaches, however, are limited in their ability to monitor regulatory contributions simultaneously, as such interactions are inferred from mutant expression patterns and it is rarely possible to obtain mutants in more than three genes. Moreover, mutant regulatory systems by definition have an incomplete or otherwise defective set of regulatory interactions. Thus, the regulatory structure of the wild-type network must be assembled on the basis of evidence from different experiments. The consistency of such an inferred mechanism can be established conclusively only by testing it in the intact and complete developmental system.
Another problem for interpretation of mutant expression patterns is to establish the uniqueness of a mechanism, i.e., to decide whether regulatory interactions are direct or indirect. At least two possible regulatory mechanisms can account for the effect of Hb on Kr. Both mechanisms are consistent with available experimental evidence. In such an ambiguous situation, independent evidence can be provided by molecular approaches. Both Hb and Kni have been shown to bind to the Kr regulatory region in vitro (Hoch et al. 1991, 1992), but the functional importance of such biochemical interactions can be established only in vivo. Ideally, this would be achieved by targeted mutation of transcription factor binding sites in the regulatory region of an endogenous gene. Such an experiment is technically difficult and has not yet been attempted. Alternative approaches involving reporter constructs are subject to two significant complications. First, it is often difficult to establish the regulatory equivalence of such constructs to the endogenous gene. For instance, in kni mutants the posterior boundary of the third stripe of even-skipped (eve) is intact (Frasch and Levine 1987), whereas the minimal enhancer for this stripe shows complete derepression between stripes three and seven (Small et al. 1996). Second, regulatory regions used in a construct may contain binding sites for multiple factors (see Kr above) or unknown binding sites, which leads to similar ambiguities in interpreting mutant expression patterns as in the case of the endogenous gene.
Finally, there is a fundamental issue concerning completeness of a proposed regulatory mechanism, which cannot be addressed by experimental approaches alone. The fact that all maternal and gap genes are necessary for correct gap gene expression does not prove that they are also sufficient. It is impossible to prove sufficiency of the inferred mechanism without reconstituting the system ab initio, using only well-defined ingredients. Such a reconstitution is obviously impossible by contemporary experimental methods and hence has to be attempted by using mathematical modeling and computer simulations.
The problems illustrated above show that to establish consistency, uniqueness, and completeness of a regulatory mechanism, we need a method that allows us to reconstitute wild-type gene expression patterns in silico, infer underlying regulatory interactions from these wild-type patterns, and keep track of all regulatory interactions in all nuclei at all times. The gene circuit method provides such an approach (Mjolsness et al. 1991; Reinitz and Sharp 1995; Reinitz et al. 1995, 1998). It is a data-driven mathematical modeling method whose main aim is to extract information about dynamical regulatory interactions between transcription factors from given gene expression patterns (Figure 2A; Reinitz and Sharp 1995). This is achieved in four steps: (1) formulation of a mathematical modeling framework, (2) collection of gene expression data, (3) fitting of the model to expression data to obtain regulatory parameters, and (4) biological analysis of the resulting gene circuits.
The Drosophila blastoderm permits exceptionally precise modeling, since pattern formation is a consequence of regulatory interactions among segmentation genes only. Segmentation gene mutations affect expression of other segmentation genes, but do not cause any morphological defects before the onset of gastrulation (Merrill et al. 1988). Thus, the internal state of each blastoderm nucleus can be described by concentration levels of transcription factors encoded by segmentation genes. Gap gene circuits include the genes bcd, cad, hb, Kr, gt, kni, and tll. We do not model RNA explicitly, since it has no known regulatory function in Drosophila segment determination. In addition, there is no tissue growth, and we do not have to consider intercellular signaling since nuclei are not yet surrounded by membranes during the syncytial blastoderm stage (Campos-Ortega and Hartenstein 1985). Finally, patterning systems along the anteroposterior (A-P) and the dorsoventral (D-V) axes are largely independent of each other in the segmented germ-band region of the blastoderm. Therefore, blastoderm nuclei, which are the basic objects of the gene circuit model, are arranged in a one-dimensional row along the A-P axis.
Gap gene circuits cover cleavage cycles 13 and 14A during the late syncytial blastoderm stage (Figure 2B; Foe and Alberts 1983), including most of embryonic stages four and five in Campos-Ortega and Hartenstein (1985). This covers the time between the first unambiguous detection of zygotically expressed Kr and Gt proteins in early cycle 13 (our own data and Gaul and Jäckle 1987; Eldon and Pirrotta 1991; Kraut and Levine 1991a) and the onset of gastrulation at the end of cycle 14A (Foe and Alberts 1983). All nuclei divide equally and simultaneously at the beginning of cycle 14A.
Change in concentrations of transcription factors within each nucleus is governed by regulated protein synthesis, protein decay, and diffusion between neighboring nuclei (Mjolsness et al. 1991; Reinitz and Sharp 1995). Due to the lack of an in vitro polymerase II assay for eukaryotic transcription that faithfully reproduces in vivo transcriptional regulation, it is currently impossible to formulate a gene network model that is based on mechanistic chemical kinetics of transcription. Instead, we use coarse-grained kinetic equations for protein concentrations, which approximate the exact biochemistry with a sigmoid regulation-expression function (Figure 2C; Mjolsness et al. 1991; Reinitz and Sharp 1995).
Note that the general modeling framework outlined above does not specify which specific regulatory interactions take place within a gap gene circuit. These interactions are determined by regulatory parameters that constitute a genetic interconnectivity matrix (the T matrix). Each regulatory effect of a specific transcription factor on a target gene is described by a single parameter in the T matrix (Figure 2D). The gene circuit method aims to determine regulatory parameters and thus regulatory interactions within a gene circuit from given gene expression data. In other words, we seek sets of regulatory parameters that cause the gene circuit model to produce expression patterns that resemble real gap gene expression patterns as closely as possible (Figure 2A). This is achieved by fitting the model to quantitative segmentation gene expression data.
The set of quantitative gene expression data used in this study contains data for bcd, cad, hb, Kr, kni, gt, and tll from wild-type embryos (Figure 1; Poustelnikova et al. 2004). Data and model can be compared by numerically calculating expression patterns for given time classes from the model and then evaluating the sum of squared differences between model output and expression data for each gene, nucleus, and time class for which we have data. We minimize this sum by using a global optimization method called parallel Lam simulated annealing (PLSA, Figure 2A; Chu et al. 1999). The optimization procedure results in a gene circuit, which is defined by a specific set of regulatory parameters. Due to the stochastic nature of PLSA, different gene circuits (i.e., different sets of parameters) may be obtained, which all show essentially correct gene expression patterns.
The last step of the gene circuit method is the analysis of gene circuits to gain biological insights. The most important aspect of the gene circuit method considered here is that it allows for very detailed analysis of direct regulatory interactions within a given gene network. This is achieved by studying the distribution of gene circuit parameters between different gene circuits and by graphical analysis of regulatory contributions to specific patterning features (see results and Reinitz and Sharp 1995). This method of analysis allows us to study quantitative regulatory contributions to gene regulation in any nucleus at any point in time during a simulation.
Here we present a dynamical analysis of the gap gene network that is based on gap gene circuits. We show that these circuits are able to reproduce gap gene expression patterns in the late Drosophila blastoderm at high accuracy and temporal resolution. We provide a detailed analysis of regulatory interactions involved in gap gene regulation and show that our results are largely consistent with existing experimental evidence. Our models extend current knowledge of the gap gene system in several important aspects. We predict an activating effect of Cad on Kr and clarify evidence on the effects of Hb on Kr, Kr on kni, and Gt on kni. Our results suggest that mutual repression by complementary gap genes is absolutely essential for correct gap gene expression. We observe spatial asymmetry with posterior dominance in repressive interactions among overlapping gap genes. Moreover, the gene circuit method can provide information on regulatory mechanisms that is difficult to obtain by current experimental methods. Control of the posterior boundaries of posterior kni and gt was found to involve a temporal succession of multiple repressive interactions. Finally, we report a correlation between regulatory input from the terminal maternal system and late formation of gap gene domain boundaries in the posterior region of the embryo.
MATERIALS AND METHODS
The gene circuit modeling framework:
The gene circuit modeling framework has been described in detail in Mjolsness et al. (1991) and Reinitz and Sharp (1995). The basic objects of the gene circuit model are blastoderm nuclei denoted by the index i. We consider a one-dimensional model in which nuclei are arranged in a row along the A-P axis where nuclei i − 1 and i + 1 are neighbors of nucleus i. The model has three rules governing the behavior of nuclei in time t: (1) interphase, (2) mitosis, and (3) division. Rules 1 and 2 are continuous and describe the dynamics of protein synthesis and decay within a nucleus and protein diffusion between nuclei. Rule 3 is discrete and describes how each nucleus is replaced by its two daughter nuclei upon division. The schedule for these rules is based on Foe and Alberts (1983) and is summarized in Figure 2B.
The internal state of nucleus i is described by concentrations vai of transcription factors encoded by segmentation genes denoted by index a. The change in transcription factor concentration over time, dvai/dt, depends on three processes during interphase: (1) protein synthesis, (2) protein diffusion, and (3) protein decay, represented by the summation terms on the right-hand side of Equation 1 below. During mitosis, protein synthesis is shut down and only diffusion and decay occur. Thus we write 1where N is the total number of zygotic genes in the model.
In Equation 1, Tab represents a matrix of regulatory coefficients where each coefficient Tab characterizes the regulatory effect of the product of gene b on the expression of gene a (Figure 2D). This matrix is independent of i, reflecting the fact that each nucleus contains a copy of the same genome. vBcdi is the concentration of Bcd in nucleus i. Bcd is exclusively maternal and its concentration is constant in time. The regulatory effect of Bcd on gene a is represented by the parameter ma. ha is a threshold parameter representing regulatory contributions of uniformly expressed maternal transcription factors. The relative rate of protein synthesis is then given by the sigmoid regulation-expression function , where is the total regulatory input on gene a (Figure 2C). The maximum synthesis rate for the product of gene a is given by Ra. The diffusion parameter Da(n) depends on the number of nuclear divisions n that have taken place before the current time t. Diffusion is assumed to vary inversely with the square of the distance between neighboring nuclei and this distance is halved upon nuclear division. λa is the decay rate of the product of gene a. It is related to the protein half-life of the product of gene a by .
Quantitative expression data:
D. melanogaster blastoderm stage embryos were fluorescently stained for Eve protein and two other gene products using antibodies described in Kosman et al. (1998). As secondary antibodies, we used FITC anti-guinea pig, Texas Red anti-rabbit, and Cy5 anti-rat. Laterally oriented embryos were scanned using the 16× oil immersion objective of a Leica TCS4D confocal laser microscope. Fluorescent dyes were excited with a single wavelength at a time to ensure no leakage between channels, using the BP-FITC filter for the 488-nm excitation line (FITC), the BP-60030 filter for 568 nm (Texas Red), and the RG665 filter for 647 nm (Cy5). Expression levels were normalized per gene to a relative fluorescence intensity range of 0–255 on the basis of the most intensely fluorescent pattern on each slide with multiple embryos. Embryo images were cropped to fit embryo size and aligned along the A-P axis as shown in Figure 1.
A detailed description of this processing step can be found at http://flyex.ams.sunysb.edu/flyex/proc_steps/dave.html. Embryo images were segmented to obtain tabulated protein concentrations per nucleus as follows: Binary nuclear masks were constructed by edge detection, and average protein concentrations were obtained by averaging pixel values covering each nucleus in the mask. Nuclear positions are based on centroids of nuclei in the binary mask.
Embryos were assigned to cleavage cycle 12 (time class C12, used for initial conditions of the model at t = 0.0), cycle 13 (C13), and eight time classes (T1–T8) in cycle 14A (Figure 2B). Time classification for C12 and C13 is based on embryo morphology and for T1–T8 on careful visual inspection of the highly dynamic eve expression pattern by two independent observers (D. Kosman and S. Surkova; cf. Myasnikova et al. 2001). Time classification was validated by membrane morphology (cf. Merrill et al. 1988), as well as automated classification of eve expression patterns by complex-valued neural networks (Aizenberg et al. 2002), and support-vector regression (Myasnikova et al. 2002).
Nonspecific background staining was approximated by a paraboloid and subsequently eliminated by a linear mapping of intensities that transforms fluorescence at or below background level to zero and transforms maximum fluorescence to itself (E. Myasnikova, unpublished results). Expression patterns were registered using fast dyadic wavelets to align expression patterns as closely as possible (Myasnikova et al. 2001). Only nuclei with positional values in the middle 10% along the D-V axis were used for further processing.
Each integrated expression profile is based on registered data from at least 10 embryos stained for a specific gene at a specific time class, with the exception of Kni at C13, which is based on only two embryos, and Tll, for which we did not have data earlier than T3. Nuclei were categorized into 25 (C12), 50 (C13), and 100 (T1–T8) equal-sized bins according to their position along the A-P axis (cf. Foe and Alberts 1983). Concentration values for all nuclei in each bin were averaged to yield the final integrated one-dimensional expression pattern (Figure 1; Poustelnikova et al. 2004). The concentration of Bcd is nearly constant with respect to time during cycles 13 and 14A and is based on averaged registered bcd expression data from T1–T7. Concentrations of Cad and Hb at the onset of cycle 13 are derived from expression data for cycle 12. Initial concentrations for Kr, Kni, Gt, and Tll are zero in all nuclei.
Optimization by parallel Lam simulated annealing:
PLSA was used as described in Reinitz and Sharp (1995) and Chu et al. (1999). The set of ordinary differential Equations 1 was solved numerically using a Bulirsch-Stoer adaptive-step-size solver scheme adapted from Press et al. (1992). Equations were solved to a relative accuracy of 0.1%, and solutions were tested for numerical stability. We minimize the following cost function by adjusting parameters Ra, Tab, ma, ha, Da, and λa in Equation 1:
Summation is performed over the total number of data points Nd, i.e., the number of protein measurements across all genes a, nuclei i, and time classes t.
Parameter search spaces were defined by explicit search limits for Ra, Da, and λa and a collective penalty function for Tab, ma, and ha as described in Reinitz and Sharp (1995). ha parameters of Kr, kni, gt, and hb were fixed to negative values representing a constitutive “off” state of the gene. This accelerated the annealing process considerably and slightly improved annealing results while not altering the overall quality of the resulting gene circuits. Optimization was performed in parallel on 10 2.4-Ghz Pentium P4 Xeon processors and took between 8 and 160 hr per optimization run.
Selection of gap gene circuits:
We use the root mean square (rms) score as a measure for the quality of a gene circuit. The rms represents the average absolute difference between protein concentrations in model and data. PLSA is a stochastic optimization method yielding gap gene circuits of varying quality. Gene circuits most faithfully reproducing gap gene expression were selected as follows: First, only circuits with an rms of <12.0 were considered (20 circuits out of 40). All gap gene circuits with an rms of >12.0 showed obvious pattern defects, some of them severe, such as displaced or missing expression boundaries. Second, each of the selected 20 circuits was carefully tested for patterning defects by visual inspection and plotting of squared differences between model and data for each protein and time class. The 10 resulting circuits are listed in Table 1. Unless noted otherwise, graphs shown below use circuit 28008 (Table 2), since it has no circuit-specific patterning defects and its regulatory parameters correspond to the gap gene network topology observed in a majority of circuits (compare Table 2A to Figure 4A).
Analysis of circuit parameters:
Parameter values Tab and ma were classified into three types of regulatory interaction: (1) repression for parameter values ≤ −0.005, (2) no interaction for parameter values between –0.005 and 0.005, or (3) activation for parameters ≥0.005 (see Figure 4A). The threshold of 0.005 for the “no interaction” category was chosen empirically. Interactions falling into the no interaction category usually had no detectable effect on pattern formation in gap gene circuits analyzed graphically (see below). The gap gene network topology observed in a majority of gap gene circuits (Figure 4A) is preserved if a threshold of 0.01 is used instead (data not shown).
Software and bioinformatics:
Simulator and optimization codes were implemented in C; data quantification tools were implemented in C and the Khoros image analysis environment; and gene circuit analysis and plotting tools were implemented in Perl and Java. Software and gene circuit files are available at http://flyex.ams.sunysb.edu/lab/gaps.html. Expression data (FlyEx database) are available at http://urchin.spbcas.ru/flyex and http://flyex.ams.sunysb.edu/flyex.
Ten gap gene circuits including bcd, cad, hb, Kr, gt, kni, and tll, and covering a range of 35–92% A-P position, were selected for analysis as described in materials and methods (Table 1). A comparison between model output and quantified expression data is shown in Figure 3. Most circuits show minor circuit-specific patterning defects consisting of small spurious domains or slight irregularities in specific domain boundaries (Table 1). Moreover, all gap gene circuits show slight defects in the establishment of the posterior borders of the posterior gt and hb domains and fail to reproduce the late parasegment 4 (PS4)-specific expression peak of hb (Figure 3). Finally, we observed slightly elevated expression levels of gap genes during early cycle 13 (data not shown).
Analysis of circuit parameters:
The distribution of parameter values between circuits can vary from parameter to parameter (Figure 4). Most parameters show a strong tendency toward a particular type of regulatory interaction, i.e., activation, repression, or no interaction. Figure 4A shows the gap gene network topology corresponding to genetic interactions observed in a majority of gap gene circuits (see Figure 9 for a schematic representation of the network). Although a gene circuit using average parameter values does not produce correct gap gene expression patterns (data not shown), we have found two circuits (26003, 28008) whose parameters exactly represent the topology of the majority of circuits (Table 2).
Some basic features of the gap gene network topology are immediately obvious from inspection of Figure 4A. First, Bcd and Cad generally activate zygotic gap gene expression. Second, hb, Kr, kni, and gt show autoactivation. Third, except for autoregulatory interactions and the effect of Gt on hb, all reciprocal interactions among gap genes are either zero or repressive. Especially strong constraints for mutual repression are present between Kr and gt, as well as kni and hb, which show complementary expression patterns in the region of 35–92% A-P position (Figure 1, G–L). Many repressive interactions between overlapping gap genes show weaker constraints toward repression, and we have found very weak or no dynamical constraints for repression of kni, gt, and hb by the products of their immediate anterior neighbors Kr, kni, and gt, respectively. Finally, the terminal gap gene product Tll represses all other gap genes except hb.
Graphical analysis of gap gene regulation:
Graphical analysis of gap gene circuits allows us to “dissect” regulatory contributions of different transcription factors on the expression of a target gene and to characterize these interactions in great detail in space and time. To achieve this, we plot individual contributions to the sum of regulatory interactions affecting a gene's expression. Thereby, we focus on regions of expression domain boundaries. We identify regulatory factors responsible for the positioning of specific boundaries by looking for regulatory inputs that change significantly and consistently over the region of an expression domain boundary (cf. Reinitz and Sharp 1995). Consistent change implies that for boundary control by activation, the activator has to show a spatial expression gradient of the same polarity as the boundary it controls. Analogously, boundary control by repression implies a gradient of repressor with opposite polarity to the boundary it controls.
We have found activation of gap genes by Bcd and Cad in broad regions of the embryo (Figure 5). Bcd contributes strong activating inputs on the anterior domains of gt (Figure 5, A and C) and hb (Figure 5, E, F, H, and I) as well as on the central domain of Kr (Figure 5, B and D). Smaller activating inputs by Bcd can be detected in the posterior domains of kni (Figure 5, G and J) and gt (Figure 5, A and C). Three circuits (28003, 25005, and 29007) show repression of kni by Bcd, suggesting that Bcd activation might not be essential for kni expression during cycle 14A (Figure 4, A and C). The predominant maternal activating input on posterior kni and gt is provided by Cad (Figure 5, C and J). Furthermore, Cad provides a relatively strong activating input to central Kr expression (Figure 5D) and even contributes significantly to early anterior expression of hb (Figure 5H). Note that a small activating contribution of Cad on anterior hb can be detected in most gap gene circuits, but the strong early activation of hb by Cad shown in Figure 5H is exceptional. Activation in the posterior hb domain is largely due to Cad and hb autoactivation (Figure 4, A and E, and data not shown), a mechanism that we consider to be an artifact of the model (see discussion).
In addition to activation by maternal genes, zygotic gap genes show a tendency toward positive autoregulation (Figure 4). Autoactivation contributes strongly to zygotic expression of Kr, hb, and kni and can become the dominant activating contribution within an expression domain during the second half of cycle 14A (Figure 5, D, I, and J). Autoactivation of gt was found to be somewhat weaker (Figure 5C) and is not present at significant levels in all circuits (Figure 4, A and D). Note that activation in the anterior hb domain is slightly special, due to the presence of maternally expressed Hb protein in the anterior half of the embryo (Figure 1, B and E), which causes exceptionally strong autoactivation of hb early in cycle 14A (Figure 5H).
Whereas activation of gap genes by maternal genes occurs in rather broad regions, repressive interactions among gap genes provide spatially specific regulatory input for boundary positioning. Note that Kr and gt have mutually exclusive expression patterns in the blastoderm (Figure 1, G–I, and Figure 6A). Kr shows repression by Gt in all circuits (Figure 4, A and B). This repressive interaction is involved in positioning both anterior and posterior boundaries of central Kr expression during cycle 14A (Figure 6C). Although repression by Gt is quite strong, the regulatory profile of Kr indicates that missing repression by Gt does not lead to significant Kr derepression outside its central domain, since total regulatory input is not elevated significantly above the 10% level of expression in the absence of Gt (Figure 6C, arrow).
Both hb and kni show overlaps of their expression domains with the central domain of Kr (Figure 1, I and L, and Figure 6B). Most circuits show repressive inputs on Kr by Hb and Kni, which are weaker than that of Gt (Figure 4, A and B). Kni is involved in setting the posterior border of the central Kr domain. Figure 6D (asterisk) shows that Kr synthesis expands posteriorly in the absence of Kni. Similarly, Hb is involved in setting the anterior border of the central Kr domain, as Kr synthesis expands anteriorly in the absence of Hb (Figure 6D, asterisk). We found one circuit (28005) in which the boundaries of the Kr domain are set exclusively by Gt. However, this caused a slight patterning defect of the posterior Kr boundary at late cycle 14A (Table 1). In addition to the repressive interactions described above, we observed strong repression of Kr by Tll in all circuits (Figure 4, A and B). This repression is not involved in setting the boundaries of the central Kr domain since it affects regulation of Kr only at the posterior pole of the embryo (data not shown).
The anterior border of the posterior kni domain (Figure 1L and Figure 7, A and B) is set by a combination of repressive inputs by Hb and Kr (Figure 7, C and D). Whereas Hb represses kni in all circuits, repression by Kr was observed in only 6 of 10 circuits (Figure 4, A and C). Gap gene circuits without repression of kni by Kr show no detectable defects in kni expression (data not shown). Regulation of the posterior border of kni reveals a dynamic succession of repressive interactions during cycle 14A (Figure 7, E and F). All circuits show diminishing repressive input on kni by Tll in the region of the posterior boundary during cycle 14A (Figure 7, E and F), as tll expression is retained only in a region posterior to 80% A-P position (compare Figure 7A with 7B). In contrast, there is increasing repression by Gt and Hb in the boundary region (Figure 7, E and F). Note that the strength of repressive inputs by Gt and Tll varies greatly between circuits (Figure 4C and Figure 7, E and F). For instance, circuit 28008 (Figure 7E) shows extraordinarily strong repression of Gt, while other circuits such as 26001 show predominant repression by Hb and Tll, with a smaller contribution by Gt (Figure 7F).
gt is expressed in two domains in the region covered by gap gene circuits (Figures 1I and 8A). The posterior boundary of the anterior domain as well as the anterior boundary of the posterior domain of gt depend almost exclusively on very strong repression by Kr (Figure 8C). We detect a small repressive contribution by Hb to the anterior gt domain. However, Hb repression is not specifically involved in positioning the posterior boundary of this domain, being uniformly distributed across it (Figure 8, E and F). In all circuits, the posterior border of posterior gt is initially established through repression by Tll (Figure 8E). During cycle 14A, repression by Tll is increasingly complemented and replaced by Hb repression (Figure 8F). We found one circuit (28002) that shows weak activation of gt by Hb. This is likely to be an artifact of this particular circuit, since its posterior domain of tll was expanded slightly anteriorly to compensate for missing repression by Hb. Only one circuit (25005) showed very weak repression of gt by Kni, whereas all other circuits showed no such interaction (Figure 4, A and D).
hb has an anterior and a posterior expression domain (Figures 1L and 8B). Regulation of hb is quite different from other gap genes in that it has only one repressive input (Figure 4, A and E). Very strong repression of hb by Kni was found in all 10 circuits (Figure 4E). This repression is involved in positioning of the posterior border of the anterior hb domain as well as the anterior border of the posterior hb domain (Figure 8D). Note that we have found no effect of Kr on hb in any gap gene circuit (Figure 4A).
Accuracy and specificity of gap gene circuits:
Some earlier models of gap gene expression did not consider the genetic nature of the underlying dynamic mechanism (Nagorcka 1988; Goodwin and Kauffman 1990; Hunding et al. 1990). Others were based on generalized genetic mechanisms, which did not consider the specific dynamics of gene regulation or details of gap gene network topology (Meinhardt 1986, 1988). As more detailed evidence became available, theoretical approaches incorporated more detailed, qualitative representations of gap gene regulation (Burstein 1995; Sanchez and Thieffry 2001; Tchuraev and Galimzyanov 2001). The gene circuit method is the only approach so far, which allows for detailed quantitative analysis of dynamic regulatory interactions among gap genes. However, earlier studies using gap gene circuits showed a high degree of variation in the distribution of regulatory parameters between circuits (Reinitz and Sharp 1995; Reinitz et al. 1995). The quantitative data set used in the current study (Poustelnikova et al. 2004) has resulted in several significant improvements. Error in gap gene expression patterns has been reduced to <5% deviation from gene expression data (Table 1), which is comparable with the experimental error in the data itself (Myasnikova et al. 2001). The dynamics of gap gene expression are now reproduced to a temporal resolution of <7 min during cycle 14A. Our models show correct timing of gap gene expression and correct extents of overlaps between neighboring gap domains, two features of gap gene expression that were not addressed in earlier studies. Moreover, gap gene circuits reproduce shifts of gap domain boundaries during cycle 14A, a phenomenon first discovered by analyzing quantitative gap gene expression data (Jaeger et al. 2004). Finally, the addition of cad and tll has allowed us to extend the region for which we obtain correct gap gene expression patterns toward the posterior pole region.
Some theoretical approaches to regulatory interactions in the segmentation gene network infer these interactions on the basis of interpretation of mutant expression patterns taken from the literature (Sanchez and Thieffry 2001; Kumar et al. 2002). A difficulty with this approach is that such models tend to reproduce the interpretations of data they are based on, rather than providing an independent interpretation. In contrast, the gene circuit method does not require any a priori assumptions about specific regulatory interactions. Instead, it attempts to reconstruct these interactions on the basis of wild-type gene expression data (Figure 2A). Given this caveat, it is noteworthy that the results of our analysis of the gap gene network are largely consistent with studies based on mutant gene expression (see below). The fact that two independent methods lead to very similar results is an important cross-validation of conclusions based on both approaches.
Fitting models with many parameters to data is always at risk of producing nonspecific results. Gap gene circuits fail to fit expression data in regions of the embryo where additional factors are required for regulation, i.e., anterior of ∼35% A-P position (cf. Reinitz et al. 1995), where gap gene regulation is known to involve head gap genes (Cohen and Jürgens 1990; Finkelstein and Perrimon 1990; Grossniklaus et al. 1994), and posterior of 92% A-P position, where activity of the terminal gap gene hkb is required (Weigel et al. 1990; Brönner and Jäckle 1991). Moreover, even though we have not obtained unique values for regulatory parameters in different circuits, we have found a strong tendency toward a specific type of regulatory interaction for most parameters (Figure 4). This suggests that gap gene circuits represent the gap gene regulatory network in a specific and reproducible way.
Mechanisms of gap gene regulation:
Although activating contributions from Bcd and Cad show some degree of localization (Figure 5), positioning of gap gene boundaries during cycle 14A is largely under the control of repressive gap-gap cross-regulatory interactions. Thereby, activation is a prerequisite for repressive boundary control, which counteracts broad activation of gap genes in a spatially specific manner (Figures 5–8). In addition, gap genes show a tendency toward autoactivation (Figure 4), which increasingly potentiates activation by Bcd and Cad during cycle 14A (Figure 5). Autoactivation is involved in maintenance of gap gene expression within given domains and sharpening of gap domain boundaries during cycle 14A. A similar, but less specific mechanism for spatially localized gene activation by maternal gradients has been proposed by Meinhardt (1988).
Regulatory loops of mutual repression create positive regulatory feedback between complementary gap genes, providing a straightforward mechanism for their mutually exclusive expression patterns. Such a mechanism of “alternating cushions” of gap domains has been proposed by Kraut and Levine (1991b) and Clyde et al. (2003). Our results suggest that this mechanism is complemented by repression among overlapping gap genes. Overlap in expression patterns of two repressors imposes a limit on the strength of repressive interactions between them. Accordingly, repression between neighboring gap genes is generally weaker than that between complementary ones (Figure 4). Moreover, repression among overlapping gap genes is asymmetric, centered on the Kr domain (see Figure 9). Posterior of this domain, only posterior neighbors contribute functional repressive inputs to gap gene expression, while anterior neighbors do not. We show elsewhere that this asymmetry is responsible for anterior shifts of posterior gap gene domains during cycle 14A (Jaeger et al. 2004).
Repression by Tll mediates regulatory input to gap gene expression by the terminal maternal system (see Introduction). Tll provides the main repressive input to early regulation of the posterior boundary of posterior gt (Figure 8E), and activation by Tll is required for posterior hb expression (Casanova 1990; Reinitz and Levine 1990; Brönner and Jäckle 1991). Note that these two features form only during cycle 13 and early cycle 14A (Figure 3), while other gap domain boundaries are already present at the transcript level during cycles 10–12 (Knipple et al. 1985; Tautz 1988; Mohler et al. 1989; Rothe et al. 1989) and largely depend on the anterior and posterior maternal systems for their initial establishment (Gaul and Jäckle 1987; Tautz 1988; Mohler et al. 1989; Rivera-Pomar et al. 1995). The delayed formation of posterior patterning features and their distinct mode of regulation are reminiscent of segment determination in primitive dipterans and intermediate germ-band insects, supporting a conserved dynamical mechanism across different insect taxa (Tautz and Sommer 1995; Davis and Patel 2002).
The set of regulatory interactions presented here provides a consistent and sufficient dynamical mechanism for gap gene expression (see Introduction). In summary, this set of interactions consists of the following five basic regulatory mechanisms (Figure 9): (1) broad activation by Bcd and/or Cad, (2) autoactivation, (3) strong repressive feedback between mutually exclusive gap genes, (4) asymmetric repression between overlapping gap genes, and (5) feed-forward repression of posterior domain boundaries by the terminal gap gene tll. In the following subsections, we discuss evidence concerning specific regulatory interactions involved in each of these basic mechanisms in some detail.
Activation by Bcd and Cad:
Activation of gap gene expression by Bcd and Cad is supported by the following. Bcd binds to the regulatory regions of hb, Kr, and kni (Driever and Nüsslein-Volhard 1989; Driever et al. 1989; Hoch et al. 1991; Rivera-Pomar et al. 1995). The kni regulatory region also contains binding sites for Cad (Rivera-Pomar et al. 1995). The anterior domains of gt and hb are absent in embryos from bcd mothers (Tautz 1988; Eldon and Pirrotta 1991). The posterior domain of gt is missing in embryos mutant for both maternal and zygotic cad, while the posterior domain of kni is absent in embryos mutant for maternal bcd plus maternal and zygotic cad (Rivera-Pomar et al. 1995). Our results suggest partial redundancy of activation of kni by Bcd, consistent with evidence from zygotic cad embryos from bcd mothers, where maternally provided Cad is sufficient to activate kni (Rivera-Pomar et al. 1995).
Kr expression expands anteriorly in embryos from bcd mothers (Gaul and Jäckle 1987), which is due to the absence of the anterior gt and hb domains (Tautz 1988; Eldon and Pirrotta 1991; Kraut and Levine 1991b). Bcd has been shown to activate expression of Kr reporter constructs (Hoch et al. 1990, 1991), supporting an activating effect of Bcd on endogenous Kr. The fact that Kr is still expressed in embryos from bcd mutant mothers has been attributed to activation by general transcription factors (Kerrigan et al. 1991) or low levels of Hb (Hülskamp et al. 1990; Struhl et al. 1992; Schulz and Tautz 1994). In contrast, our models predict that this activation is provided by Cad (Figure 4, A and B, and Figure 5D). Although Kr expression is normal in embryos overexpressing cad (Mlodzik et al. 1990), repressive control of Kr boundaries could account for the lack of expansion of the Kr domain in such embryos.
The activating effect of Cad on hb found in gap gene circuits is likely to be spurious. The anterior hb domain is absent in embryos from bcd mutant mothers (Tautz 1988), which show uniformly high levels of Cad (Mlodzik and Gehring 1987). Moreover, the complete absence of the posterior hb domain in tll mutants (Casanova 1990; Reinitz and Levine 1990; Brönner and Jäckle 1991) suggests activation of posterior hb by Tll rather than by Cad. We believe that this spurious activation of hb by Cad is due to the absence of hkb in gap gene circuits. The posterior hb domain fails to retract from the posterior pole in hkb mutants (Casanova 1990; Brönner and Jäckle 1991), suggesting a repressive role of Hkb in regulation of the posterior hb border. Consistent with this, the posterior boundary of the posterior hb domain never fully forms in any of our circuits (Figure 3). Moreover, Tll is constrained to a very small or no interaction with hb (Figure 4E) due to the absence of the posterior repressor Hkb, since activation of hb by Tll would lead to increasing hb expression extending to the posterior pole.
A role for autoactivation in the late phase of hb regulation (Schröder et al. 1988; Hülskamp et al. 1994) is supported by the fact that the posterior border of anterior hb is shifted anteriorly in a concentration-dependent manner in embryos with decreasing doses of zygotic Hb (Simpson-Brose et al. 1994). Weakened and narrowed expression of Kr in mutants encoding a functionally defective Kr protein (Warrior and Levine 1990) suggests Kr autoactivation. Similarly, a delay in the expression of gt in mutants encoding a defective Gt protein (Eldon and Pirrotta 1991) indicates gt autoactivation. However, our results suggest that gt autoactivation is not essential. It is generally weaker than autoactivation of other gap genes (Figure 4, B–E), and circuits lacking gt autoactivation show no specific defects in gt expression. Finally, in the case of kni, there is no experimental evidence for autoactivation, while some authors have even suggested kni autorepression (Howard 1990; Rothe et al. 1994). We have not been able to detect such autorepression in any gap gene circuit (Figure 4, A and C).
Repression between complementary gap genes:
Mutual repression of gt and Kr is supported by the following. gt expression expands into the region of the central Kr domain in Kr embryos (Eldon and Pirrotta 1991; Kraut and Levine 1991a). In contrast, Kr expression is not altered in gt mutants before germ-band extension (Gaul and Jäckle 1987; Reinitz and Levine 1990; Eldon and Pirrotta 1991). However, Gt binds to the Kr regulatory region (Capovilla et al. 1992), and the central domain of Kr is absent in embryos overexpressing gt (Kraut and Levine 1991b). Moreover, Kr expression extends further anterior in hb gt double mutants than in hb mutants alone (Kraut and Levine 1991b). The above is consistent with our analysis, which shows no significant derepression of Kr in the absence of Gt even though repression of Kr by Gt is quite strong (Figure 6C).
Hb binds to the kni regulatory region, and the posterior kni domain expands anteriorly in hb mutants (Hülskamp et al. 1990; Rothe et al. 1994; Clyde et al. 2003). Embryos overexpressing hb show no kni expression at all (Nauber et al. 1988; Rothe et al. 1989; Kraut and Levine 1991b), and embryos misexpressing hb show spatially specific repression of kni expression (Clyde et al. 2003). There is no clear posterior expansion of kni in hb mutants (Hülskamp et al. 1990; Clyde et al. 2003). This could be due to the relatively weak and late repressive contribution of Hb on the posterior kni boundary or due to partial redundancy with repression by Gt and Tll (Figure 7, E and F). The posterior hb domain expands anteriorly in kni mutants, but anterior hb expression is not altered in these embryos (Jäckle et al. 1986; Clyde et al. 2003). Nevertheless, a role of Kni in positioning the anterior hb domain is suggested by the fact that misexpression of kni leads to spatially specific repression of both anterior and posterior hb domains (Kosman and Small 1997; Wu et al. 2001; Clyde et al. 2003). Moreover, only slight posterior expansion of anterior hb is observed in Kr mutants, while hb is completely derepressed between its anterior and posterior domains in Kr kni double mutants (Clyde et al. 2003).
Repression between overlapping gap genes:
gt, kni, and Kr show repression by their immediate posterior neighbors hb, gt, and kni, respectively (Figure 4). Retraction of posterior Gt from the posterior pole during midcycle 14A fails to occur in hb mutants (Mohler et al. 1989; Eldon and Pirrotta 1991; Kraut and Levine 1991a), and no gt expression is observed in embryos overexpressing hb (Eldon and Pirrotta 1991; Kraut and Levine 1991b). The posterior kni boundary is shifted posteriorly in gt mutant embryos (Eldon and Pirrotta 1991), and kni expression is reduced in embryos overexpressing gt (Capovilla et al. 1992). Note that these effects are very subtle and were not reported in similar studies by different authors (Kraut and Levine 1991b; Rothe et al. 1994). A weak but functional interaction of Gt with kni is consistent with our results. This interaction was found to be essential even in a circuit (29007) where it was deemed below significance level (Figure 4, A and C, and data not shown). Finally, Kni has been shown to bind to the Kr regulatory region (Hoch et al. 1992), and the central Kr domain expands posteriorly in kni mutants (Jäckle et al. 1986; Gaul and Jäckle 1987).
In contrast, we have been unable to detect any effect of Kr on hb (Figure 4, A and B). However, hb expression expands posteriorly in Kr mutants (Jäckle et al. 1986; Gaul and Jäckle 1989; Clyde et al. 2003). This effect is likely to involve repression of hb by Kni. Kni levels are reduced in Kr embryos (Pankratz et al. 1989). hb is completely derepressed between its anterior and posterior domains in Kr kni double mutants, whereas anterior hb does not expand at all in kni mutants alone (Clyde et al. 2003). Taken together with our results, this suggests that there is direct repression of hb by Kr in the embryo, but it is at least partially redundant with repression of hb by Kni.
Unlike repression by posterior neighbors, we have found no or only weak repression of posterior kni, gt, and hb by their anterior neighbors Kr, kni, and gt, respectively (Figure 4). Most gap gene circuits show weak activation of hb by Gt (Figure 4, A and E). Graphical analysis failed to reveal any functional role for such activation (Figure 5, H and I). Moreover, we have found no functional interaction between gt and Kni (Figure 4, A and D). Although relatively weak repression of kni by Kr was found in 6 out of 10 circuits (Figure 4, A and C), no specific patterning defects could be detected in the other 4. Consistent with the above, expression of posterior hb is normal in gt mutants, and both the anterior boundaries of posterior gt and kni are positioned correctly in kni and Kr mutant embryos, respectively (Mohler et al. 1989; Pankratz et al. 1989; Eldon and Pirrotta 1991; Rothe et al. 1994).
Note that we have never observed activation of kni by Kr (Figure 4, A and C), which has been proposed to explain decreased expression levels of kni in Kr mutants (Pankratz et al. 1989; Rothe et al. 1994). Our results strongly support the view that this interaction is indirect through Gt, which is further corroborated by the fact that kni expression is completely restored in Kr gt double mutants compared to that in Kr mutants alone (Capovilla et al. 1992).
We have found a significant repressive effect of Hb on Kr (Figure 4, A and B). Consistent with this, Hb has been shown to bind to the Kr regulatory region (Hoch et al. 1991), and the central Kr domain expands anteriorly in hb mutants (Jäckle et al. 1986; Gaul and Jäckle 1987). However, partial redundancy of this interaction is suggested by correct positioning and shape of the anterior Kr domain in a circuit (28005) that does not show repression of Kr by Hb (Table 1).
It has been proposed that Hb plays a dual role as both activator and repressor of Kr (see Introduction). In the framework of the gene circuit model, concentration-dependent switching of regulative action could be implemented by allowing genetic interconnection parameters to switch sign at certain regulator concentration thresholds. Our current model explicitly does not include such a possibility. Nevertheless, we have been able to obtain circuits that reproduce Kr expression faithfully (Figure 3), suggesting that a dual role of Hb is not required for proper Kr expression. Moreover, we have never observed activation of Kr by Hb in any of the circuits (Figure 4, A and B). Therefore, our results support a mechanism in which the activation of Kr by Hb is indirect through derepression of kni.
Repression by Tll:
Only a few earlier theoretical approaches have considered terminal gap genes (Meinhardt 1986; Tchuraev and Galimzyanov 2001). Gap gene circuits accurately reproduce tll expression (data not shown). However, in gene circuits, tll is subject to regulation by other gap genes, which is inconsistent with experimental evidence (Brönner and Jäckle 1991). In contrast, the correct expression pattern of tll in gap gene circuits allows us to study its effect on other gap genes in great detail. We have found strong repressive effects of Tll on Kr, kni, and gt (Figure 4). Tll binding sites have been found in the regulatory regions of Kr (Hoch et al. 1992) and kni (Pankratz et al. 1992). In tll mutants, Kr expression is normal (Gaul and Jäckle 1987; Reinitz and Levine 1990), whereas expression of kni expands posteriorly (Pankratz et al. 1989), and the posterior gt domain fails to retract from the posterior pole (Eldon and Pirrotta 1991; Kraut and Levine 1991a). No expression of Kr, kni, or gt can be detected in embryos overexpressing tll under a heat-shock promoter (Kraut and Levine 1991b; Steingrimsson et al. 1991).
Comparison to logical analysis:
The logical analysis by Sanchez and Thieffry (2001) is the only other theoretical study of the gap gene system that achieves a level of detail comparable to the analysis presented here. Our results largely agree with Sanchez and Thieffry (2001), with the following notable exceptions. tll and the posterior hb domain were not considered in the logical analysis. The absence of posterior hb could explain why Sanchez and Thieffry (2001) did not report a repressive feedback loop between hb and kni, which we have found to be essential for gap gene regulation. Moreover, a difficulty with logical analysis is that functional thresholds must be assigned to continuous protein concentrations prior to the analysis. This leads to assigning functional borders of expression domains in the embryo, which may not coincide with observable expression borders. In the case of Sanchez and Thieffry (2001), a priori assignment of thresholds implicitly results in the posterior borders of the anterior hb domain, central Kr domain, and central kni domain coinciding, while the posterior domains of kni and gt show no overlap (cf. Figure 1, I and L). The authors conclude that a dual role of Hb in Kr regulation is required to account for the large overlap between the two respective expression domains. Our expression data indicate that the posterior hb boundary (Figure 1, I and L) lies in the middle of the Kr domain, and our analysis suggests that a dual role of Hb is not required for correct expression of Kr. Finally, the discrete logical approach failed to reveal the role of autoactivation in sharpening gap domain boundaries during cycle 14A. The thresholds selected by Sanchez and Thieffry (2001) divided the embryo into four discrete zones along the A-P axis, but modeling boundary sharpening requires an approach with a larger spatial resolution.
Limitations of the model:
We observe artificially high levels of gap proteins during early cycle 13 (data not shown) and earlier cleavage cycles if included in the model (Reinitz et al. 1995, and data not shown). This is a serious problem for analysis of early gap gene regulation, since premature accumulation of gap proteins causes premature gap-gap regulatory interactions that rapidly dominate early inputs from maternal genes. In the embryo, production delays between the time when a transcription factor binds to a regulatory region and the completion of subsequent protein synthesis have a significant influence on the timing of gene expression (Rothe et al. 1992). Cleavage cycles 10–12 are only ∼7–13 min long, which is significantly shorter than cycles 13 and 14A (Foe and Alberts 1983). A production delay on a scale of 5–15 min combined with transcript degradation during mitosis (Shermoen and O'Farrell 1991) can account for the absence of zygotic gap proteins before cycle 13. Therefore, production delays will have to be incorporated into gap gene circuits to obtain correct early gap gene expression and regulation.
Gene circuits can be used for prediction of expression patterns in mutants (Sharp and Reinitz 1998). Mild changes in genotype, such as varying Bcd dosage, led to successful prediction of mutant gap gene expression patterns using gap gene circuits (Simpson-Brose et al. 1994; Reinitz et al. 1995). In contrast, we have not been able to predict gap gene expression patterns in null mutants. This could be due to spurious early gap gene regulation (see above). Alternatively, it might be due to scaling indeterminacy in our quantitative expression data. We currently do not know the proportionality constant, different for each protein, that relates fluorescence levels with absolute protein concentrations. Just as improvements in the data used in the present study improved results over previous studies, we expect that further improvements in data quantification will lead to further improvement in the predictive capacity of our models.
Our analysis yields a much more dynamic picture of gap gene expression than previously thought. During the late blastoderm stage, gap gene expression patterns and their regulatory interactions change on a very rapid scale. Many open questions remain about how or if the transient and highly dynamic nature of these patterns and interactions affects the establishment of the stable segmentation prepattern of segment polarity gene expression. To address these questions, future gene circuit models will have to include more downstream layers of the segmentation gene network, namely pair-rule and segment-polarity genes. Just as the gap gene system is only the first step in the regulatory hierarchy of the segmentation gene system, our current models are only the first step toward more comprehensive gene circuit models of segment determination in D. melanogaster.
We thank Hilde Janssens, Nick Monk, Julia Dallman, and James D. Baker for discussions and comments on the manuscript. This work was supported by the National Institutes of Health, grants RR07801 and TW01147, and by Civilian Research and Development Foundation Grant Assistance Program Award RBO-1286.
Communicating editor: T. Schüpbach
- Received February 6, 2004.
- Accepted April 20, 2004.
- Genetics Society of America