Canalization of gene expression and domain shifts in the Drosophila blastoderm by dynamical attractors

Abstract

The variation in the expression patterns of the gap genes in the blastoderm of the fruit fly Drosophila melanogaster reduces over time as a result of cross regulation between these genes, a fact that we have demonstrated in an accompanying article in PLoS Biology (see Manu et al., doi:10.1371/journal.pbio.1000049). This biologically essential process is an example of the phenomenon known as canalization. It has been suggested that the developmental trajectory of a wild-type organism is inherently stable, and that canalization is a manifestation of this property. Although the role of gap genes in the canalization process was established by correctly predicting the response of the system to particular perturbations, the stability of the developmental trajectory remains to be investigated. For many years, it has been speculated that stability against perturbations during development can be described by dynamical systems having attracting sets that drive reductions of volume in phase space. In this paper, we show that both the reduction in variability of gap gene expression as well as shifts in the position of posterior gap gene domains are the result of the actions of attractors in the gap gene dynamical system. Two biologically distinct dynamical regions exist in the early embryo, separated by a bifurcation at 53% egg length. In the anterior region, reduction in variation occurs because of stability induced by point attractors, while in the posterior, the stability of the developmental trajectory arises from a one-dimensional attracting manifold. This manifold also controls a previously characterized anterior shift of posterior region gap domains. Our analysis shows that the complex phenomena of canalization and pattern formation in the Drosophila blastoderm can be understood in terms of the qualitative features of the dynamical system. The result confirms the idea that attractors are important for developmental stability and shows a richer variety of dynamical attractors in developmental systems than has been previously recognized.

Figure 1. Reduction of variation in segmentation gene expression patterns over time.

Kr and gt expression patterns 47 min (A, early) and 3 min (B, late) before gastrulation. ftz expression patterns 34 min (C, early) and 3 min (D, late) before gastrulation. The expression patterns shown here and subsequently are from the middle 10% of dorsoventral positional values; 0% egg length (EL) is the anterior pole. The standard deviation () of expression level at the peak of the central Kr domain is 33.0 early ( embryos), and reduces to 17.4 late (). The position of the domain peak has early and late. The expression level at the peak of the gt posterior domain has early (), and late (). The position of the domain peak has early and late. The expression pattern of ftz has extensive qualitative variation early (panel C, ), but is very reproducible just before gastrulation (panel D, ). Embryos with the most diverse patterns of ftz were chosen for panel C, while they were randomly chosen for all other panels.

Figure 2. Gap gene data and expression patterns in the gene circuit chosen for analysis.

(A,B) Time class T8 (∼3 mins before gastrulation; Table S1) embryos immunostained for Kr and Gt (A) and Hb and Kni (B). Anterior is to the left and dorsal is above. Bars indicate the modeled region. (C) Data for the maternal protein gradients Bcd (cleavage cycle 13), Hb (cycle 12), and Cad (cycle 12). (D) Average gap gene data in time class T8. (E,F) Cleavage cycle 13 (E) and time class T8 (F) gap gene expression patterns produced by the gene circuit. The arrow shows the main patterning defect, which is related to experimental noise in Tll data . (G,H) Gap gene expression patterns produced by the same circuit in cleavage cycle C13 (G) and time class T8 (H) in the absence of diffusion ( for all proteins). (I,J) Expression patterns produced by the circuit in cleavage cycle 13 (I) and time class T8 (J) in the absence of diffusion and tll. The dashed vertical line shows the region (35%–71% EL) in which the expression patterns of the circuit excluding tll (J) agree with the circuit that has tll (H). The anterior and posterior regions identified in the stability analysis (Mechanisms of Canalization and Pattern Formation section) are highlighted in panel J in blue and red respectively. (K) The topology of the gap gene network determined by the gene circuit method.

Figure 3. Equilibria determined by the continuous analysis.

The points are equilibria calculated by Newton-Raphson at discrete nuclear positions using Cad data (red curve in Fig. S3), while the solid lines are equilibria determined by the continuous analysis using the interpolated Cad profile (black curve in Fig. S3). Point attractors are blue or cyan and saddle equilibria having one or two eigenvalues with positive real part are red or brown respectively. The -axis is the projection of equilibria positions on the Kr axis (A and B) or the Gt axis (B). The -axis is the bifurcation parameter, the A–P position . (A) The bifurcations observed in both the analysis at discrete positions and the continuous analysis are encircled. The bifurcations only seen in the continuous analysis are highlighted in boxes. The arrow highlights the annihilation between at 53% EL that divides the anterior and posterior regions. See Table S4 for bifurcation values of the A–P position. (B) The point attractor (right -axis), showing its continuous movement from the hb,gt-on to the hb-on state. is plotted on the left -axis.

Figure 4. Two distinct dynamical regimes that control canalization.

(A,B) Three-dimensional projections of four-dimensional phase portraits. Point attractors are blue and saddle equilibria having one eigenvalue with positive real part are red. Time is represented as a color gradient along the trajectories, with start of cycle 13 as green, and gastrulation as red; trajectories are blue after gastrulation. The midpoints of time classes (Table S1) T1, T3, T5, and T7 are indicated with cyan points. The sharp bends in trajectories are mitoses, when synthesis shuts down. See caption of Fig. S5 for details of rendering the phase portraits. (A) Hb-Kr-Gt projection of the phase portrait at 37% EL, highlighting the anterior dynamical regime. 10 trajectories are plotted. (B) Hb-Kr-Kni projection of the phase portrait at 57% EL, highlighting the posterior dynamical regime. 25 trajectories are shown. (C,D) Reduction of initial variation. (C) The time evolution of the volume of the box (52,100)×(0,1)×(0,60)×(0,1) representing initial variation in the anterior region nucleus at 37% EL. Dashed line is volume in log-scale; it reduces by a factor ∼10⁸ by gastrulation. (solid line) gives the average shrinkage in a dimension. By gastrulation, each dimension shrinks by a factor of ∼20. are shown normalized to 1. (D) The time evolution of the volume of the box (0,20)×(0,80)×(0,80)×(0,80) representing initial variation in the posterior region nucleus at 57% EL. reduces by a factor of ∼10⁶ by gastrulation, and each dimension shrinks by an average factor of ∼10. (E) Kr-Gt-Kni projection of the phase portrait at 57% EL showing that the manifold traverses the gap gene states in the posterior region, Kr-on, kni-on, and gt-on.

Figure 5. Maternal Hb is the morphogen in the posterior region.

Gap gene concentrations at the midpoint of T8 in a single nucleus as a function of initial Hb concentration. The initial Hb concentration was varied uniformly in the nucleus at 63% EL (peak of the kni domain at gastrulation), keeping the Bcd and Cad inputs constant. The nucleus produces all gap gene states in the posterior region from the Kr peak (53% EL) to the gt peak (71% EL) as initial Hb concentration is decreased from 40 to 0. The shapes of the “domains” are distorted since maternal Hb has faster than linear decay with position; as a consequence anterior “domains” are exaggerated here. The -axis on the top shows A–P positions determined from the values of maternal Hb, showing that the domains are in correct proportions spatially. Posterior region nuclei form domains by responding to maternal Hb without any instruction from Bcd.

Figure 6. Shifts due to attraction by the manifold .

(B,D) Dynamics of protein concentrations in the nucleus at 59% EL, through which the Kr posterior and Kni anterior boundaries pass as they shift to the anterior. The nucleus at 59% EL is indicated with a dashed vertical line. (A) The Hb-Kr-Kni projection of the phase portrait of the nucleus. The invariant manifold is shown as a magenta tube. The trajectory in the nucleus is plotted in a continuous color gradient from green (t = 0 min) to red (t = 71.1 min, gastrulation). Times after gastrulation are depicted as blue. The nucleus passes through intermediate states (indicated with an arrow) with high Kr concentrations before reaching a state with high Kni concentration by the onset of gastrulation. This registers as an anterior shift in the posterior Kr and anterior kni borders. (C) A two dimensional projection of the trajectory in the Kr, Kni plane. The trajectory (red) starts at the origin. It attains a high Kr value at t = 45 min (arrow) before approaching high Kni values. The temporary reversal in the trajectory is a mitosis, during which the trajectory moves toward the origin. Time after gastrulation is shown in blue.

Figure 7. Summary of the dynamical mechanisms of canalization and pattern formation.

(A) Maternal gradients in the analyzed region 35%–71% EL. (B) Dynamical mechanisms. 2D phase portraits and trajectories of highlighted nuclei are shown. Dotted lines connect the highlighted nuclei and their phase portraits with the gap gene state shown in panel C. (C) Gap gene expression patterns (T8) in the absence of diffusion and tll. (D) Dynamical mechanisms of canalization.