Modeling cell intercalation during Drosophila germband extension

Our kinetic model is based on experimental observations of the biochemical pathways leading to PCP [9, 13, 14, 18], the mechanisms of AJ remodeling [9, 11, 13, 19] and a tension-myosin positive feedback [10, 16, 24]. Figure 3 depicts these interactions schematically. By tracking the amount of Shroom (S), Rok (R), Myo-II (m) and Baz (B) on each cell edge, the model predicts polarization of those proteins and in turn an anisotropic contractile force on AP edges. On each edge, the evolution of these species are governed by the following differential equations:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \frac{dR}{dt} = q_r - k_r e^{-S} R, \label{eq:drdt} \nonumber \end{align} \tag{ 1 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \frac{dB}{dt} = q_b - k_bRB, \label{eq:dbdt} \nonumber \end{align} \tag{ 2 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \frac{dm}{dt} = k_1R - k_2 mB - k^-m, \label{eq:dmdt} \nonumber \end{align} \tag{ 3 }$$

where for S we impose an anisotropic, time-independent level on the AP, DV and slanted ('shoulder') edges [14, 24]. The terms and the parameters are as explained in [24], with each unit of S, m, R and B consisting of 162 molecules of the respective protein. In particular, the e^−S form in equation (1) is chosen for mathematical convenience, the real kinetics of the Shroom-Rok chemical pathway being unclear at present. Note that $k^- = k_3 e^{-k_4 f }$ is a myosin detachment rate that is suppressed by tensile force f on the edge, a well-known positive feedback between junctional myosin and tension [10].

Zoom In Zoom Out Reset image size

We model the pulsating medial actomyosin network by imposing a sinusoidal myosin signal on the centroid of each cell:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \label{eq:myosin_spokes} m_m = A_m\sin^2 \left(\frac{t\pi}{T} + \phi \right). \nonumber \end{align} \tag{ 4 }$$

Although the cyclic contraction of medial actomyosin is an interesting phenomenon in its own right [15], it is not a central element in the current model, and is therefore modeled ad hoc. The phase angle ϕ is specified for each cell in the tissue (see figure 2) such that the neighbors exhibit anti-correlation on average, in agreement with experimental observations [15]. See the online supplemental information (SI), section S2 for more details on the neighbor-neighbor correlation for medial myosin.

As in previous vertex models, the cells deform and move through the over-damped motion of its vertices:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \eta\frac{d \boldsymbol{x}_i}{d t} =\sum_j f_{ij}\frac{\boldsymbol{x}_j - \boldsymbol{x}_i}{|\boldsymbol{x}_j - \boldsymbol{x}_i|} + \sum_k \boldsymbol{f}_{ik} + \sum_k \boldsymbol {p}_{ik}, \label{eq:xi} \nonumber \end{align} \tag{ 5 }$$

where η is an effective viscosity, $\boldsymbol{x}_i$ is the position of the ith vertex, f_ij is the total force on the edge connecting vertices i and j, $\boldsymbol{f}_{ik}$ is the medial myosin force from neighboring cell k pointing toward its centroid, and $\boldsymbol{p}_{ik}$ is a pressure force in cell k along the outward 'normal direction', which bisects the polygonal corner at i. Each cell edge exhibits a passive elastic force and a myosin-based active contractile force:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle f_{ij} = \mu (l_{ij} - l_0) + \beta m_{ij} , \label{eq:fij} \nonumber \end{align} \tag{ 6 }$$

where μ is the elastic modulus of the edge, l_ij and l₀ being its current and resting length, m_ij is the Myo-II level on the edge, and β is the myosin force constant [24]. The magnitude of the medial myosin force is simply $|\boldsymbol{\,f}_{ik}| =\beta m_{mk}$, m_mk being the medial myosin of cell k, and that of the pressure force is $|\boldsymbol{p}_{ik}|=\alpha (A_0-A_k)$, A₀ and A_k being the resting and current cell area of neighbor k, and α the bulk modulus for the cell.

Besides mechanical forcing, remodelling of AJs is another factor in the contraction of AP edges and elongation of DV edges. The exact mechanisms for the remodeling are not clear at present, but several plausible pathways have been suggested, involving Rok, Abl, Baz, β-catenin and E-cad [11, 13, 18]. Among these, only Rok and Baz are accounted for in our model, and they localize to complementary cell borders upon onset of PCP [24]. We have chosen Baz to be the variable through which to encode AJ remodeling, as Baz exhibits similar polarization patterns to β-catenin and serves as a good surrogate for the latter [13]. In turn, the strength of AJs is modeled by passive elasticity on each edge, with a resting length l₀(B) that varies with Baz. On AP edges, therefore, l₀ shortens with declining Baz and thus facilitates contraction of the edge, while on DV edges, l₀ grows with Baz to facilitate DV elongation. More details can be found in section S1 of the online SI.

Finally, the geometric, kinetic and mechanical parameters are evaluated from experimental data where possible. These are tabulated in SI (section S3) with sources and justifications. The numerical techniques for the computation are also discussed in SI.

The 52-cell assembly of figure 2 represents only a portion of the real germband. Unlike in whole-embryo models [28], we need to treat the outer boundaries of the model tissue so as to mimic the interaction with the surrounding tissue in reality. As the germband does not exhibit spatial periodicity, it would be erroneous to impose periodic boundary conditions. Modeling the surrounding tissue by elastic springs attached to a fixed circular boundary produced poor results in Lan et al [24]. The springs hindered the contraction along the DV axis and extension in the AP axis and skewed the convergent extension of the tissue. In reality, the surrounding tissue would be undergoing convergent extension as well, and its deformation would accommodate that of the model tissue. To represent this interaction, we have introduced 'cables' emanating from the boundary nodes more or less isotropically (figure 2), each carrying a constant tension (to be defined in the next subsection). As the boundary vertices move, the cables shift as well, without changing their orientation.

In addition, PMI exerts a stretching force on the germband along the AP axis. Lye et al [5] quantified this stretching not as external forces on the germband, but as shape change of the cells. Collinet et al [6] represented it as movement of the boundary along the AP axis. Following [5], we view the germband extension as comprising the effects of cell intercalation and cell-shape change due to PMI. Using their measurements, we convert the rate of cell-shape elongation into a rate of AP extension for the tissue, and add that, at each time step, onto the movement of the boundary vertices that results from the intercalary forces and external cables. Thus, the global stretching due to PMI amounts to an additional horizontal outward movement on our boundary vertices. Details can be found in SI (section S4).

It is convenient to first imagine an equilibrium state in which the tissue is in chemical and mechanical balance. Imposing a uniform Shroom level of S_e  =  0.5 on all edges, we find a chemical equilibrium with R_e  =  0.659, B_e  =  0.607 and m_e  =  0.912 on all edges. Force balance is achieved by having perfectly uniform hexagonal cells of edge length $l_e=l^{eq}_0 + \tau_{\rm pre}/\mu - \beta m_e/\mu$, with each edge carrying a uniform pre-tension $\tau_{\rm pre}$, $l^{eq}_0=l_0 (B_e)$ being the resting length at equilibrium level of Baz B_e. On the boundary vertices, the pre-tension is balanced by the external cables. See SI (sections S1, S3 and S4) for detailed discussion on choosing the values of l₀, $\tau_{\rm pre}$ and the cable forces F_cable to ensure mechanical equilibrium prior to the onset of cell intercalation.

The procedure for simulation is as follows. For the initial period, up to t  =  400 s, we maintain the uniform S_e  =  0.5, R_e  =  0.659 and B_e  =  0.607 on all edges. The medial myosin pulse is applied with phase lags among neighboring cells as prescribed in figure 2. Thus, the cells oscillate gently during this period, with no cumulative changes or polarization. Because of the positive feedback of junctional tension on myosin, the junctional myosin will oscillate about a mean value slightly below m_e thanks to the additional compression on the cell edge by the medial myosin.

At t  =  400 s, we initiate the PCP by randomly selecting an AP edge by the stochastic scheme described in the next subsection, and raising S to 1 on it to initiate polarization and AP contraction. In the mean time, we assign an intermediate S  =  0.625 on shoulders connected to the contracting AP edge. As new DV edges are created, from T1 transition or rosette resolution, S  =  0.5 is assigned on them. This anisotropic Shroom pattern is based on experimental observations [14], and is assumed to be time-independent throughout the simulation.

During germband extension, the contraction of AP edges appears at apparently random times and locations. Here we model temporal stochasticity as a Poisson process. The frequency parameter is estimated from experimental data [13] to be roughly 0.05 event of AP contraction per cell per minute (see details in SI, section S5). Note that the Poisson distribution thus designed yields the expected frequency of AP contraction only at the start of our simulation. The geometric regularity of our model tissue, with relatively few AP edges, means that the actual rate of AP contraction drops in time and fewer T1's and rosettes will have taken place in the end than indicated by the target frequency.

Spatial stochasticity is implemented by selecting a random available AP edge for each contraction. The term available refers to AP edges that are connected to at least 3 cells and have yet to begin contraction. An interior AP edge, with neither end vertex exposed to the outside, are connected to 4 cells, while those at the boundaries are connected to 3 or fewer cells. T1 swap is impossible with AP edges that are connected to 2 cells or 1 cell. We also impose a constraint that no boundary AP edge be selected if an interior AP edge is available. The reason is that a sequence of early contractions on the boundaries tend to make some cells concave and cause that part of the tissue to collapse. AP contraction typically culminates in neighbor swaps through T1 transition or rosette formation, to be discussed next.

We have adapted the scheme for T1 transition of Lan et al [24]. When an AP edge contracts to a critical length of $l_c = 0.02 l_e$, i.e. 2% of its equilibrium length, a neighbor swap is implemented by rotating the AP edge counter-clockwise into a new horizontal DV edge and reconnecting its vertices to the closest neighbors. Shroom is set to S  =  0.5 on the nascent DV edge, in keeping with the observed anisotropy [14]. Different from Lan et al [24], however, we do not set m  =  0, R  =  0 and B  =  0 on the nascent DV edge. Instead, the DV edge carries over all the protein levels from before the rotation. This is because the model of Lan et al [24] lacked any treatment of the establishment of cell-cell junctions along the growing DV edge. As our model captures this effect via l₀(B), we no longer need to make these drastic assumptions.

Rosettes form after two adjacent AP edges and the shoulder edge in between contract to a single vertex (figure S3). We model the initiation of rosettes as a stochastic process. If, by chance, two AP edges connected to the same shoulder edge are selected successively to contract, we use a second Poisson distribution to determine whether a rosette should form. If an affirmative answer is returned, Shroom is elevated to S_m  =  1.4 on the two AP edges as well as on the shoulder sandwiched in between such that the 3 edges form a 'linked AP edge' [9, 10]. Its contraction to a single vertex then produces a six-cell rosette. To our knowledge, the only previous vertex model that accounted for rosettes is Trichas et al [29], for the mouse visceral endoderm. In experiments, rosettes made of 5–11 cells have been reported [9]. For simplicity, however, we restrict our model to six-cell rosettes only. The frequency parameter in the second Poisson process is estimated from in vivo data of Tamada et al [13], at 0.026 rosette per cell per minute (see details in SI, section S5), about half the rate of formation for T1 transitions. An additional constraint is that a cell is barred from participating in two or more rosettes at the same time. This avoids the appearance of overlapping or neighboring rosettes, which the model fails to resolve properly owing to its geometric simplicity. Note that in vivo the germband cells undergo multiple rounds of rosette formation, and neighboring and overlapping rosettes do occur [9].

Rosette resolution follows a similar scheme to the T1 transition described above. Once the linked AP edge has shrunk to a critical length of $l_c = 0.06 l_e$, we replace it by a new horizontal DV edge, with neighboring vertices separated by a uniform distance of 0.02l_e. The critical length is three times that for T1 transition since 3 edges are involved in the rosette formation. After the transformation, the cell edges are reconnected onto the new DV edge as explained in SI, section S6, so as to ensure a unique topological outcome. On these new DV edges we set the Shroom level to S_e  =  0.5 according to PCP [14], and let Rok, Baz and Myo-II evolve according to the kinetic model.

We start by simulating the intercalation in the 52-cell tissue that undergoes T1 transition only, without rosettes. Because of the medial myosin pulse (equation (4)), the cells exhibit a gentle oscillation for the initial 400 s of the simulation, before PCP and intercalation set in. This period corresponds to the early Stage 6 of experimental observations (e.g. figure 1 in Blankenship et al [9]), and is similar to what Lan et al [24] described using the four-cell assembly.

Starting at t  =  400 s, we randomly select AP edges to start PCP, at the frequency mentioned above, until all available AP edges have undergone contraction and T1 transition. The whole process is shown in movie 1 of the SI. Most features of the T1 transition are similar to those predicted by Lan et al [24], and some details are given in section S7 of the SI. For example, the development of PCP produces a greater contractile force on the AP edge (figure S4). The AP contraction is driven mainly by the junctional myosin force against resistance of the shoulder forces, while DV elongation is driven by the medial myosin forces in the two cells anterior and posterior to the new DV edge as well as the shoulder forces, mostly against the myosin force on the DV edge (figure S5). A new feature is the disassembly of the AJs on the contracting AP edges and the strengthening of AJs on the growing DV edge. Modeled through the resting length l₀(B) of the cell edges, the declining elastic resistance on the AP edge promotes its contraction, while the lengthening l₀ on the DV edge relieves the elastic tension and facilitates DV elongation.

Over the course of the simulation, T1 transitions occur stochastically in space and time. To probe this stochasticity, we have run 10 independent simulations under identical conditions. Because of the regularity of the geometric setup, a total of 49 T1's take place in each run, and the tissue elongates to the same regular pattern (see movie 1). An interesting indicator of the stochasticity of the process is the times for AP contraction T_AP and DV elongation $T_{DV}$. For each AP contraction, we measure the time between the onset of PCP and the moment that the AP edge length shrinks to the threshold of 0.02l_e. The duration of DV elongation is measured between the appearance of the new DV edge and the moment when the DV length first attains the equilibrium length of l_e  =  4 μm. This treatment differs from that of Lan et al [24], whose DV edge length oscillated more or less regularly, and its third peak was taken to be the end of the elongation process. In the current model, the elongation of the DV edge is less regular and a third peak cannot be easily defined. Compiling T_AP and $T_{DV}$ data for all the T1 transitions in the 10 simulations, figure 4 shows that both times exhibit a rather wide distribution. This is because AP contraction and DV elongation are affected by what the neighboring cells happen to be doing, and the latter varies stochastically thanks to the random selection of AP edges for contraction. After averaging all 490 T1's, we have mean values T_AP  =  493 s and $T_{DV}=304$ s. Experiments have recorded T_AP ranging from 450 to 600 s and $T_{DV}$ ranging from 420 to 600 s [8, 9]. The predicted mean values are somewhat shorter than the measurements. Potential causes include the uncertainties in evaluating the parameters and in defining the start and end points of contraction and elongation. Besides, the broader distribution of figure 4 suggests that our model tissue is too small and sensitive to noise.

Zoom In Zoom Out Reset image size

Next, we will examine the degree of anisotropy in terms of the distribution of the proteins and the contractile forces. Following prior experimental and modeling work [10, 14, 15, 18, 22, 24, 30], we focus on the ratio of these quantities between the AP and DV edges. In experiments, the level of various proteins are measured by averaging vertical or horizontal edges of differing lengths [21, 30]. These edges are presumably at different stages of T1 transition. To approximate the averaging in the experiments, we have taken the time average on AP edges from the onset of AP contraction to the T1 swap, and on DV edges from the T1 swap to the completion of DV elongation. Additionally, we have averaged among the T1 transitions for the 10 realizations. The results are compared with prior experiments and modeling in table 1.

Table 1. AP/DV ratios of various proteins and the edge force during T1 transition: comparison between model predictions and experimental measurements.

	Experiments	Model of Lan et al [24]	This model
Rok	1.5–1.9 [14]	1.65	$1.38\pm 0.08$
Baz	0.33–0.5 [14, 18, 31]	0.607	$0.81\pm 0.07$
Myo-II	1.3–2.1 [14, 18, 27, 31]	2.74	$1.80\pm 0.40$
Edge force	1.7–2.55 [10, 30]	7.5	$1.60\pm 0.70$

The AP/DV contrasts in Rok and Baz levels are somewhat weaker in our model than in the experiments, while the Myo-II contrast falls in the experimental range. Our force ratio is slightly below the experimental data, but is a marked improvement over the prediction of Lan et al [24]. Experimentally, the force ratio is estimated from the ratio of end retraction speeds after laser ablation of the cell edges. What is measured, therefore, is the total force on the edges, which in our model comprises the myosin and elastic contributions. Since Lan et al [24] did not consider the myosin tension on the nascent DV edges, they have over-predicted the ratio by a wide margin. In addition, a more precise comparison can be made of the magnitude of the edge forces. Bambardekar et al [30] have reported an average force of 0.324 nN on AP edges and 0.127 nN on DV edges. In comparison, our model predicts 0.68 nN and 0.46 nN, respectively, more than twice as high as the measured values. In the experiment, the force values were averaged over edges oriented within 45° of the AP or DV axis, and may have included various shoulder edges. In our model, the AP and DV edges are clearly defined. This may partially explain the discrepancy.

Finally, let us turn to the overall deformation of the tissue. Following Collinet et al [6], we define a relative length L_r as the current AP length of the tissue, measured between the centroids of the leftmost and rightmost cells, divided by the initial length before the onset of axis elongation. Figure 5 shows the temporal evolution of L_r averaged over the 10 runs. The relative length of the model tissue (large box) rises to a value of L_r  =  1.74 at the end of the simulation, roughly 40 min after the onset of intercalation. This compares with experimental measurements of L_r  =  2–2.5 [1, 6]. In vivo, the germband undergoes many rounds of T1 and rosette processes. In our model tissue, the simple geometry allows fewer T1 transitions, and no rosettes in this case. These have contributed to the under-prediction of tissue elongation.

Zoom In Zoom Out Reset image size

By comparing the relative deformation of cell clusters of different sizes, Collinet et al [6] observed that 'tissue extension is spatially homogeneous and scale-independent'. That is to say, the strain is uniform throughout the tissue. To test for spatial homogeneity in our model tissue, we compare the relative length L_r of the three tissue clusters delineated in the inset of figure 5. The relative length agrees to within 4% at the end of the simulation among the three boxes. The earlier differences are mainly due to spatial variations of where a T1 happens first. Thus, our model has successfully recapitulated the spatial homogeneity that Collinet et al [6] demonstrated in vivo. Note further that the larger boxes show slightly larger elongation in figure 5, a feature that is also present in the experimental data (figure 1(e) of [6]). This may be related to the stretching along the AP axis due to midgut invagination, which, in our model at least, favors greater deformation of cells at the left and right edge of the tissue, as is evident from images near the end of movie 1. Interestingly, the effect of PMI on cell shape change is actually inhomogeneous in vivo [5]; cell deformation is more pronounced closer to the posterior tip of the embryo. As explained in section S4 of the SI, this feature is not included by our model.

Now we allow both T1 and rosettes in the simulation. In this subsection, we examine closely the formation and resolution of the first rosette as a representative of all rosettes. The purpose is to discuss the various mechanisms that can affect the rosette formation and resolution, and explore its local features. The global features of the full simulation, with both T1 and rosettes appearing and resolving with spatial and temporal stochasticity, will be discussed in the next subsection.

At the end of the initial oscillation, t  =  400 s, we start to select randomly AP edges to initiate PCP and contraction, as described in the previous subsection. In the representative simulation that we will describe now, the first rosette appears when two nearby AP edges in the middle of the tissue, edges 19–41 and 42–59 in figure 2, are selected at t  =  499 s and 501 s for contraction, respectively. Thus, S is raised to 1.4 on them and to 0.625 on the shoulders connected to them on the outside. At t  =  538 s, the sandwiched shoulder (41–42) is selected for rosette formation by a Poisson probability and S is raised to 1.4 on the shoulder as well. This causes contraction of the three connected edges, referred to as the linked AP edge following [9, 10] (t  =  650 s in figure 6(a)). Note that when the two AP edges start to contract, the shoulder edge in-between is first elongated. With the development of its own Rok and myosin, however, it contracts as well until the linked AP edge shrinks down to the critical length of $l_c=0.06 l_e$ at t  =  1158 s. After rotating the linked AP edge into the horizontal direction and reconnecting the vertices according to figure S3, the three nascent edges elongate (figure 6(b)). Interestingly, it is the middle edge that elongates most rapidly. The two horizontal edges elongate more slowly, and DV elongation is complete when the shortest of the three reaches l_e  =  4 μm at t  =  1882 s. The forces driving and resisting the rosette transformation are the same ones as in T1 transition, but larger in magnitude (see section S7 and figure S6 of SI). Note that the topology of the cells after rosette resolution is identical to what would have resulted from two separate T1 transitions involving the two AP edges initially chosen to contract. The rosette formation and resolution are also shown in movie 2 in SI. During this process, the surrounding tissue undergoes stochastic T1 transition, and the second rosette has started at t  =  550 s involving edges 34–51, 51–50, and 50–67.

Zoom In Zoom Out Reset image size

To effect the remarkable contraction of the three linked edges into a vertex, our model has adopted the ad hoc device of raising the Shroom level on them. Experimental observations in vivo have identified several factors that may have contributed to this contraction. First, tension is known to suppress the detachment rate of myosin motors on actin filaments [10], and this amounts to a positive feedback between tension and myosin on a cell border. This factor has been incorporated into our model. Second, the protein Abl tends to intensify β-catenin turnover on AP edges, causing disassembly of cell junctions along the linked AP edge and facilitating its contraction [13]. This has been modeled via the resting length l₀ of the linked AP edges, which decreases as Baz becomes depleted. Without this effect, AP contraction would not be able to proceed to the point of producing a T1 transition or a rosette. The third is the alignment of several cell borders into a linked edge oriented perpendicular to the AP axis [9, 10]. This aligns the contractile forces on the connected segments so as to reinforce each other's effect. Experimental images and videos suggest that the fluid-like deformation of the germband and the irregular shapes of the cells frequently bring the edges to near alignment. In our model tissue, however, the regular hexagonal shape is such that it is difficult to align the shoulder edge within $30^\circ$ of the DV axis. As a result, this factor is not well represented in our model. Fourth and finally, recent experiments suggest that stretching of a cell edge to a large strain or strain rate may predispose the edge to a strong surge of myosin recruitment afterwards [10, 32, 33]. We have tested this factor in our model by making the rates of myosin attachment and retention functions of the strain or strain rate (see SI, section S8 for more details). This at most provides an initial boost for myosin that disappears later with the AP contraction, as the strain and strain rates both decline. In reality, the reaction of a stretched cell edge may possess a certain memory or relaxation time. This remains to be explored using more sophisticated mathematical formulations.

Is the rosette more effective than T1 at elongating the tissue in the AP direction? Zallen and Blakenship [11] reported that 'rosette formation amplifies the effects of neighbor exchange, producing a 5.5-fold average change in the AP:DV aspect ratio, compared to a 2.5-fold average change for T1 processes'. When Collinet et al [6] measured AP extension associated with T1 and rosettes, they 'did not find a significant difference (figures 1g, h)'. To probe this question using our model, we first note the different ways of measuring tissue elongation in the literature. Some researchers [9, 11, 24] use the change in the AP:DV aspect ratio, while others [6] use the relative length, as in figure 5 of this paper. As the extensional strain along the AP axis and the compressive strain along the DV axis have roughly the same magnitude [6, figure 1k], one can compare the relative length with the square root of the change in AP:DV aspect ratio. Besides, measuring the relative length from the outline of the cell cluster yields a smaller ratio than using the centroid of the outmost cells (see figure S1 of [6]), especially for smaller clusters. With these clarified, we can compare our model predictions with the experimental data.

Using the cell outlines, we have measured the evolving relative length L_r of four-cell clusters during T1 and six-cell clusters during rosette transformation. After shifting the start of such events to t  =  0, we compute an average L_r for the T1's and for the rosettes, and compare them with the measurements of Collinet et al [6] in figure 7. On average, the T1 clusters achieve a relative length of $L^{\rm T1}_r=1.39$ at the end, whereas the rosettes attains a somewhat stronger stretching at $L^{\rm ro}_r=1.48$. The data of Collinet et al [6] show little difference between the two structures, both attaining a relative length of L_r  =  1.44 (figure 1(h) therein). Thus, our model indicates a moderate advantage to the rosettes in terms of tissue elongation.

Zoom In Zoom Out Reset image size

How to reconcile the above with the apparently conflicting trends in experimental observations? Blakenship et al [9] offered a clue when they noted 'a change in the AP/DV aspect ratio ... by a factor of 2.1 for rosettes of 6–7 cells (n  =  74) and by a factor of 3.0 for rosettes of 8–10 cells (n  =  9)'. The factor of 2.1 translates to a relative length of $L_r=\sqrt{2.1}=1.45$, close to that of our six-cell rosette and that of Collinet et al [6], whose figure 1(g) shows rosettes of 6 cells. Therefore, it appears that rosettes involving a larger number of cells may produce significantly more stretching than T1's, while six-cell rosettes are only moderately more effective. Incidentally, the larger numbers reported by Zallen and Blakenship [11] are based on cell centroid positions, as opposed to the cell outlines. If one uses their figure 3 to estimate the conversion factor between these two measures, then their centroid-based data—5.5-fold for rosettes and 2.5-fold for T1—translate to 2.94 and 1.57, respectively, for the outline-based measure of change in aspect ratio. These are consistent with the data of Blakenship et al [9].

Having analyzed the T1 transition and rosette transformation separately, we turn now to tissue-wide features of the simulation with both processes occurring. Again, we have carried out 10 realizations from identical initial conditions. Their outcome are qualitatively similar but quantitatively different thanks to the temporal and spatial stochasticity built into the model. We will average among the 10 runs to measure the various aspects of tissue elongation.

Movie 3 in the SI depicts a representative run of the full simulation, which lasts 3000 s with 41 T1 transitions and 4 rosettes during the simulation. On average, the AP contraction phase of a T1 transition takes 438 s, while the DV elongation phase 270 s. These are close to the times in the simulations of figure 4 with only T1 transitions. Rosette formation and resolution take considerably longer times (figure 8). The formation and resolution times of 643 s and 523 s compare with experimental values of 648 s and 600 s [9]. As noted earlier, the broad distributions of the contraction and resolution times reflect the stochastic disturbance from the neighboring cells that may delay or precipitate cell intercalation.

Zoom In Zoom Out Reset image size

Depending on how nearby AP edges are chosen randomly to form rosettes, the total numbers of rosettes and T1 transitions vary among our 10 simulations, the number of rosettes ranging from 4 to 6 while that of T1 from 37 to 41. These are much lower than those recorded in vivo [13, 31]; based on the more recent experimental data of [31] one would expect some 35 rosettes and 63 T1s for 52 cells. Furthermore, the final configuration is always the same in the simulation, with hexagonal cells in neat rows. Both features are symptomatic of the simple and regular configuration of our cells in the model tissue. The total numbers of T1 and rosettes are limited by the number of AP edges available. In vivo, cells exhibit considerable variations in shape, size and connectivity, and the germband flows and deforms in time to create new AP edges. Thus, the tissue undergoes multiple rounds of intercalation during germband extension [9]. This should be kept in mind when comparing our model predictions with experiments in the following.

Anisotropy in the distribution of Rok, Baz and Myo-II and the edge tension are tabulated in table 2, in comparison with experimental data in the literature. In the full simulation, the AP:DV ratios in T1 transitions are close to those in the T1-only simulations of table 1. What is of more interest here is the ratio between linked AP edges and DV edges in rosettes, computed by averaging the three linked edges. The Rok and Baz ratios are more strongly polarized than those in T1 transitions. They are also comparable with experimental data for AP:DV ratios, although the experimental sources [14, 18, 31] did not distinguish between T1 and rosettes, and probably contained both. In our model, the linked AP edges exhibit a Myo-II level more than twice that of isolated AP edges. Partly this is due to the elevated Shroom level on linked AP edges. Moreoever, Myo-II enjoys a positive feedback with tension [10]. The strong tension on the linked AP edges has given Myo-II an additional boost. Again, most experimental reports did not distinguish the linked and isolated AP edges. For example, Farrell et al [31] reported an AP:DV myosin ratio of 2.1 for both types of AP edges. Besides, this ratio includes edges within 15° of the vertical and horizontal axes. Considering these two factors, the real linked AP:DV ratio for myo-II in vivo is likely much higher than 2.1. The only experimental study that distinguishes myosin levels on linked and isolated AP edges gives the following ratios [10]: linked AP:AP:DV  =  1.14:1.05:1. Though the ratios are lower than in the other reports [14, 18, 27, 31], the linked AP edge still exhibits a clear advantage over the isolated AP edge in terms of myosin enrichment relative to the DV edge. Note that the strain-activated myosin recruitment mechanism may have contributed to the myosin enrichment on linked AP edges in the experiments [33]. Finally, the model successfully captures the strong tension force on linked AP edges, relative to isolated AP edges, and the tension ratio is not far above data gathered from laser ablation experiments [10, 15].

Table 2. Anisotropy of various protein species and the edge force: comparison between model predictions and experimental measurements. The AP:DV ratio corresponds to edges in T1 transitions, while Linked AP:DV refers to edges in rosettes.

	Experiments AP:DV	Experiments Linked AP:DV	This model AP:DV	This model Linked AP:DV
Rok	1.5–1.9 [14]	(No data)	$1.35\pm 0.09$	$1.99\pm 0.15$
Baz	0.33–0.50 [14, 18, 31]	(No data)	$0.83\pm 0.08$	$0.54\pm 0.05$
Myo-II	1.05–2.1 [10, 14, 18, 27, 31]	1.14 [10]	$1.60\pm 0.40$	$3.70\pm 1.0$
Edge force	1.7–2.55 [10, 30]	2.9–3.1 [10, 15]	$1.40\pm 0.80$	$3.60\pm 1.30$

Figure 9 compares the model predictions of the tissue stretching for the T1 only simulations of section 3.1 and the full simulations of this section. The in vivo data of Collinet et al [6] are also shown for comparison. First, the inclusion of rosettes in our model leads to a slight increase in the germband extension, consistent with the observations in figure 7. The effect is weak because the simulation sees much fewer rosettes than T1's. Second, our model predictions of relative length fall some 15% below the in vivo data. In the end, our relative length L_r reaches about 1.7, whereas in vivo L_r rises above 2. The reason for this discrepancy, as alluded to already, lies in the fact that the simple geometry of our model tissue allows just one round of intercalation, after which no vertical AP edge remains and no further intercalation can occur. In vivo, on the other hand, the germband undergoes complex overall deformation that allows movement and rotation of cell edges and multiple rounds of T1 and rosette processes. This is clear from videos of germband extension in Blankenship et al [9]. Therefore, the simplicity of our model tissue has limited its ability to capture the true magnitude of germband extension.

Zoom In Zoom Out Reset image size

To test the spatial homogeneity of stretching, we again compare the strain (or relative length) of the 3 rectangular clusters depicted in figure 5. Results show that such homogeneity remains in the full simulation; the relative length agrees to about $4\%$ among the 3 boxes. This is perhaps no surprise given that the T1 transitions alone incur uniform deformation throughout the tissue (figure 5), and that rosettes produce only slightly larger tissue strain than T1, and are relatively few in our model predictions.

Finally, we discuss the effect of the 'extrinsic force' due to posterior midgut invagination (PMI) on germband extension [5, 6]. As explained in the SI (section S4), the PMI effect is imposed as a strain rate along the AP axis. If we use too low a strain rate in our model, the interior cells tend to develop concave shapes, possibly another symptom of the overly simple topological setup of the model tissue. This limits us to a relatively narrow range of variations for the PMI strain rate. Figure 10 shows that PMI contributes modestly to the total germband extension. With the PMI strain rate at 70% or 150% of the baseline value, the relative length at the end goes down by 0.8% or up by 1.5%, respectively. Thus, the model prediction is consistent with the consensus in the literature that germband extension is mostly due to cell intercalation, but tissue-level external forces also contribute to it [5, 6, 10].

Zoom In Zoom Out Reset image size