Global morphogenetic flow is accurately predicted by the spatial distribution of myosin motors

1Animal development is characterized by highly dynamic rearrangements of mechanically coupled cells. Such rearrangements must be tightly coordinated across the embryo to achieve normal morphology and organogenesis. During gastrulation of Drosophila, for example, the embryonic blastoderm – an epithelial monolayer of about 6000 cells on the surface of the embryo – undergoes a dramatic deformation that changes tissue topology and gives rise to the three germ layers. These processes involve a coherent flow of cells along the surface of the epithelial monolayer, which in turn drives folding and defines future shape of the embryo. The most prominent aspects of gastrulation are the formation of the ventral furrow which initiates the invagination and internalization of the mesoderm, and germ-band extension which involves convergent extension of the lateral ectoderm and the flow of the ventral germ-band onto the dorsal side of the embryo. Both of these processes have been extensively studied, leading to the identification of developmental patterning genes specifically required for each process. Live imaging has also uncovered process-specific cell behaviors such as apical constriction of presumptive mesoderm cells during ventral furrow formation and intercalation of neighboring cells in the lateral ectoderm during convergent extension. These behaviors are associated with localized activity of the force generating non-muscle myosin II. However, despite considerable understanding of the local processes involved in such cellular rearrangements, a coherent picture of global morphogenetic flows has remained elusive.

2Understanding how cell flows are coordinated across different cell populations requires distinguishing the roles of local cell behavior and long-range intercellular interactions. To what extent is the transformation of tissue driven locally by the processes associated with cells at that position? How important is the long-range interaction between different regions? In the context of the fly embryo, VF formation seems well explained locally by the apical area contraction of ventral mesoderm cells. On the other hand, the non-local interactions between the VF constriction (or posterior midgut invagination) and the convergent extension of lateral ectoderm remain a subject of active study which requires quantitative multi-scale analysis.

3There are two complementary approaches towards quantitative analysis of tissue flow. One approach focuses on cell-scale behavior aiming to decompose tissue flow into specific cellular processes such as cell-shape change and intercalation Alternatively one can ‘zoom out’, taking a continuum mechanics approach that aims to describe tissue flow on the whole organ scale. This coarse-grained mesoscopic perspective captures correlations in cell behavior which stem from intercellular interactions and the supra-cellular organization of the cytoskeleton in epithelial tissues. In biophysics, the continuum mechanics approach has been developed to understand the behavior of active gels during myosin driven viscoelastic flow and has been successfully used to model cortical flows in C. elegans zygotes at the first-cleavage state. Recent theoretical work provides a bridge between cell-based and meso-scale continuum descriptions, focusing on the non-trivial consequences of stress-dependent active cytoskeletal processes. Here, we shall use continuum mechanics approach to set up a framework for predicting global tissue flow at the whole organ level.

4The main advantage of the continuum mechanics approach is its ability to capture key aspects of force balance associated with local deformation and flow. It allows to describe quantitatively, with only a few parameters, how the effect of local forcing spreads across a tissue. The tendency of cells to stick together and resist deformation results in a non-local relation between the myosin activity that drives the flow and actual flow velocities. The continuum mechanics approach therefore enables one to test different hypotheses helping to identify key contributing processes. For example, the question of local versus non-local response in the continuum mechanics approach translates into specific hypotheses regarded force balance: are the myosin-generated forces balanced locally by traction relative to a substrate or do they propagate within the epithelium layer through cell deformation and viscous coupling? Quantitative analysis can then be used to build, starting from the simplest model, a sequence of approximations that capture biological reality in increasing detail.

5We shall describe below a novel image analysis-based approach that will use continuum mechanics to quantitatively relate different observables and will show that myosin distribution and anisotropy on mesoscopic scale is a fully adequate proxy of physical stress, thereby enabling a surprisingly predictive description of global flow. Specifically, we shall show that embryo-scale tissue transformations during Drosophila gastrulation are represented by a temporal sequence of three topologically distinct flow field configurations. Each phase is accompanied by a characteristic spatial distribution of myosin molecular motors both on the basal as well as apical cell surface, which we quantify in terms of a coarse-grained ‘myosin tensor’ that captures both myosin concentration and anisotropy. To relate the observed global flow fields to myosin apical and basal distributions we assume that tissue flow is driven by stress proportional to the myosin tensor, and is effectively viscous with two parameters: effective shear and bulk viscosities, the latter controlling the compressible component of the flow. With a total of just three global parameters (only one of them time dependent), this simple model achieves remarkable agreement between predicted and measured spatio-temporal pattern of the flow. The analysis uncovers the importance of a) spatial modulation of myosin distribution and b) the long-range spreading of its effect due to mechanical interaction of cells. In particular, we find that transition to the germband extension phase of the flow is associated with the onset of effective areal incompressibility of the epithelium, which makes the relation of the flow and myosin forcing strongly non-local. Our quantitative analysis also reveals a new function for basal myosin in generating a dorsally directed flow and, combined with mutant analysis, points to an unconventional control mechanism of this function through twist -dependent reduction of basal myosin levels on the ventral side. Finally, we shall argue that the ability to quantitatively describe the relation between the flow and its myosin-generated forcing provides a new approach to the study of the processes that control morphogenetic transformations, which can be used to disentangle novel control mechanisms such as mechanical feedback from the effects of gene expression patterning.

6To enable our study, we generated a pipeline that combines in toto light sheet microscopy (Figure 1a), tissue cartography (Figure 1b), and segmentation-free anisotropy detection to quantify global tissue flows, and myosin activation patterns (Figure 1—figure supplements 1–4). Using optical flow velocimetry applied to cylinder projections of the Surface of Interest (SOI) passing through cells below the apical cell surface (see SI for details Figure 1—figure supplement 5c), we find that tissue remodeling during Drosophila gastrulation is characterized by three simple flow field configurations (Figure 1c–e). The earliest flows start well before the ventral furrow (VF) forms, and are characterized by a dorsal sink and ventral source (Figure 1c). In contrast to the VF, no cells are internalized during this flow, but rather cells reduce cross section on the dorsal side (Figure 1—figure supplement 1c). As the VF forms, source and sink swap sides and a large group of cells internalize on the ventral side, as mesoderm precursors leave the surface of the blastoderm (Figure 1d). During germband extension (GBE), the flow pattern exhibits two saddles arranged on the dorsal and ventral sides as well as four vortices, two in the posterior and two in the anterior end (Figure 1e). Each of the three flow fields is accompanied by a typical spatial myosin configuration. The pre-VF flow associates with basal myosin that exhibits a pronounced Dorso Ventral (DV) asymmetry, with high levels of myosin on the dorsal and low levels on the ventral side (Figure 1f), while the apical pool appears uniform across the surface (Figure 1i, Figure 1—figure supplement 5d–g). The basal pool remains asymmetric during VF flow (Figure 1g), but the apical pool now also develops DV asymmetry in reversed orientation (Figure 1j). The asymmetry on the apical surface becomes further pronounced in the GBE-phase (Figure 1h,k, Figure 1—figure supplements 6 and 7).

7Global changes in myosin pools are a hallmark for transitions in flow field configuration (Figure 2a, Figure 1—figure supplements 1, 2, 3 and 4). Myosin is initially enriched in the basal pool, and as sink and source swap position, it begins to accumulate on the apical side. While the basal pool is isotropic (Figure 2—figure supplement 4a), cortical myosin on the apical cell surface is known to polarize during convergent extension. To quantify this effect at the tissue level, we developed an automated segmentation-free anisotropy detection algorithm (Figure 2b, Figure 2—figure supplement 3a,b). Available methods for anisotropy detection mostly operate at the single cell level and construct a nematic tensor by integrating signal intensities along cell outlines. At the organismal scale membrane segmentation is costly, and often fails to define closed outlines of cells using only a polarized membrane marker. We overcome the need for fiduciary markers that increase experimental complexity by shifting the perspective to cell edges and using the Radon transform to implement a robust and rapid segmentation-free algorithm for detecting course-grained anisotropy (Figure 2b). Radon transforms integrate signal along lines of given orientation and normal distance from the origin. In this way, edges are mapped to peaks that reflect the total intensity along the length of an edge (Figure 2b) (see Appendix 1 for detail). Edge orientation and average myosin intensity are described by a 2 × 2 symmetric matrix (of rank one) defining the local ‘myosin tensor’ (Figure 2b). By averaging the resulting tensors in a given region, we obtain a quantitative description of local tissue anisotropy and overall levels that reflects the intensity-weighted average of cell edges. The resulting course-grained tensor has a non-zero trace, and thus can be separated into an isotropic and a traceless anisotropic part (Figure 3a, Figure 2—figure supplement 3b), from which we construct a measure for anisotropy (which is low in the basal pool Figure 2—figure supplement 4a). The anisotropic signal in the apical pool starts out low, but increases from about 8 min corresponding to late stage 6 (Figure 2a, Figure 2—figure supplement 4b,c). The anisotropy axis, readily computed by the eigenvectors of the myosin tensor, aligns well with local tangent to pair rule gene expression boundaries (Figure 2a,d, Figure 2—figure supplement 4d). This is the expected result given that anisotropies are thought to be driven by the patterned juxtaposition of pair-rule gene expression.

8We have examined and quantified tissue flow and myosin distribution in multiple (N = 22) wild-type embryos and found it highly reproducible (Figure 1—figure supplement 5). For the purpose of quantitative analysis presented below, we shall use the (suitably aligned) ‘ensemble’-averaged flow and myosin distribution (see SI).

9To relate myosin to stress, we assume signal intensity is proportional to myosin motor concentration and its local activity. The latter – pulling on cytoskeletal actin filaments – generates local force dipoles, which can be explicitly described in terms of local stress tensor (see Appendix for details). On the coarse grained level, resulting stress would be defined by the activity weighted average over filament orientations and hence proportional to the myosin tensor as we define it. The resulting force per unit area of the epithelial layer is then proportional to the divergence of the myosin tensor. Note that the isotropic component of the myosin distribution (observed both in the apical and the basal pool) also generates a force that is proportional to the gradient of the measured concentration intensity profile (Figure 1f–k).

10To relate myosin generated stress to morphogenetic flow, we assume that on the mesoscopic scale tissue flow is governed by effective viscoelasticity which arises from the mechanical properties of the underlying cytoskeletal network within the two dimensional epithelial layer of cells. This model assumes that on short time scales tissues respond elastically to mechanical perturbations, yet on longer time scales elastic stress is relaxed through active rearrangement of the cytoskeleton as cells adapt to the imposed deformation. On the longer time scale tissue dynamics can be described by a two-dimensional effective viscous flow equation with two effective viscosity parameters that (see Appendix 2) are directly related to the two elastic constants: shear modulus (controlling ‘sliding’ of cells relative to each other) and the planar bulk modulus (controlling areal compression or dilation) (Figure 3a, Figure 3—figure supplement 1b). We note that effective viscosity spreads the impact of local forcing, generating a non-local response so the flow at any given point integrates the influence of forces acting all over the embryo. Inverting the equation using the finite element method, we obtain a quantitative prediction for the flow field generated by measured myosin localization patterns (see SI for details). Our model has only three global parameters: the ratio of effective viscosities, and the conversion factors relating normalized apical and basal myosin intensity to stress (Figure 3—figure supplement 1c). To keep the model as simple as possible, we do not allow spatial dependence of these parameters and keep conversion factors constant, leaving the ratio of viscosities as the only time-dependent global fitting parameter (Figure 3—figure supplement 1b).

11Even without spatial modulation of the parameters, the model achieves about 90% accurate description of the flow pattern before and after VF invagination (see Figure 3b). The main discrepancy of model predictions for pre-VF flow (see Figure 3c) is a displacement of sink and source positions along the AP axis by less than 30 µm. Prediction of GBE flow essentially agrees with measurements across the entire embryo, with the exception of a domain close to the vortices on the posterior end, due to a mismatch of fixed-point location (Figure 3d). Remarkably, our model is even able to correctly predict subtle differences between anterior and posterior fixed points along the DV axis (Figure 3d). Measured flow is first dominated by sources and sinks that disappear later during GBE, suggesting that before and during VF invagination cells are less resistant to surface area compression than during GBE. Indeed, quantitative comparison with an independently measured flow field (Figure 1c–e) shows that the ν20 /ν27 ratio increases dramatically at the start of GBE phase (corresponding to the relative increase of the underlying 2D bulk modulus, see Appendix 2, Figure 3—figure supplement 1b) resulting in effective incompressibility of apical surface of cells. The temporal coincidence between completion of cellularization and increase of the bulk modulus provides an intriguing possible explanation of how the continuous transition in our time dependent variable might be realized. Poor agreement during VF invagination is due to a significant fraction of cells internalizing and thus leaving the surface. To account for this effect, we extend the model to allow a ‘cut’ in the lattice along ventral midline with an imposed in-plane boundary force (perpendicular to the cut) representing the pulling effect of the VF (see SI for detail, Figure 3—figure supplement 1a). This relatively simple extension allows to recover ~90% accuracy (Figure 3be), illustrating how regional inhomogeneity associated with particular morphogenetic events could be quantitatively captured by suitable generalizations.

12To evaluate the fit obtained in wild type embryos, we examined flows in mutant embryos in which the distribution of myosin is altered. Analysis based on tissue tectonics has shown that strain rates in twist (twi) embryos, which lack the VF, exhibit slower kinetics compared to WT, however, the cause of this remains a subject of debate. We have quantified the flow field and myosin activity patterns in twi mutants (Figure 4—figure supplement 1), and find that our model is able to accurately predict the flow profiles (Figure 4a). During early flow phases – corresponding to times of pre-VF flow in WT – DV asymmetry of the basal myosin pool is strongly reduced in comparison to WT, as is tissue movement towards the dorsal pole (Figure 4b, Figure 1—figure supplement 1, Figure 4—figure supplement 1a,d). Moreover, anisotropy of the apical myosin pool increases at a slower rate as compared to WT. As previously reported for strain rates, this is most pronounced for the first 20 min (Figure 4c, compare Figure 4 – fig. supplement = 0 with Figure 2—figure supplement 3e,f). In bcd nos tsl (bnt) embryos lacking all AP patterning, the early basal DV asymmetry is similar to WT, with only slightly reduced myosin asymmetries and dorsal movement (Figure 4b, Figure 2—figure supplement 1). At later stages, however, anisotropy of the apical myosin pool remains low and comparable to pre-VF WT levels. This result is expected given the uniform expression of pair-rule genes in a bnt genetic background (Figure 4c). Consistent with these myosin distributions, we see the early dorsal flow associated with basal myosin asymmetry but a failure to produce the more complex later flow patterns with their characteristic saddles and vortices. On a quantitative level our model’s predictive power for AP patterning deficient bnt mutant embryos is comparable to WT and twi mutants (Figure 4a).

13In summary, we have presented a simple biophysical model of morphogenetic flow that quantitatively describes complex tissue motion in terms of a hydrodynamic equation parameterized by two effective viscosities. The flow is driven by the stress defined by a linear superposition of two myosin tensors describing the apical and basal myosin pools. We propose that the basal myosin pool forms an isotropic and contiguous network, contracting in a similar fashion as purified actomyosin gels in vitro. Imbalance within this network, caused by the twi dependent depletion of myosin on the ventral side, drives global dorsal-ward flow in the pre-VF phase, which continues to contribute until early GBE (Figure 4d). Interestingly, in silico perturbations indicate the local depletion (on the ventral side) has a global effect, most evidently manifested by a ‘sink’ on the dorsal side, which is lost in a simulation using the same parameters as WT, but no DV modulation of the basal myosin pool (Figure 4d, left panels). The apical pool decomposes into isotropic and anisotropic components. In addition to previously described accumulation of isotropic myosin in ventral regions, we observe a striking gradient of anisotropic apical myosin along the DV axis, reaching highest levels in lateral ectoderm and lowest levels in amnioserosa tissue at the dorsal pole. Because the force driving the flow arises from the non-uniformity of the stress, this modulation of myosin distribution is critical for the dynamics. While local rate of cell intercalation is often interpreted in terms of local myosin distribution on cellular and sub-cellular scales, our model shows that the local rate of strain is a result of the tissue-wide distribution of forces generated by the spatial non-uniformity of myosin (mathematically described by the divergence of the myosin tensor) (Figure 4d, right panels). The importance of spatial modulation (Figure 4d) suggests a novel role of the dorsal signaling pathway in generation of GBE flow. Surprisingly, in twi mutants both the rate of increase as well as the peak myosin anisotropy are significantly reduced in the first 20 min of GBE flow (Figure 4c, Figure 4—figure supplement 1). The reduced intercalation and strain rates observed in these mutants has been previously reported, and interpreted in terms of possible generation of AP forces by the internalized VF (absent in twi mutants). Our model accounts for the reduced rate of strain in terms of the changes in spatial distribution and the reduced level of myosin anisotropy. This however brings up the question of how elimination of twi expression in the ventral mesoderm affects myosin anisotropy in the lateral ectoderm. We suggest that this effect may be due to mechanical feedback on myosin recruitment, which relates the later to local rate of strain. Through this ‘dynamic recruitment’ effect, changes in the ventral region that modify global flow patterns can affect myosin distribution and anisotropy in the lateral region. In this way, local modification of the myosin pattern can produce not only a non-local perturbation of the flow, but also a non-local perturbation of myosin distribution. The global nature of the flow is reinforced by the observed transition towards areal incompressibility at the onset of GBE-flow, which together with reduced polarization kinetics and reduced strain rates observed in twi mutants, indicates that non-local consequences of stress generated largely in lateral ectoderm can account for the dorsal movement of the posterior midgut.

14Taken together, our observations show that morphogenetic flow is a global response to local forcing which arises from the spatial modulation of myosin density and anisotropy. The latter is derived from the spatial patterns of developmental transcription factors, but we suggest may also involve mechanical feedback affecting recruitment of myosin. Our quantitative approach provides a framework for integrating the effect of local factors in the description of the global flow.

15His2Av-mCherry, bcde1nosbntsl4/TM3, halo twiID96/Cyo (twiID96 is also known as twi), sqh-GFP klar, OregonR. Embryos where dechorionated following standard procedures, and mounted in agarose gels as previously described.

16Fluorescence-based live imaging was carried out on a MuVI SPIM. Briefly, the optics consisted of two detection and illumination arms. Each detection arm forms a water-dipping epifluorescence microscope, consisting of an objective (Apo LWD 25x, NA 1.1, Nikon Instruments Inc.), a filter wheel (HS-1032, Finger Lakes Instrumentation LLC), with emission filters (BLP01-488R-25, BLP02-561R-25, Semrock Inc.), tube lens (200 mm, Nikon Instruments Inc.), and an sCMOS camera (Zyla 4.2, Andor Technology plc.), with an effective pixel size of 0.26 µm. Each illumination arm consisted of a water-dipping objective (CFI Plan Fluor 10x, NA 0.3), a tube lens (200 mm, both Nikon Instruments Inc.), a scan lens (S4LFT0061/065, Sill optics GmbH and Co. KG), and a galvanometric scanner (6215 hr, Cambridge Technology Inc.), fed by lasers (06-MLD 488 nm, Cobolt AB, and 561LS OBIS 561 nm, Coherent Inc.). Optical sectioning is achieved by translating the sample using a linear piezo stage (P-629.1cd with E-753 controller) sample rotation is performed with a rotational piezo stage (U-628.03 with C-867 controller) and a linear actuator (M-231.17 with C-863 controller, all Physik Instrumente GmbH and Co. KG).

17Stages and cameras are controlled using Micro Manager, to coordinate time-lapse experiments, running on a Super Micro 7047GR-TF Server, with 12 Core Intel Xeon 2.5 GHz, 64 GB PC3 RAM, and hardware Raid 0 with 7 2.0 TB SATA hard drives. Samples were recorded from two, by 900 rotated views, at a typical optical sectioning of 1 µm, and temporal resolution of 75 s. As previously described, MuVI SPIM optical stability allows a fusion strategy based on a diagnostic specimen. Recorded once per experiment, the diagnostic specimen is used to determine an initial guess for an affine transformation, which we feed into a rigid image registration algorithm, to fuse individual views, resulting in an isotropic resolution of. 26 µm in the registered image.

18We used tissue cartography to extract surfaces of interest (SOI) from embryos. Briefly, we identify the outline of the sample using the Ilastik detector to determine a point cloud for SOI construction. In a fitting step (implemented in the sphere-like fitter), we create a smooth description of the SOI in terms of cylinder coordinates defined by AP axis and azimuth. Image intensity data are then projected onto a nested group of 5 layers (two normally evolved layers above and below the SOI), each three pixels apart, defining a 3–4 µm thick 'curved image stack'. For analysis the maximum intensity projection of nested layers was used. Although the shape of embryos of the same genotype is highly reproducible, small differences in the underlying point cloud can result in small differences of the SOI passing through the apical cell surface. To simplify comparison between embryos, we create a standard projection on a cylinder grid of fixed size, with the embryo surface oriented such that apical is left, posterior right, dorsal in the center, and ventral on top and bottom (Figure 1). Systematic distortions of measurements due to projecting the curved embryo surface to the plane are corrected using the metric tensor.

19The apical surface is static, while the dynamic basal cell surface moves with the cellularization front. Projections of the latter could be created by reading signal on a surface obtained by evolving the apical SOI along its normal basal wards. However, small differences in cell height (<10% of a typical cell height, and ~1% of the embryo diameter), could result in small but systematic bias of the SOI around the cellularization front and impair projection quality. We avoid this problem by determining a new point cloud for each time point, for which we focus the ilastik detector on the interface between basal myosin and yolk. Our model approximates the embryo as a thin shell (see below), and hence as 2D surface. Therefore, we map the dynamic basal cell surface onto the cylinder grid of the static apical SOI.

20We measure the flow field using the particle image velocimetry (PIV) method, that identifies local displacements between two time points. Briefly, we implemented the phase correlation method that leverages favorable execution times of fast fourier transforms, to estimate local flow in a region of interest on the projections. To minimize effect from systematic distortions towards polar regions on cylinder projections, we adjust the size of the ROI according to the local metric strain, which we define as the deflation of the metric from flat space.

21Although gastrulation in Drosophila is highly reproducible from embryo to embryo, in practical terms experiments are subject to a constant time shift, depending on the developmental stage of the sample at the start of imaging. Thus, we developed an automated routine that allowed us to identify a common time frame that the 36 (WT: N = 22, twist:7, bcd nos tsl:7) live-imaged embryos are registered to. Specifically, we introduced a constant time shift for a given flow field that minimizes the squared difference with respect to a reference flow field averaged over the embryo surface:where denotes averaging across the embryo, is an arbitrary chosen reference from the ensemble, and, denote the i-th flow field and offset time respectively. In this way, we align samples to a chosen reference, in which we use the first occurrence of the cephalic furrow (CF) as a landmark indicating our choice for min. Within a given genotype, we could automatically determine the offset time. However, to align mutants to WT, we first aligned all mutant datasets, and then used landmarks such as the CF (twist), or the VF (bcd nos tsl), to define a common time frame as best as possible.

22Time shifted accordingly, we created an ensemble average flow field for each genotype:

23The magnitude of ensemble average is highly reproducible from embryo to embryo (note the small standard deviation Figure 1—figure supplement 5a). Flow trajectories during cellularization point towards the dorsal side (Figure 1—figure supplement 5b), showing persistent movement towards dorsal regions during pre-CF flow. This is accompanied by reduction of apical cell area in these regions (Figure 1—figure supplement 5c), as measured using confocal microscopy. While the length of pre-VF flow lines peaks on the anterior and posterior poles in WT, it is substantially reduced near poles in twist, and only mildly reduced in bcd nos tsl (Figure 1—figure supplement 5d–f). Together with the loss of basal DV asymmetry, this suggests that twist mediated reduction of basal myosin levels on the ventral side is responsible for dorsal-ward flow.

24Using the imaging and pre-processing procedure as outlined above with samples expressing sqhGFP, we created projections of the apical and basal cell surfaces, with the goal of establishing a quantitative measurement of global myosin patterns (Figure 1—figure supplements 6 and 7). Ideally, quantification of signal intensities is carried out using identical conditions for each sample in the pool used for statistics, to minimize variability across samples. However, when performing in toto live-imaging, it is difficult to image more than one sample at a time and keep a high recording frequency. To minimize variability in a sequential recording scheme, we keep imaging conditions constant, but there are still possible variabilities in recorded signal intensity for biological but also technical reasons.

25To account for such variability between experiments, we normalize recorded data (Figure 1—figure supplement 6b,b’,c). Signal intensity of all time points in a given experiment are summarized by normalizing the intensity distribution: upper and lower range are determined according to the and -percentile; normalization is done by subtracting the and dividing by, yielding a dimensionless normalized signal distribution (compare Figure 1—figure supplement 6b and b’). This strategy should not only allow for comparison on the same microscope, but also across microscopes, allowing for validation of in toto live-imaging from sequential experiments against fixed batches imaged e.g. on a confocal.

26Here, we briefly outline the results for the DV asymmetry in the basal myosin pool that we reported in the main text. First, we time align intensity normalized basal projections as described above. Next we convolve each pullback with a Gaussian of width cell diameters to obtain basal myosin at the mesoscale (see discussion in model section below for definition). The results are then ensemble averaged to obtain ensemble myosin distribution as shown in Figure 1—figure supplement 6. To assess DV asymmetry, we focus on the region outlined by the black/white dashed line, where we first take an average along the AP axis and then compute average signal on dorsal side, and subtract from it the average signal on the ventral side. Repeating the outlined routine for all time points, we obtain the plot show in main text Figure 4b.

27Figure 1—figure supplement 7c shows DV cross sections of fixed embryos stained for rb anti zipper and mouse anti dorsal, cut along the AP axis, and imaged on a confocal microscope. DV orientation of the samples is automatically determined based on the dorsal signal. To estimate the age of fixed embryos in relation to live-imaging data, we constructed a calibration curve for cell apico-basal height shown in Figure 1—figure supplement 7a. Using the known monotonic relation between cell height and age, which we find lasts until about after CF formation, we obtain estimate for the age of a fixed embryo based on measuring cell height. By segmenting the outline of basal myosin, we can then measure DV asymmetry in the same way as described for live-imaging data above.

28A direct comparison between live imaging-based DV asymmetry measurement, and based on fixed DV cross sections from confocal shows that after applying normalization routine as described above, we find similar estimates for the DV asymmetry using both light sheet and confocal imaging (see Figure 1—figure supplement 7b). Note that uncertainties in the calibration curve propagate to exact age determination in the embryo, and thus increased fluctuations in DV asymmetry determined using confocal imaging.

29Inversion of the continuum equation of state relating the coarse-grained myosin tensor and cellular flow-field was achieved using Finite Element Methods (FEM) in the weak formulation implemented within the FELICITY toolbox for MATLAB. Equations were inverted on a static triangular mesh representing the 'canonical' embryo surface produced via a point cloud (described above) subsequently turned into a smooth triangulation using MeshLab. As such, all objects within the equation of state had to be parameterized within the 3D embedding space of the mesh, which can be done by using the direction relation between projections and SOI. The only dynamic input to this inversion algorithm is the divergence of the myosin tensor. This was computed by interpolating the gradient of each Cartesian component of the tensor onto triangular faces of our mesh, producing a 3 × 3 × 3 object on each face. The partial trace of this object over directions within the tangent plane of the face result in the estimated divergence. This operation was repeated for all myosin pools used. All equations within the FEM software are projected onto the surface of our 3D mesh to manually ensure solutions only exist within the tangent plane.

30To benchmark the quality of our FEM solver, we implemented the equation of motion on a sphere and tested results for known solutions, which confirmed our solver works within the expected numerical accuracy. To test dependence on discretization used, we compared our results using different meshes at varying mesh sizes, and found good agreement with all test cases.

31In order to model internalization of the VF, we introduced the ability to add a 'cut' within the triangular mesh. Contraction of the tissue was modeled by manually introducing local force dipoles on edges within the mesh pointing along the bond. All vertices along the cut were given zero bulk modulus to allow for local tissue compression needed to simulate invagination. The location of the cut was estimated from the PIV flow fields. Specifically the ratio of the divergence of the velocity field to the velocity field's magnitude was used to estimate the spatial extent of the cut over times during ventral furrow formation.

32Predictions for the flow field obtained via inversion are subject to an overall scale factor, that can't be determined by the model. To compare ensemble averaged flow field measurement to model predictions in a quantitative fashion, we define a global measure for the spatial residual that is insensitive to such a scale factor. With the short hand notation to define overall magnitude of the field u across the surface of the embryo ( denotes averaging the space dependent field across the embryo surface, so is not space dependent.), the residual is defined asprovides a spatial discrepancy map, indicating the prediction quality as a function of location on the embryo, that is in-sensitive to noise dominated fluctuations in domains of no flow (i.e. fixed points), as opposed to e.g. inner product.