1. Developmental Biology
  2. Physics of Living Systems
Download icon

Theoretical tool bridging cell polarities with development of robust morphologies

  1. Silas Boye Nissen
  2. Steven Rønhild
  3. Ala Trusina  Is a corresponding author
  4. Kim Sneppen  Is a corresponding author
  1. University of Copenhagen, Denmark
Research Article
  • Cited 0
  • Views 1,080
  • Annotations
Cite this article as: eLife 2018;7:e38407 doi: 10.7554/eLife.38407

Abstract

Despite continual renewal and damages, a multicellular organism is able to maintain its complex morphology. How is this stability compatible with the complexity and diversity of living forms? Looking for answers at protein level may be limiting as diverging protein sequences can result in similar morphologies. Inspired by the progressive role of apical-basal and planar cell polarity in development, we propose that stability, complexity, and diversity are emergent properties in populations of proliferating polarized cells. We support our hypothesis by a theoretical approach, developed to effectively capture both types of polar cell adhesions. When applied to specific cases of development – gastrulation and the origins of folds and tubes – our theoretical tool suggests experimentally testable predictions pointing to the strength of polar adhesion, restricted directions of cell polarities, and the rate of cell proliferation to be major determinants of morphological diversity and stability.

https://doi.org/10.7554/eLife.38407.001

eLife digest

Cells have the power to organise themselves to form complex and stable structures, whether it is to create a fully shaped baby from a single egg, or to allow adult salamanders to grow a new limb after losing a leg. This ability has been scrutinised at many different levels. For example, researchers have looked at the chemical messages exchanged by cells, or they have recorded the different shapes an embryo goes through during development. However, it is still difficult to reconcile the information from these approaches into a description that makes sense at multiple scales.

When an embryo develops, sheets of cells fold and unfold to create complex 3D shapes, like the tubes that make our lungs. Moulding sheets into tubes relies on interactions between cells that are not the same in all directions. In fact, two types of asymmetry (or polarity) guide these interactions. Apical-basal polarity runs across a sheet of cells, which means that the top surface of the sheet differs from the bottom. Planar cell polarity runs along the sheet and distinguishes one end from the other. For instance, apical-basal polarity marks the inner and outer surfaces of our skin, while planar cell polarity controls the direction in which our hair grows.

Nissen et al. set out to investigate how these polarities help cells in an embryo organise themselves to form complicated folds and tubes. To do this, simple mathematical representations of both apical-basal and planar cell polarities were designed. The representations were then combined to create computer simulations of groups of cells as these divide and interact with each other.

Simulations of ‘cells’ with only apical-basal polarity were able to generate different shapes in the ‘tissues’ produced, including many found in living organisms. External conditions, such as how cells were arranged to start with, determined the resulting shape. With both apical-basal and planar cell polarities, the simulations reproduced an important change that occurs during early development. They also replicated how the tubes that transport nutrients and oxygen form.

These results show that simple properties of individual cells, such as polarities, can produce different shapes in developing tissues and organs, without the need for a complicated overarching program. Abnormal changes in cell polarity are also associated with diseases such as cancer. The mathematical model developed by Nissen et al. could therefore be a useful tool to study these events.

https://doi.org/10.7554/eLife.38407.002

Introduction

Multicellular organisms are amazing in their ability to maintain complex morphology in face of continuous cell renewal and damages. Adult salamander can regenerate entire limbs (Eguchi et al., 2011), and, during development, some regions can maintain patterning when moved to different parts of an embryo or if the size is varied (Lyons et al., 2012). Given the vast complexity and diversity of living shapes, how can we reconcile the robustness to perturbations with flexibility to diversify? While undoubtedly the end result is encoded in the DNA and protein networks, looking for an answer at this level is challenging. Examples of phenotypic plasticity (Libby and Rainey, 2011), convergent evolution, and contrasting rates of morphological and protein evolution (Cherry et al., 1979) show that morphological similarity may not couple to the protein sequence similarity (Gould, 1970). Inspired by the unfolding of morphological complexity in development, we propose that cellular polarity may be the key for reconciling complexity, robustness, and diversity of organismal morphologies.

During early development, the increase in morphological complexity coincides with the progressive polarization of cells – first apical-basal (AB) polarity and then planar cell polarity (PCP) (Müller and Bossinger, 2003; Roignot et al., 2013; Andrew and Ewald, 2010; Li and Bowerman, 2010). This theme is ubiquitous across vertebrates and invertebrates (Figure 1): starting from the single fertilized egg cell, first the morula formed by non-polarized cells turns into the blastula – a hollow sphere of cells with AB polarity. Then, as cells acquire additional PCP, primary head-tail axis forms and elongates during gastrulation and neurulation (Loh et al., 2016). Because of the optical transparency, these stages are particularly prominent in sea urchin. At the morula stage, a lumen in the center is formed and is gradually expanding as cells proliferate and rearrange into the hollow sphere. Next, during gastrulation, a group of cells invaginate and rearrange into a tube that narrows and elongates primarily by cell rearrangement and convergent extension movements (Martik and McClay, 2017). The tube then merges with the sphere at the side opposite to invagination, and as a result, the sphere transforms into a torus. Emerging data suggest that PCP drives both invagination and tube elongation (Nishimura et al., 2012; Croce et al., 2006; Long et al., 2015) – a recurring theme across species.

Figure 1 with 1 supplement see all
Two symmetry-breaking events, gain of apical-basal (AB) polarity and planar cell polarity (PCP), on cellular level coincide with the appearance of a rich set of morphologies.

Starting from an aggregate of non-polarized cells (globular symmetry), individual cells can gain AB polarity and form one or multiple lumens (spherical symmetry). Additional, gain of PCP allows for tube formation (axial symmetry). Complex morphologies can be formed by combining cells with none, one, or two polarities. In Figure 1—figure supplement 1, we schematically illustrate how existing models capture different elements of development.

https://doi.org/10.7554/eLife.38407.003

Mutations in PCP pathways produce shorter and wider tubes (Ochoa-Espinosa et al., 2012; Saburi et al., 2008; Kunimoto et al., 2017), somites (Song et al., 2010), and embryos (Gong et al., 2004). Formation and elongation of the tubes can proceed without cell division and cell death by cells rearranging along the tube’s axis, termed convergent extension (CE) (Andrew and Ewald, 2010; Tanimizu et al., 2009; Martik and McClay, 2017). While the importance of PCP in gastrulation and tubulogenesis is well established (Andrew and Ewald, 2010; Tanimizu et al., 2009; Martik and McClay, 2017; Kunimoto et al., 2017), it is unclear how polarity may control tube morphology.

The bulk-lumens-folds/tubes transition seen across animal species in early embryogenesis, is also a key feature of the later organ formation. The early stages of organogenesis in the liver, kidney, brain, gut, and pancreas are apparently so robust, that they can be recapitulated in vitro, allowing for advanced quantification and manipulation (Little, 2017). The case of pancreatic organoids is interesting as it illustrates an increase of morphological complexity from spheres to folds. Cells in in vitro pancreatic organoids first grow as a bulk and later acquire AB polarity and develop lumens. Depending on the growth conditions, organoids develop into a hollow sphere or acquire a complex folded shape (Greggio et al., 2013). It is currently unknown what drives the transition from sphere to folded state; the two possible hypotheses are rapid proliferation or physical pressure from growing into a stiff matrigel.

Is the apparent link between cellular polarity and morphological complexity accidental? Or, could it be that morphological transitions, stability, and diversity are emergent features in a population of proliferating polarized cells? If true, can we identify what drives the transition from lumens to folds and tubes? Why are these stable? Can we predict what controls fold depth, and tube length and width? To answer these questions, we lack a unified approach that could bridge polar interactions between single cells to the global features emerging on the scale of thousands of cells in 3D.

Starting with D’Arcy Thompson’s seminal contribution (Thompson, 1942), quantitative models aided in understanding specific morphogenetic events (Figure 1—figure supplement 1). Among these are invagination (Odell et al., 1981; Rauzi et al., 2015; Polyakov et al., 2014; Hočevar Brezavšček et al., 2012), primitive streak formation (Newman, 2008), convergent extension (Collinet et al., 2015; Belmonte et al., 2016), epithelial folding (Buske et al., 2012; Osterfield et al., 2013; Monier et al., 2015; Murisic et al., 2015), emergence of global PCP alignment from local cell–cell coupling (Amonlirdviman et al., 2005; Le Garrec et al., 2006; Burak and Shraiman, 2009), origins of tubulogenesis (Engelberg et al., 2008), and recently statistical properties of branching morphogenesis (Hannezo et al., 2017). However, they are often on either of the two ends of the spectra: those modeling single cells explicitly, often rely on vertex-based approaches and are limited to dozens of cells (Alt et al., 2017; Misra et al., 2016; Aigouy et al., 2010; Le Garrec et al., 2006). To capture the large features spanning thousands of cells, one typically turns to elastic models where AB polarity is implicit and epithelia is presented as a 2D elastic sheet (Hannezo et al., 2014; Etournay et al., 2015; Hufnagel et al., 2007; Nagai and Honda, 2009; Aliee et al., 2012; Nagai and Honda, 2001).

We developed a theoretical approach that, with only a few parameters, bridges cellular and organ scales by integrating both types of polarity. A main difference to earlier approaches is that a cell’s movement is coupled to how its AB polarity and PCP are oriented relative to each other and relative to neighbor cell polarities. In other words, in our approach, the adhesion strength between neighbor cells is modulated by the orientation of their polarities. We find that polarity enables complex shapes robust to noise but sensitive to changes in initial and boundary constraints, thus supporting that morphological stability and diversity are emergent properties of polarized cell populations. Lumens, folds, and stable tubes emerge as a result of energy minimization. We make testable predictions on morphological transitions in pancreatic organoids, tubulogenesis, and sea urchin gastrulation. Our approach illustrates the evolutionary flexibility in the regulatory proteins and networks, and suggests that despite differences in proteins between organisms, the same core principles may apply.

Model

There are three key elements that allow us to bridge the scale from cellular level to macroscopic stable morphologies.

(1) Cells are approximated by point particles

Cell–cell adhesion is modeled by repulsive and attractive forces acting between cell centers. This allows a substantially gain in computation time compared to vertex-based models where cell surface adhesion is explicitly considered (Alt et al., 2017). The potential for pairwise interaction between two interacting neighbors, i and j, separated by distance rij is

(1) Vij=e-rij- S e-rij/β,

where the first term corresponds to repulsion, and the second term to attraction (see Figure 2A). For a pair of non-polar cells, the strength of attraction S = 1. β >1 is the parameter that sets how much longer the attraction range is compared to repulsion. We set β = 5 throughout the paper, but our results and conclusions are consistent for smaller β. The main results are also not sensitive to the exact choice of the potential, thus for example the higher power in the exponential,

(2) Vij=e-rij4-S e-rij/β4,

give qualitatively similar results (see Figure 2—figure supplement 1A). The potential energy of a cell is the sum of pairwise neighbor interactions

(3) Vi= jVij.
Figure 2 with 4 supplements see all
Cells are modeled as interacting particles with a polarity-dependent potential.

(A) Potential between two interacting cells with apical-basal polarity (see Equation 6). Cells repulse when polarities are antiparallel (top/green part) and attract when they are parallel (orange/bottom part). (B–C) Two cells interact only if no other cells block the line of sight between them. (B) Cell i and j do not interact if ij’s midpoint (black dot) is inside of the Voronoi diagram for cell k (shaded in grey). (C) Cell i and j interact because cell k is further away than the distance ij/2 and ij’s midpoint therefore lie outside of cell k’s Voronoi diagram. In the related Figure 2—figure supplement 1A–D, we test the sensitivity of our model to the details of the potential and neighborhood assignments. In Figure 2—figure supplement 2, we relate changes in cell shapes to the model components, and in Figure 2—figure supplement 3 (Figure 2—video 1), we illustrate how altering polarity affects the dynamics of the systems with two and six cells.

https://doi.org/10.7554/eLife.38407.005

(2) Cells interact with (a subset of Voronoi) neighbors

Interacting neighbors of cell i are selected from a subset of cells sharing a Voronoi surface. The subset is limited to the nearest neighbors j which are closest to the midpoint between i and j (Figure 2B–C). This constrain effectively corrects for the finite volume associated with point particles and assures that two cells will not interact if the line of sight between their centers is separated by a surface of a third cell. Without this constraint, the macroscopic morphologies collapse. However, our results are robust to replacing the line of sight constraint with full Voronoi and a cut-off distance for attraction force (Figure 2—figure supplement 1B).

(3) Cell–cell adhesion depends on the orientation of polarity

To capture directional adhesion, we set the strength of attraction, S, to be dependent on the relative orientation of the polarities in each of the cells. We assume that AB polarity and PCP are orthogonal, and that the polarities of one cell align with the polarities of its neighbor cells. Mathematically, we introduce unit vectors, p^i and p^j representing AB polarity, and q^i and q^j representing PCP for cell i and j, respectively. We set

(4) S=λ1S1+λ2S2+λ3S3

where λ1, λ2, and λ3 are the strengths of the different polarity terms. We require that

(5) λ1+λ2+λ3=1,

to satisfy the constraint that perfectly aligned cells always have a steady state distance of 2 cell radii (for β = 5). To capture that in an epithelial sheet, AB polarities align parallel to each other and tight adherens junctions form in the plane perpendicular to AB polarity, we introduce the quadruple product

(6) S1=(p^i× r^ij)  (p^j× r^ij).

This makes two interacting cells with AB polarity maximally attracted (S1 = 1) if the two apical sides are next to each other. On the other hand, if apical side of one cell is next to the basal side of another cell, the two cells will be maximally repulsing (S1 = −1), see Figure 2A.

In case of planar polarization, we define

(7) S2=(p^i× q^i)  (p^j× q^j),

which makes the attraction maximal if the PCP of two cells are parallel to each other and perpendicular to their AB polarities. In addition, we assume that similarly to AB polarity, two cells with PCP are maximally attracted if their PCP are parallel and cells have the same kind of pole (e.g. Vangl-enriched) next to each other,

(8) S3=(q^i× r^ij)  (q^j× r^ij).

We later show that this assumption makes neighbor exchange on a sheet possible and results in CE. However, unlike with tight junctions, preferred directional adhesion with PCP is not as well established. Cells adhere to each other by membrane proteins assembled in adherens junctions just below the apical surface. Both proteins regulating adherens junctions, for example Smash (Beati et al., 2018), as well as adherence proteins forming adherens junctions can be planar polarized, for example Bazooka, E-cadherins (Simões et al., 2010; Tamada et al., 2017; Levayer and Lecuit, 2013; Warrington et al., 2013; Aigouy and Le Bivic, 2016). These indirectly support our assumption of anisotropic, planar polarized adhesion.

The motion of the cells and their polarities are calculated assuming overdamped dynamics

(9) dr¯idt=dVidr¯i+η,
(10) dp¯idt=dVidp¯i+η,

and

(11) dq¯idt=dVidq¯i+η,

where the p¯i and q¯i differentiation takes into account the rotation of polarity vectors, and η is a random uncorrelated Gaussian noise. In practice, we implemented the model in a MATLAB script (available in the Materials and methods section), where we use the Euler method. We perform the differentiation along the polarity by differentiating along all three cartesian coordinates (see Model details in the Materials and methods section). After each time step, we normalize the updated polarity vectors. The above differentiation does not include the change in partners when neighborhood changes. This is treated as a non-equilibrium step where potential energy can increase (Equation 3). Biologically, this is similar to cells spending biochemical energy as they rearrange their neighborhood.

The point particle approximation has been utilized earlier for modeling non-polar cell adhesion in early blastocyst (Krupinski et al., 2011), slug formation in amebae (Dallon and Othmer, 2004), and PCP organization in primitive streak formation (Newman, 2008). The main novelty of our approach is the dynamical coupling of cell positions and polarity orientations (Equation 6–8).

Model implications

One of the implications of the coupling between position and polarity is that in a sheet of cells, turning AB polarity in one cell will cause a force on its neighbors. In case of two cells (Figure 2—figure supplement 3 and Figure 2—video 1), the pair relaxes the imposed stress by rotating both the polarity and their positions. In biological terms, turning AB polarity in one cell (e.g. by apical constriction, illustrated in Figure 2—figure supplement 2) of an epithelial sheet will induce bending of the sheet as is the case with bottle cells in invagination.

The present formulation of PCP has several implications. First, we restrict the effects of PCP to directed (anisotropic) cell–cell adhesion and do not consider its other possible roles, in for example asymmetric cell differentiation, thus primarily focusing on its role on CE. Second, in our current formulation, AB polarity and PCP influence each other’s orientation on equal footing (Equation 7). PCP, however, is typically constrained to the apical plane and thus is expected not to influence the orientation of AB polarity. Disabling PCP’s effect on AB polarity (see Materials and methods) does not influence our main results on tube formation and gastrulation. However, the symmetry in polarities is appealing for its simplicity and is indirectly supported by the following experimental observations: First, cells can acquire PCP without AB polarity present (Baer et al., 2009; Zorn and Wells, 2009). Second, proteins required for AB polarity can be planar-polarized (Warrington et al., 2013; Aigouy and Le Bivic, 2016; Beati et al., 2018; Choi and Sokol, 2009; Dollar et al., 2005; Kaplan and Tolwinski, 2010). Third, changes in cell shapes during invagination (e.g. sliding of adherens junctions and formation of bottle cells) are regulated by PCP in neural tube closure (Ossipova et al., 2014; Nishimura et al., 2012; Kinoshita et al., 2008), gastrulation in C. elegans (Lee et al., 2006), sea urchin (Croce et al., 2006), and Xenopus (Choi and Sokol, 2009). These changes in cell shape effectively reorient AB polarity (Figure 2—figure supplement 2).

Results

We have recently introduced effective representation of AB polarity, and showed that it is sufficient for capturing spherical trophectoderm in the early blastocyst (Nissen et al., 2017). Expanding on that work, we here explore how AB polarity supports diverse yet stable and complex morphologies.

Stable complex shapes emerge from randomly polarized cell aggregates

Adult organismal shapes are stable over long time, maintaining sizes and relative positions of lumens and folds, despite continual local damages and cell renewal. To test if cellular polarity could enable such stability in time and to random local perturbations, we first performed a series of tests with AB polarized cells (Figure 3 and Figure 3—video 1).

Figure 3 with 3 supplements see all
Development of 8000 cells from a compact aggregate starting at time 0.

(A) Cells are assigned random apical-basal polarity directions and attract each other through polar interactions (see Equation 6). (A–D) Cross-section of the system at different time points with red and blue marking two opposite sides of the polar cells. Cells closest to the viewer are marked red/blue, whereas cells furthest away are yellow/white. (E) Full system at the time point shown in (D). (F) Development of the number of neighbors per cell (red) and the energy per cell (blue), as defined by the potential between neighbor cells in Figure 2. Dark colors show the mean over all cells while light-shaded regions show the cell–cell variations. The yellow dot marks the energy for a hollow sphere with the same number of cells. See Figure 3—video 1 for full time series. In Figure 2—figure supplement 1E–G and Figure 3—figure supplement 1, we study how the final morphology depends on noise. In Figure 3—figure supplement 2, we show how the outer surface self-seals, and that the shape is maintained when cells divide.

https://doi.org/10.7554/eLife.38407.010

When starting a bulk of cells with AB polarities pointing randomly, an initial rapid expansion (Figure 3A–C) stabilizes into a complex morphology of interconnected channels (Figure 3C–E). The shape remains unchanged for at least 10 times longer than the initial expansion (Figure 3C–E). The stability of the shape is illustrated by the time evolution of the average energy per cell (Figure 3F) that after an initial fast drop converges to a constant value. As expected, this value is higher than the energy of a hollow sphere (yellow dot in Figure 3F) – a configuration obtained if we start with radially, instead of random, polarized cells and preserve radial polarization at all times. The observed behavior is not sensitive to the shape of the potential (Figure 2—figure supplement 1A) but is sensitive to how the neighborhood is defined (Figure 2—figure supplement 1B–D). Rerunning the simulation in Figure 3 with different initial conditions results in a different stable shape (Figure 2—figure supplement 1E and Figure 3—figure supplement 1).

The macroscopic features of the shapes are robust to noise (Figure 2—figure supplement 1E–F and Figure 3—figure supplement 1). While the shapes emerging under high and low noise are not identical, the relative position and sizes of the majority of channels and lumens are preserved. The changes caused by noise stem from perturbations during initial expansion stage. If the same level of noise is applied after the system reached stable state, after time t = 10000, noise does not cause any major macroscopic changes (Figure 2—figure supplement 1G). The obtained shapes have self-sealing features, as an initial cut and unwrapping of a section of a surface refolds and seals back into the original morphology (Figure 3—figure supplement 2A–C). Furthermore, the shapes (Figure 3D–E) are also robust to overall growth (Figure 3—figure supplement 2D–F) retaining the same macroscopic features, just scaled to a larger size. Robustness to noise and cell proliferation further support the link between polarity and stability of morphologies, for example organ shapes, as they expand from infant to adult.

The final shapes are robust to noise but sensitive to initial and boundary conditions

The orientation of polarities in a subpopulation of cells may be set by the environment that the cells are embedded in, for example signaling molecules deposited into extracellular matrix can influence orientation of the AB polarity (Overeem et al., 2015) or signals from neighboring cells of a different type can orient PCP (Chu and Sokol, 2016). We will refer to these constraints as boundary conditions.

To investigate sensitivity to boundary conditions, we consider three cases where polarities are fixed at all times and point either radially out from the center of mass (Figure 4A), radially out from a central axis (Figure 4B), or pointing away from a central plane (Figure 4C). As anticipated, the difference in symmetries of boundary conditions results in a sphere, a cylinder, or two parallel planes (Figure 4—video 1). At the same time, in these symmetric cases, the differences in initial conditions but without imposed boundary conditions are not sufficient to generate different structures; they all converge to the nested ‘Russian doll’-like hollow spheres (Figure 4D). In development, this highlights the importance of the neighboring tissues for defining boundary conditions.

Figure 4 with 1 supplement see all
Different morphologies can be obtained by varying boundary conditions (Figure 4—video 1).

(A) A hollow sphere emerges if polarities are fixed and initially point radially out from the center of mass. (B) A hollow tube is obtained if polarities point radially out from a central axis. (C) Two flat planes pointing in opposite directions are obtained if polarities point away from a central plane. (D) For all three initial conditions (A–C), if the polarities are allowed to change dynamically and the noise is high (η = 100 compared to η = 10−1 in A–C), the resulting shape consists of three nested ‘Russian doll’-like hollow spheres that will never merge due to opposing polarities. In contrast to the random initial condition in Figure 3, the initial conditions in (D) are symmetric.

https://doi.org/10.7554/eLife.38407.014

Our results thus support the idea that polar adhesion enables stable and robust macroscopic shapes. The closest biological parallels would be the complex luminal morphologies emerging in reaggregation experiments on for example Hydra (Seybold et al., 2016) or in vitro culture of purjunkie brain cells (Muguruma et al., 2015). Together with our simulations, these experiments highlight how stable and complex morphologies can develop in non-proliferating populations from cell rearrangements alone.

Folding by pressure or rapid proliferation result in different fold-morphologies

Transitions from spheres to folded shapes are ubiquitous in development. Folds are an important part of in vivo organ development, and the composition of cell types in the folded organoids is closer to that in real organs (Greggio et al., 2013). To date, it is unclear what drives the transition from spheres to folded lumens. One possibility is that it is driven by the mechanical properties of the matrigel that effectively may place the growing organoid under pressure. Alternatively, data from 3D brain organoids suggests that the rapid cell proliferation leads to the emergence of surface folding (Li et al., 2017).

The simplicity of our tool allows to explore both of these scenarios. To model dividing cells, we pick a random cell from the entire population and introduce a new daughter cell with inherited polarity direction placed in a random location a half cell radius away from the mother cell. This event introduces dynamic perturbation by locally increasing cell density and requires some time to relax back to equilibrium. If proliferation is slow, and the time between two cell divisions anywhere in the system is longer than relaxation time of the whole system (the time it takes to reach equilibrium), the system approaches global equilibrium and will expand as a sphere. However, if proliferation is increased, the system will be pushed out of equilibrium and folds will emerge (Figure 5A, Materials and methods). In more quantitative terms, our forces are such that a single cell can move up to 0.2 cell diameter per time unit. For cells that divide every 1000 time units, the transition to non-equilibrium buckling happens when the system has grown to about 5000 cells (Figure 5—video 1).

Figure 5 with 2 supplements see all
The number of complex folds in a growing organoid depends on the generation time and the pressure from the surrounding medium (Figure 5—video 1).

(A) Number of local minima as a function of 1/(generation time), tG−1. In silico organoids grow from 200 cells up to 8000, 12,000, or 16,000 cells with different generation times and no outer pressure. (B) Number of local minima as a function of pressure, P. In silico organoids grow to the same size with the same 1/(generation time), tG−1 = 1.4⋅10−4 but different outer pressure. The images illustrate the 16,000 cells stage. Blue dots mark the average, while light shaded regions show the SEM based on triplicates. See also Figure 5—figure supplement 1 for additional measurements on the differences between rapid growth and pressure.

https://doi.org/10.7554/eLife.38407.016

As cells divide faster, our simulations predict a transition from a smooth spherical shell to an increasingly folded structure with multiple pronounced folds, in line with the observation of brain organoids proliferating at different rates (Li et al., 2017). In comparison with the model for cortical convolutions by Tallinen et al., 2016 in which folding is a result of expanding cortical sheet adhered to the non-expanding white matter core, our mechanism does not require a bulk core. Instead folds emerge in a fast-expanding sheet when the growth is faster than the global relaxation to dynamical equilibrium.

While we find that the external pressure is not necessary for folding, pressure alone can also drive folding (Figure 5B, Materials and methods). However, this scenario contradicts the observation that pancreatic organoids can grow as spheres or folded morphologies in gels with the same stiffness but different media composition (Greggio et al., 2013).

In principle, both scenarios may contribute to folding, but visually the fold morphologies are different. To differentiate between the two, we have quantified the final folded structures in terms of their local minima (Figure 5, Materials and methods, see also Figure 5—video 1). Our simulations predict that in the pressure-driven case, the number of local minima will reach an upper limit as organoids increase in size (Figure 5B). In the case of out-of-equilibrium proliferation, new folds can continue forming as organoids grow (Figure 5A). Increased proliferation causes more and shallower folds. These folds are different than obtained with pressure which causes fewer but deeper minima. Quantitatively, both the depth and the horizontal extension of the folds are about double as large with pressure than with growth-induced folding (Figure 5—figure supplement 1).

PCP enables convergent extension and robust tubulogenesis

Despite the numerous evidences supporting the role of PCP in tubulogenesis, it remains unresolved whether oriented cell division or the extent of CE controls tube length and width (Karner et al., 2009; Carroll and Yu, 2012). It is also debated if it is important for the tubes to maintain regular shape, or if it is only important for tube initiation and growth (Kunimoto et al., 2017).

The simplicity of our approach allows us to address these questions by introducing cell–cell interactions through PCP. This term favors front-rear cell alignment in the interaction potential with only two additional parameters: the strength of the orientational constraint of AB polarity with respect to PCP, λ2, and the strength of PCP, λ3 (see Equations 4–8). For simplicity, we focus on the stability (ability to maintain regular diameter over time) and tube morphogenesis in systems without cell division.

Inducing PCP in a spherical lumen leads to two significant events. First, independent of initial orientation, after some transient time PCPs becomes globally ordered and point in direction parallel to an emerging equator that self-organizes around the sphere (inset in Figure 6). This arrangement has the lowest energy. Second, cells start intercalating along the axis perpendicular to PCP orientation, gradually elongating the lumen (Figure 6—video 1). During intercalations, cells exchange their neighbors through T1-like transitions as reported experimentally, Figure 6—figure supplement 3 (Nishimura et al., 2012; Sanchez-Corrales et al., 2018). The intercalations along the axis continue until the force balance between AB polarity and PCP is restored at a new equilibrium. Thus, our model predicts that the strength of PCP (λ3) relative to AB polarity (λ1) determines the width and the length of the tube (Figure 6). We obtain similar results if we constrain PCP to always remain in the apical plane and thus does not allow PCP to reorient AB polarity (Figure 6—figure supplement 1, Materials and methods). Note, that this result is very different in nature from the tube presented in Figure 4B as both AB polarity and PCP can now reorient in each cell at any time.

Figure 6 with 4 supplements see all
The length and width of tubes are set by the strength of planar cell polarity (PCP, λ3).

For each value of λ3, we initialize 1000 cells on a hollow sphere with PCP whirling around an internal axis (PCP orientation marked by cyan arrows in the top-left inset). Semi-major axis (dark blue) and semi-minor axis (light blue) are measured at the final stage (Materials and methods). Images show the final state. Throughout the figure, λ2 = 0.5 and λ1 = 1 - λ2 - λ3. The animated evolution from sphere to tube is shown in Figure 6—video 1. See also Figure 6—figure supplement 1 where we show that tubes also form when we disable the direct influence of PCP on apical-basal polarity, and Figure 6—figure supplement 2 where we vary the degree of PCP along the axis of the tube. In Figure 6—figure supplement 3, we show that cell intercalations result in experimentally reported T1 neighbor exchanges during convergent extension.

https://doi.org/10.7554/eLife.38407.019

These results support the observations that stable tubes can emerge without cell proliferation. In addition, when first the tube is formed, loss of PCP does not lead to cyst formation as recently shown by Kunimoto et al. (2017). However, localized cysts could result if the lumen is initialized with varying strength of PCP along the axis (Figure 6—figure supplement 2).

Two polarities are sufficient to explain major features of sea urchin gastrulation

Currently invagination in neurulation and gastrulation is understood and quantitatively modeled as a process driven by changes in cell shapes or the mechanical properties of cells with AB polarity (Rauzi et al., 2013; Tamulonis et al., 2011; Misra et al., 2016; Hočevar Brezavšček et al., 2012). This process is often assumed to be driven by apical constriction and decoupled from the eventual tube formation and elongation. However, emerging data suggests that PCP drives both invagination and tube elongation (Nishimura et al., 2012; Croce et al., 2006; Long et al., 2015) likely because apical constriction is controlled by PCP (Ossipova et al., 2014; Nishimura et al., 2012; Croce et al., 2006) and cell intercalations, similar to those in CE, contribute to invagination (Nishimura et al., 2012; Rembold et al., 2006; Sanchez-Corrales et al., 2018).

To probe the limits of our approach, we investigated if AB polarity, PCP, and boundary conditions reminiscent of posterior organizer (Loh et al., 2016) are sufficient to recapitulate the main stages of sea urchin gastrulation: invagination, tube formation, elongation (by CE), and finally merging of the gastrula tube with the pole opposite to invagination site.

The current understanding is that sea urchin gastrulation consists of primary invagination, driven by swelling of the inner layer of extracellular matrix beneath invaginating cells (Lane et al., 1993), and formation of a ring of bottle cells due to apical constriction (Kimberly and Hardin, 1998), and secondary invagination where tube elongates due to CE (Lyons et al., 2012). PCP is necessary for both invagination, possibly through its effect on apical constriction in bottle cells (Nishimura et al., 2012), and tube extension (Croce et al., 2006).

Motivated by these observations, we set boundary condition such that PCP of the invaginating cells are oriented around the anterior-posterior (top-bottom) axis, and are always in the apical plane. This constrain on PCP orientation allows for CE. While this particular configuration is not documented, it is consistent with observed effects of WNT orienting PCP within the apical plane (Humphries and Mlodzik, 2018). Second, we simulate the combined effect of bottle cells (Figure 2—figure supplement 2) and bending by swelling of the extracellular matrix by applying an external force, F, on AB polarity (see Materials and methods). This force gradually reorients AB polarity away from the anterior-posterior axis, thus leading to bending of the epithelial sheet. The effect is maximal for cells closest to the anterior-posterior axis. This external force is a phenomenological description aiming at capturing the observed effects of how change in AB polarity results in tissue bending and does not aim at capturing mechanisms driving the reorientation of AB polarity.

As a result, cells start to rearrange, the bottom flattens (Figure 7B) and bends inward (Figure 7C). Subsequently, the CE-driven by PCP causes the invaginated cells to rearrange, tube elongates, and merges with the top of the sphere (Figure 7F and Figure 7—video 1). In line with experimental observations, the tube elongates due to cells moving in the tube (Martik and McClay, 2017).

Figure 7 with 2 supplements see all
External constraints on apical-basal (AB) polarity and planar cell polarity (PCP) can initiate invagination and drive gastrulation in sea urchin.

(A) The lower third of the cells in a blastula with AB polarity (apical is blue–white, basal is red–orange) pointing radially out acquire PCP (cyan–green) in apical plane pointing around the anterior-posterior (top-bottom) axis (as in the inset to Figure 6). (B) Flattening of the blastula and (C) invagination occur due to external force reorienting AB polarity (Materials and methods). (D–E) Tube elongation is due to PCP-driven convergent extension and (F) merging with the top of the blastula happens when the tube approaches the top. Throughout the simulation, λ1 = 0.5, λ2 = 0.4, and λ3 = 0.1 for the lower cells while the top cells have λ1 = 1 and λ2 = λ3 = 0. For full time dynamics see Figure 7—video 1. In Figure 7—figure supplement 1, we consider alternative scenarios of sea urchin gastrulation and and neurulation.

https://doi.org/10.7554/eLife.38407.024

The contribution of extracellular matrix swelling is specific to sea urchin invagination. As bottle cells alone are sufficient to drive invagination in other systems, we have tested a cell-intrinsic scenario where AB polarities of bottle cell neighbors prefer to be tilted towards each other (by e.g. modifying the potential in Equation 6 to capture AB polarities as in Figure 2—figure supplement 2, unpublished results). This also resulted in successful invagination.

Several observations suggest that radial intercalation movements towards the center of invagination may drive tissue bending (Sanchez-Corrales et al., 2018; Panousopoulou and Green, 2016; Rembold et al., 2006; Chung et al., 2017). When tested in our model, PCP alone (omitting external or cell-intrinsic reorientation of AB polarity) results in exvagination and tube elongating outside of the sphere (this is more energetically favorable than extending inwards). A similar, exogastrulating, phenotype was observed in a PCP mutant (Long et al., 2015). While the similarity may be accidental, it is possible that this mutation abrogates PCP-driven apical constriction (Ossipova et al., 2014; Nishimura et al., 2012) and since cellular apices face outside, abrogating apical constriction eliminates the bias in direction of tube formation. We have also tested the consequences of apical constriction (reorienting AB polarity) in the model without PCP and while this was sufficient for invagination, the cavity remained spherically symmetric and failed to form a tube.

Discussion

Despite the stunning diversity and complexity of morphologies, the same concepts seem to emerge across organismal development. One of them is the link between local, cellular, and global, organs/whole organism, symmetry breaking. We know, from experimental and, to a lesser degree, theoretical work that cellular polarity is essential for forming axis, complex folded sheets, and interconnected tubes (Figure 1—figure supplement 1). What we do not know is why are these shapes so stable, and where do the differences between species and organs come from. To understand the differences, we typically compare genes or gene regulatory networks, thus limiting our understanding to analog processes within a few related species.

Phenomenological description bridges cell polarity to macroscopic morphologies

To address the origins of morphological diversity and stability across species and organs, we focused on a phenomenological description of polarized cell–cell interactions. This allowed us to bridge local single-cell symmetry breaking events to global changes in morphologies spanning tens of thousands of interacting cells. With this tool at hand, we find that with only a few parameters, we can recapitulate the two global symmetry breaking events: formation of epithelial sheets and folds by cells with AB polarity, and emergence of global axial symmetry (tubes) among cells with PCP.

Remarkably, our results show that interactions among AB polarized cells lead to stable morphologies, that after initial relaxation remain indefinitely in their final configuration. The morphologies are robust to noise, growth, and local damage (Figure 3—figure supplement 2). These results may explain how organs and embryos preserve their architecture while growing. Polar cell–cell interactions not only provide clue to the morphological stability, but also point to a simple explanation to the origin of the diversity. We find that the exact morphological details are defined by initial conditions, for example initial positions and orientation of polarities, and boundary conditions, for example polarities restricted to certain direction for a fraction of the cells. It is thus tempting to speculate that diverse shapes do not require multiple interacting morphogen gradients, but can be a result of differences in initial and/or boundary conditions: as for example presence of yolk cells at start and boundary constraints by vitelline membrane (Schierenberg and Junkersdorf, 1992; Wu et al., 2010).

The diversity of shapes and forms is further enriched by a second symmetry breaking event, PCP, oriented perpendicular to AB polarity. Within our phenomenological framework, addition of PCP component is simple, and requires only two additional parameters: one favoring perpendicular orientation of AB polarity and PCP within a cell, and another, favoring parallel PCP alignments between neighbor cells. These constraints are the coarse-grained representation of the well-established experimental and computational results on intracellular symmetry breaking events and global ordering of planar polarities mediated by cell–cell coupling (Le Garrec et al., 2006; Amonlirdviman et al., 2005; Wang et al., 2006). The first constraint allowed formation of axial symmetry and in combination with AB polarity, stable tubes, with length and diameter remaining constant with time. The second constraint resulted in cell rearrangements and intercalations consistent with the cell-autonomous CE typically associated with PCP. The patterns of neighbor exchanges during cell-rearrangements are in line with the ubiquitous T1 exchanges through formation and resolution of four cell vertex, Figure 6—figure supplement 3 (Nishimura et al., 2012; Sanchez-Corrales et al., 2018). The mechanism of the CE in our model is in line with the results by ‘filopodia tension model’ where elongated structures of many cells emerge from local cell–cell interactions in a direction defined by PCP (Belmonte et al., 2016). The presented formulations of our model captures only some of the known events contributing to CE and does not include PCP-driven changes in cell shape mediated by for example apical constriction (Nishimura et al., 2012) or contribution of external forces (Lye and Sanson, 2011).

Combining AB polarization and a local induction of PCP in a subpopulation of cells was sufficient to obtain main stages of sea urchin gastrulation: invagination, tube formation, and elongation through CE as well as merging of the tube with the animal pole at the top of the blastula. It is important to notice that the model for gastrulation uses an external force (Equation 25) acting on AB polarity, and the model for neurulation uses an imposed external constraint on PCP. This is not fully satisfying, and suggests an extension of the model to capture tissue bending induced by changes in cell shapes at the edge of the region that invaginates.

The existing in silico models treat invagination and CE-driven tube elongation as independent processes (Figure 1—figure supplement 1). Recent data, however, suggests that multiple mechanisms (intracellular apical constriction, intercellular directed cell division and cell intercalations, and supracellular actomyosin cables) act simultaneously and contribute to both invagination and tubulogenesis (Chung et al., 2017; Nishimura et al., 2012; Ossipova et al., 2014). Within our approach apical constriction (modeled as reorientation of AB polarity) and CE can act in parallel (Figure 7). It will be interesting to parallel recent experimental work (Chung et al., 2017) and computationally investigate how a combination of intra-, inter- and supracellular mechanisms contribute to the robustness of tubulogenesis, and to what extent the model can capture the range of observed phenotypes.

It has been proposed that the above mechanisms may all be coordinated by PCP (Nishimura et al., 2012). Besides the reported molecular links, a simple logic suggests that these mechanisms cannot be isotropic as in this case the initial bending will result in spherical structures. Thus, apical constriction, cell intercalations, and actomyosin cables have to be anisotropic (planar polarized) in directions consistent with the eventual tube orientation. This anisotropy is reported for both ‘wrapping’ tubes forming parallel to the epithelial plane, for example neurulation (Nishimura et al., 2012), and ‘budding’ tubes forming orthogonally to the epithelial plane, for example salivary glands (Chung et al., 2017; Sanchez-Corrales et al., 2018).

As organizing signals such as WNT can induce and orient PCP (Chu and Sokol, 2016) within the apical plane, we asked if it is in principle possible to design PCP constraints (not limited to the apical plane) that would result in ‘wrapping’ and ‘budding’. First, pointing PCP out of the plane was sufficient for a sheet of cells to bend (Figure 7—figure supplement 1). This is because in our formulation, PCP can drive reorientation of AB polarity and that in its turn is able to bend the sheet (Figure 7). Both ‘budding’ and ‘wrapping’ were qualitatively captured by the model when the axial and radial anisotropy were set by constraining orientation of PCP for cells within a circle (‘budding’ in sea urchin) and two stripes of cells (mimicking hinge points in neurulation ‘wrapping’). While CE was needed for proper tube forming in sea urchin example, the tube formed without CE in neurulation (Figure 7—figure supplement 1).

Thus, the simulations suggest that ‘wrapping’ in neurulation and gastrulation in Drosophila (Figure 7—figure supplement 1) vs. ‘budding’ in sea urchin and organogenesis (Andrew and Ewald, 2010; Zegers, 2014) may be outcomes of different constraints imposed on PCP. While it is intriguing to speculate that PCP may be oriented out of epithelial plane directly by organizing signals, this may also be an indirect effect of a sequence of intermediate steps. As the organizing signals not only induce and orient PCP but also drive apical constriction (effectively reorienting AB polarity) the PCP may be gradually oriented out of the original plane of epithelium by the following sequence of events: PCP → apical constriction → tilt in AB polarity → tilt in apical plane → PCP out of original epithelial plane. This less precise but simpler interpretation highlights the fact that PCP may drive many of the alternative mechanisms of tubulogenesis and shifts the focus from the differences in mechanisms driving tubulogenesis to the differences in boundary conditions – a set of constraints imposed on cell polarities by neighboring tissue (e.g. notochord in neurulation and organizer in gastrulation).

Testable predictions

In addition to our conceptual findings, we propose three testable predictions. First, we predict that two potential mechanisms behind the emergence of folds in pancreatic organoids – matrigel resistance and rapid, out-of-equilibrium, cell proliferation – will result in distinct morphologies. Our results suggest that in case of rapid proliferation, the growing structure will develop many shallow folds close to the surface which later tend to deepen. In contrast, external pressure causes fewer but deeper and longer folds (Figure 5—figure supplement 1). And further, as organoids grow in size, the number of folds will reach an upper limit when under pressure, however, in case of rapid proliferation, the number of folds will keep growing (Figure 5). Visual inspection of published morphologies seems to support the out-of-equilibrium growth (Greggio et al., 2013; Li et al., 2017). To assess if the growth is out-of-equilibrium in 3D organoids, one can quantify the distributions of cell shapes (Cerruti et al., 2013). Our model thus predicts that quantitative counting of folds and measurements of the fold depth and length relative to the size of the growing structure may discriminate between the alternative hypothesis. Quantification of the folds can be done in in vitro organoids by either phase or confocal fluorescence microscopy of whole-mount immunostained samples (Greggio et al., 2013; Li et al., 2017). The fold depth and length can be quantified with the same approach as applied to simulated shapes (Materials and methods) in binarized images of the 3D organoid surfaces. As oganoids cannot be cultured without gel supporting 3D growth, it will be necessary to vary both gel stiffness and generation time to uncouple their respective contribution to the folding. The work by Little, 2017 shows that in brain organoids generation time can be both slowed down and speeded up by either genetic manipulation or by adding small molecule inhibitors of pathways regulating cell proliferation. Unfortunately, changing stiffness of matrigel also changes its biochemical composition and may affect cell proliferation and differentiation. One will have to turn to synthetic hydrogels, where it is now possible to uncouple mechanical and biochemical clues (Gjorevski et al., 2016). To illustrate possible applications of our approach, we have only focused on two out of several possible mechanisms that may contribute to folding. The other likely alternative is that folding may result from differences in biomechanical interactions or generation times characteristic to the different cell types. These alternative scenarios are straightforward to consider in our model and will be an exciting venue to explore when more quantitative data on differences in organoid morphologies is available.

Our second prediction is that in case of tubes formed by non-proliferating cells, the length and width of the tubes are controlled by the relative strength of AB polarity and PCP. This result calls for quantification of adhesion proteins along the AB polarity and PCP axes. In PCP mutants with shorter and wider tubes, one would expect less planar polarization in adherens junctions and actomyosin, for example larger spread compared to wild type in their orientation quantified relative to the tube axis (Nishimura et al., 2012). Alternatively, the balance between AB polarity and PCP can be altered by weakening AB polarity, for example mutating tight junction proteins should result in longer tubes. A similar phenotype has already been reported for Drosophila tracheal tube (Laprise et al., 2010). With the recent advancements in in vitro systems of tubulogenesis, allowing for easy genetic manipulations and more amenable for quantitative imaging, it may in principle be possible to relate the extent of planar anisotropy in PCP mutants and strength of AB adhesion in tight junction mutants with tube length and diameter. The existing coupling between PCP and AB polarity may, however, make it challenging to tweak one polarity at a time.

Our third, and probably most challenging to test, prediction is on the conditions differentiating between tubes forming perpendicular (e.g. sea urchin gastrulation) or parallel (as in Drosophila gastrulation or neurulation) to the plane of epithelium. We predict that the outcome will be defined by the orientation of PCP in the invaginating region and the geometry of the boundaries (circular for budding and axial for wrapping) set by for example WNT organizing signals. Recent development in imaging localization of PCP complexes in single cells (Wu et al., 2013; Chu and Sokol, 2016; Minegishi et al., 2017; Habib et al., 2013) allows monitoring localization of PCP complexes, and thus PCP orientation, in individual cells. By placing WNT-soaked (Habib et al., 2013) beads or WNT-secreting cells (Chu and Sokol, 2016) one can vary PCP orientation in the cells at the epithelial boundary facing WNT and test for the direction of the epithelial bending and possibly tube formation.

Our approach is by purpose phenomenological and by its nature cannot make predictions about specific molecular details. In all cases, we do not see our simulations as finalized predictions, but rather as pointing in the most promising direction for further exploration of these complex developmental processes. Our setup easily allows for changes as we learn more. The proposed tool should be used in close collaboration with gained experimental knowledge on initial conditions, cell generation times and differentiation processes where polarities play a central role.

Our results open for a series of biological generalizations both in development and diseases. On one hand, we now may be able to explain and unify the apparently very distinct morphological transitions during gastrulation in flies, frogs, fish, mice, and humans by accounting for different initial and boundary conditions. Our model suggests how a moderate change in expression of polarities during some critical evolutionary stages could lead to widely different final morphologies. Thereby, development driven by cell–cell polarity interactions could provide major morphological transitions from local and transient modulations in polarity.

On the other hand, it becomes possible to think of gastrulation, neurulation, tubulogenesis, and organogenesis as the same class of phenomena, where the orientation of the tube is guided by local organizers, and lengths/widths of the tubes are determined by the relative strength of AB polarity and PCP. At the same time, there is an emerging view that wound healing and cancer are local perturbations – for example local loss of cells, dysregulation of cell polarities (Martin-Belmonte and Perez-Moreno, 2012), proliferation, or autonomously induced organizing signals – of otherwise conceptually the same developmental processes (Humphries and Mlodzik, 2018). The power of our model is that it allows to address these hypotheses through predictive models for the dynamics of many cells that interact through combinations of AB polarity and PCP.

Materials and methods

Throughout the paper, we use the Euler method to integrate the ordinary differential equations stated in the Model section with dt set to 0.1 or 0.2. Lower dt values will give qualitatively similar results but with increased simulation time. Higher dt values, will result in a collapse of the presented morphologies. The time unit is arbitrary, and the same throughout the paper.

The noise parameter η = 10−4 where nothing else is stated. Lower noise values will give more smooth simulations, while η on the order of 100 will result in collapsing shapes (see also Figure 2—figure supplement 1E–F).

Model details

In our model, we use the following potential to describe the pairwise interaction between cells

(12) Vij=e-rij-S e-rij/β,

 where rij is the center–center distance between cell i and cell j, and S is the polarity factor

(13) S=λ1S1+λ2S2+λ3S3

Here, λ1, λ2, and λ3 are the strengths of the respective polarity terms which are given as

(14) S1=(p^i× r^ij)  (p^j× r^ij),
(15) S2=(p^i× q^i)  (p^j× q^j),
(16) S3=(q^i× r^ij)  (q^j× r^ij).

The unit vectors p^i,p^j and q^i,q^j represent the apical-basal polarity and planar cell polarity of cell i and j. Throughout the paper, β is a constant which we set to 5. In order to use the Euler method, we need the gradient of Vij differentiated with respect to position, r¯i, and the two polarities, p¯i and q¯i:

(17) dVijdr¯i=erij/β{ γ r^ijλ1rij[(r^ijp^j)p^i+(r^ijp^i)p^j]λ3rij[(r^ijq^j)q^i+(r^ijq^i)q^j]},
(18) dVijdp¯i=erij/β{ λ1[S1p^ip^j+(r^ijp^j)r^ij]+λ2[S2p^i(q^iq^j)p^j+(q^ip^j)q^j]},
(19) dVijdq¯i=erij/β{ λ2[S2q^i(p^ip^j)q^j+(p^iq^j)p^j]+λ3[S3q^iq^j+(r^ijq^j)r^ij]}.

In order to derive Equations 17, 18, and 19, we have used the following:

(20) ddr¯ierij/β=1β erij/β r^ij,
(21) γ=erij(β1)/βSβ+2rij[λ1(r^ijp^i)(r^ijp^j)+λ3(r^ijq^i)(r^ijq^j)],
(22) dS1dr¯i=1rij[(r^ijp^i)p^j+(r^ijp^j)p^i]2[(r^ijp^i)(r^ijp^j)r^ij],
(23) dS1dp¯i=1pi[p^j(r^ijp^j)r^ijS1p^i],

where pi is the length of the polarity of cell i which is equal to one at all times.

Generation time in growing organoids

In Figure 5, the number of cells, N, at a given time, t, is N = 200 exp(ln(2) t/tG) where tG is the generation time. In these simulations, the AB potential (Equation 6) between cells is set to zero when the angle between their polarity is larger than π/2.

Modeling resistance from matrigel

In Figure 5, we model resistance from the matrigel by imposing a surface force pointing toward the center of mass. The potential of the pressure in the growth medium is given by

(24) VM=Pr22rmax

where P is the stiffness of the medium, r is distance from the center of mass, and rmax is the distance to the cell that is the furthest away from the center of mass. The resulting force will be constant in time at the periphery. Thus, all cells on a growing sphere will be exposed to a force of equal size. However, cells that end up deep inside a folded morphology will experience weaker resistance.

Quantification of the local minima

In Figure 5, the number of local minima is defined as the number of cells that do not have any neighbor cells that are closer to the center of mass than themselves, and at the same time have an average angle between their AB polarity and their neighbor cells displacement vector that is less than π/2.

Measuring the tube length and width

In Figure 6, the semi-minor and semi-major axes correspond to the half-width and half-length of the tubes, respectively. As cells on opposite sides of a tube have AB polarity pointing in opposite directions, we approximate the semi-major and semi-minor axes, by finding the half of the maximum and minimum distance between two cells with AB polarity pointing in opposite directions.

Modeling cells with different polarities

In the gastrulation simulation (Figure 7), each cell is assigned a specific value of polarity strengths (λ1,i, λ2,i, and λ3,i). We define the mutual interaction strength between a pair, i and j, of cells in Equation 4 with different polarity strengths by setting λ1 = mean(λ1,i, λ1,j), and λ2 = mean(λ2,i, λ2,j) as well as λ3 = mean(λ3,i, λ3,j). This choice makes sure that two neighbor cells interact with a force with equal magnitude but opposite sign. Furthermore, it makes sure that λ1 + λ2 + λ3 = 1 holds for all cells.

Disabling PCP’s effect on AB polarity

In the model described by Equations 6–9, AB polarity and PCP can influence each other’s orientation. To constrain PCP to the apical plane and thus disable its influence on AB polarity, we set λ2 = 0 when updating AB polarity in time (during the numerical integration of Equation 10).

Modeling bottle cells and apical constriction by an external force

The invagination in gastrulation is implemented by adding an external force, F, that act on the AB polarity in addition to our usual intrinsic forces from Equation 10

(25) F=-k r^ e(-x2-y2)/σ2,

The two parameters are k which is the strength of the force (in Figure 7, k = 0.02), and σ which defines the decline of the gaussian force (in Figure 7, σ = 10). r^is the unit vector pointing from x = y = 0 to the cells’ position, and x and y are the respective coordinates of the cells. This force is applied on the AB polarity and bends the orientation of the polarity away from a z-axis (the anterior-posterior axis).

MATLAB script

MATLAB script to generate and visualize data. MATLAB R2016b or newer is required together with the Statistics and Machine Learning Toolbox. In addition, the Parallel Computing Toolbox is required if the PAR parameter in the basic script is set to 1. The input folder contains initial conditions for three standardized systems. Bulk systems have neither apical-basal (AB) polarity nor planar cell polarity. Plane and shell systems have only AB polarity. ‘N’ in the file names gives the number of cells in the system. The initial polarity directions can be modified on line 4–5 in the basic.m file. Inside this file, it is also possible to set the degree of noise (η), the size of the time steps (dt), and the relation between the polarity strengths (λ1, λ2, and λ3). The parameter inc is used to speed up the simulations by only applying the neighborhood function to the nearest 100 neighbors. Generated data is saved in the output folder, and the visualization script is in a separate folder.

References

  1. 1
  2. 2
  3. 3
  4. 4
    Vertex models: from cell mechanics to tissue morphogenesis
    1. S Alt
    2. P Ganguly
    3. G Salbreux
    (2017)
    Philosophical Transactions of the Royal Society B: Biological Sciences 372:20150520.
    https://doi.org/10.1098/rstb.2015.0520
  5. 5
  6. 6
  7. 7
  8. 8
  9. 9
  10. 10
  11. 11
  12. 12
    Chapter Eight - The Kidney and Planar Cell Polarity.”
    1. TJ Carroll
    2. J Yu
    (2012)
    In: Y. ingzi Yang, editors. Current Topics in Developmental Biology. Academic Press. pp. 185–212.
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
    Regenerative capacity in newts is not altered by repeated regeneration and ageing
    1. G Eguchi
    2. Y Eguchi
    3. K Nakamura
    4. MC Yadav
    5. JL Millán
    6. PA Tsonis
    (2011)
    Nature Communications, 2, 10.1038/ncomms1389.
  23. 23
    In Silico Simulation of Epithelial Cell Tubulogenesis
    1. JA Engelberg
    2. M Kim
    3. KE Mostov
    4. C Anthony Hunt
    (2008)
    30th Annual International Conference of the IEEE Engineering in Medicine and Biology Society. pp. 1036–1039.
  24. 24
  25. 25
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
  42. 42
    A Role for Regulated Secretion of Apical Extracellular Matrix during Epithelial Invagination in the Sea Urchin
    1. MC Lane
    2. MA Koehl
    3. F Wilt
    4. R Keller
    (1993)
    Development 117:1049–1060.
  43. 43
  44. 44
  45. 45
  46. 46
  47. 47
  48. 48
  49. 49
  50. 50
  51. 51
  52. 52
  53. 53
    "Chapter Five - Tension and Epithelial Morphogenesis in Drosophila Early Embryos"
    1. CM Lye
    2. B Sanson
    (2011)
    In: Michel Labouesse, editors. Current Topics in Developmental Biology, 95. Academic Press. pp. 145–187.
  54. 54
  55. 55
  56. 56
  57. 57
  58. 58
  59. 59
  60. 60
  61. 61
  62. 62
  63. 63
  64. 64
    Computer simulation of wound closure in epithelial tissues: Cell–basal-lamina adhesion
    1. T Nagai
    2. H Honda
    (2009)
    Physical Review E, 80, 10.1103/PhysRevE.80.061903.
  65. 65
  66. 66
  67. 67
  68. 68
  69. 69
  70. 70
    The mechanical basis of morphogenesis. I. Epithelial folding and invagination
    1. GM Odell
    2. G Oster
    3. P Alberch
    4. B Burnside
    (1981)
    Developmental Biology 85:446–462.
  71. 71
  72. 72
  73. 73
  74. 74
  75. 75
  76. 76
  77. 77
    Embryo-scale tissue mechanics during Drosophila gastrulation movements
    1. M Rauzi
    2. U Krzic
    3. TE Saunders
    4. M Krajnc
    5. P Ziherl
    6. L Hufnagel
    7. M Leptin
    (2015)
    Nature Communications, 6, 10.1038/ncomms9677.
  78. 78
  79. 79
  80. 80
  81. 81
  82. 82
  83. 83
  84. 84
  85. 85
  86. 86
    On the growth and form of cortical convolutions
    1. T Tallinen
    2. JY Chung
    3. F Rousseau
    4. N Girard
    5. J Lefèvre
    6. L Mahadevan
    (2016)
    Nature Physics 12588:588–593.
  87. 87
  88. 88
  89. 89
  90. 90
    On Growth and Form
    1. DW Thompson
    (1942)
    Cambridge: Cambridge University Press.
  91. 91
  92. 92
  93. 93
  94. 94
  95. 95
  96. 96

Decision letter

  1. Aleksandra M Walczak
    Reviewing Editor; École Normale Supérieure, France
  2. Naama Barkai
    Senior Editor; Weizmann Institute of Science, Israel

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

Thank you for submitting your article "Theory bridging cell polarities with development of robust complex morphologies" for consideration by eLife. Your article has been reviewed by Naama Barkai as the Senior Editor, a Reviewing Editor, and three reviewers. The reviewers have opted to remain anonymous.

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.

The reviewers and Editor greatly appreciated the novelty of your theory that will hopefully stimulate discussion, experiments and new work. They found the paper a creative, and thought-provoking contribution to the field of development. It puts the emphasis on local and rather generic interactions between cells' polarities. It is a worthwhile goal to illustrate the range of morphogenetic events that could occur based on the assumed cellular behaviors, even if they reflect biological reality to greater or lesser extents. Nevertheless, the reviewers raise a number of concerns, summarised below. Please consider all of the reviewers' points. The reviewers often suggest presentation ideas since they do think the work is worthwhile, even if it is sometimes controversial it offers a good starting point for new models and discussion.

- The reviewers think, if the model is to inform our understanding of morphogenesis, that the mechanical forces in the model reflect, at least loosely, a physically plausible state of forces in the tissue. This is a very serious concern. Specifically, reviewer 1 raises serious concerns about the origins of tissue rotation that absolutely needs to be addressed. reviewer 2, from a very different perspective, seconds these concerns, worrying about the ability of PCP axes to influence AB axes. These points are better explained in the reviewer comments.

- All reviewers raised concerns about the sea urchin example and it would probably be useful to find a good way of either discussing this example. One suggestion is to use the current manuscript to show what is possible. The sea urchin example is (probably – unless we are misunderstanding things) implausible, but is a possible outcome given the input assumptions, and we'd therefore suggest not to eliminate it but discuss the plausibility. In the experts’ opinion, the initial conditions for neurulation are equally implausible, but still demonstrate the possibilities of the model.

- In short, a clear description of assumptions and their influence on the results would be much appreciated.

- In general, as you can easily see some of the reviewers are physicists, and some are biologists. It is a hard exercise but please try and rewrite the manuscript to appeal to both (e.g. define terms that are obvious to physicists, be critical in discussing biological reality for the ignorant physicist), as it is clear it does have the potential to influence both groups of researchers. There is no need to oversell.

- A discussion of testable predictions needs to be improved by making more explicit what manipulations would be performed and how outcomes would be analyzed.

Reviewer #1:

This manuscript presents a cell-based model for the morphogenesis of 3D epithelial structures. Cells, described as point particles, are endowed with two-unit vectors representing apico-basal (AB) polarity and planar cell polarity (PCP). Cells interact through a pair potential that is modulated according to the orientation of the polarity vectors. The positions of the cells and the orientations of the polarity vectors follow a simple relaxational dynamics, with a stochastic term.

The introduction of a minimal model to study the collective dynamics of polarized cells, in the context of epithelial morphogenesis, is of potential interest. However, to be of value, such a model should present a physically consistent description of the phenomena of interest. I am unclear that this is the case of the present model. Although its energy function is invariant under simultaneous rotation of the cell positions and polarity vectors, the dynamics as written here (below equation 8) can give rise to a net torque and rotation of the system, which should not be if motion is driven by internal forces. To rephrase my concern, if I understand the model correctly, the orientation-dependence of the energy means that the tissue as whole can rotate to orient in a specific direction, like a compass in a magnetic field (I think that's part of how flat sheets are made in one example). But for a freely-floating tissue there isn't an extrinsic force like the magnetic field to rotate it around, and it's not a motion that can occur through intercellular forces.

This is particularly evident when the polarities are fixed: with a fixed AB polarity alone (as in Figure 4), a doublet of cells minimizing its energy would rotate to align perpendicular to their polarity vectors (cf. the coupling term in equation 6). As far as I understand it, the same effect is at play in the application of the model to sea urchin gastrulation (Figure 7): the imposition of a fixed PCP orientation, through a displacement of AB polarity, gives rise to an extrinsic torque that bends the tissue. If this interpretation of the model is correct, if there's an extrinsic force imposing the orientation of the surface at every point, the model is not a serious contender to explain how an urchin is made. It's maybe less critical for other examples but leaves me unclear how seriously the possibilities illustrated by the model can be taken. That being said, the model certainly stands as a stimulant for further thinking along the same lines.

Similarly, in the case of fold formation on a sphere (Figure 5), which the authors suggest can occur if cell proliferation is on a faster time scale than relaxation, the physical origin of a slow relaxation time scale, and the mechanism of the observed instability, would merit discussion. Although assuming relaxational dynamics is standard and justified when the local minima of the energy are of interest, it might be questioned when turning to dynamical phenomena. Do folds develop through buckling, because the sphere grows (in number of cells) faster than it can expand in the surrounding medium? In many contexts, mechanical equilibration is much faster (seconds or minutes) than proliferation (hours). What would be the origin of a longer time scale for relaxation? Resistance from the surrounding medium might be one, but the authors present fast proliferation as an alternative mechanism to pressure from growth in a stiff gel.

I also have some issues of form with the manuscript (e.g. I do not find the connection that the authors draw between the model and evolution very substantive), but the above questions on the physical relevance of the model are the main issues that I would like to see clarified before the manuscript is given further consideration.

Reviewer #2:

This manuscript describes a mathematical model and associated simulations of large clusters of cells that can be assigned properties of apical basal polarity and planar cell polarity and may proliferate, to follow their 3D morphological evolution. The aim is to test the hypothesis that "stability, complexity, and diversity are emergent properties in populations of proliferating polarized cells." The authors conclude that "polarity enables complex shapes robust to noise but sensitive to changes in initial and boundary constraints." The work is of potential interest, but a number of substantial issues should be addressed. Most important, some essential definitions are not provided, the writing needs to be made more accessible to biologists who are not modelers, some key assumptions should be made explicit, and others.

Essential revisions:

Definitions needed: What precisely is meant by "stability?" Are the authors referring to feature constancy over time? Do they mean that models will converge to a specified set of features? What features? For example, they suggest that the trajectory in Figure 3 is stable, but noise seems to qualitatively change the shape (Figure 2—figure supplement 1), though it still produces a complex folded shape. It's not clear what, exactly, the authors are claiming.

I struggled with the idea of "boundary conditions." I eventually concluded that, here, the term appears to mean restrictions on polarity orientations, and maybe other parameters. However, to a biologist, "boundary conditions" means something very different. Difficulty with this passage was compounded by the challenging sentence "At the same time, in these symmetric cases the differences in initial conditions (no boundary conditions) are not sufficient to generate different structures." This implies that differences in initial conditions means no boundary conditions. I think the authors mean to say "differences in initial conditions but without imposed boundary conditions." Finally, in the Discussion section the authors say, "boundary conditions, e.g. fixed polarities for a fraction of the cells in some time." A precise definition should be provided when the idea of boundary conditions is first introduced.

"Minima" is defined in the Materials and methods section, but the authors also use "maxima." By analogy, one can guess what this means, but please define.

Improve writing to make accessible to biologists: I'm a biologist with ample experience collaborating with modelers. Still, I found the manuscript difficult to follow and needed several readings to understand. If the goal is to make this work accessible to a broad audience (and it should be if published in eLife), the authors should carefully revise the text to replace jargon with precise simple language.

Describe implications of model structure: Several features of the model are not as explicit as they should be. One is the conflation of PCP with convergent extension behavior. PCP defines a vectorial molecular asymmetry, and PCP can have a number of outcomes. One is CE, but others include the asymmetric elaboration of cellular structures, or the asymmetric induction of gene expression that results in different cell fates, for example. Here, the model structure equates PCP to forces that drive CE. Furthermore, the descriptions of this implementation are misleading. In the Introduction, the authors, in describing the implementation of AB and PCP polarity, say "In other words, the adhesion strength between neighbor cells is modulated by the orientation of their polarities." Are the authors describing how AB and PCP (=CE) are implemented in the model, or how they function biologically? There is no evidence that molecular PCP functions by regulation of adhesion (involvement of Celsr proteins notwithstanding), and differential adhesion is only one part of a complex and only partially understood set of mechanisms that drive CE.

Another implicit part of the model as I understand it, and one that is probably very important to the outcomes, is the notion that altering PCP vectors can influence the AB vectors of cells. I know of no biological precedent to justify this behavior. This input assumption is inconsistent with my understanding of the biological reality. AB polarity is an initial condition that must be met for the subsequent establishment of PCP, which is then by definition orthogonal to the AB axis. There are examples in which PCP determines edges of groups of cells that then increase tension, thereby inducing folds, but this is a different mechanism from the one implied here. I think it would be interesting for the authors to explore how their model would operate if AB polarities orient PCP, but PCP cannot directly influence AB. I suspect that the examples of tube formation would fail. At a minimum, they should discuss the need for this assumption explicitly, noting that it is not currently supported by biological data (or citing the data if there is some that I'm not aware of).

Some unlikely assumptions: While it may be interesting that complex forms such as the sea urchin-like gastrulation and the folded neural tube CAN be produced within the confines of the model, the conditions needed stretch the bounds of plausibility. The model simulating sea urchin gastrulation requires a precise choreography of PCP parameters that could happen but seem unlikely. Even more difficult to believe is the initial conditions invoked to simulate neural plate folding. Specifically, invoking "two rows where the PCP points parallel and antiparallel to AB polarity" invokes conditions for which I know of no precedent at all.

Additional discussion of testable predictions needed: Re-emergence of folds in pancreatic organoids, pressure causes fewer and much deeper folds that form early during growth. This seems to be a qualitative measure. How would biological observations be assigned to one or the other category? Isn't the key difference that the pressure hypothesis predicts that the number of minima maxes out, whereas the proliferation model predicts that the number keeps growing with increase in size (Figure 5)?

Re tubes formed by nonproliferating cells, "in PCP mutants with shorter and wider tubes one would expect higher fraction of "AB adhesion" relative to "PCP adhesion". However, as noted above, PCP is not directly related to any known adhesion events.

Re the conditions differentiating between tubes forming perpendicular (e.g. sea urchin gastrulation) or parallel (as in Drosophila gastrulation or neurulation) to the plane of epithelium, the measurements could in principle be made, but as noted above, the conditions seem to me to be highly implausible, and ask us to ignore the likely contribution of other signals.

Reviewer #3:

This is a creative and innovative article on development. The authors introduce a simple theoretical model based on the idea that apical-basal (AB) polarity and planar cell polarity (PCP) are (the) main driving forces relevant for the formation of complex tissue morphologies. The model idealises cells as point particles with two independent axes (polarities), an AB and a PCP axis. Pairwise interactions between cells contain an attractive and a repulsive component where the magnitude of the attractive component depends on three additive factors: alignment of the AB axis, alignment between AB and PCP axis, and alignment of PCP axis. Noise comes into the model at various levels: Langevin noise in the equations of motion, noise at the level of choosing the axes orientation at cell proliferation, and noise at the level of initial conditions. In summary, this seems to be a computationally efficient model to explore the question: How can cell-cell adhesion that depends on the relative orientation of polarisation axes alone drive different morphologies during development?

Indeed, the authors simulations show that their model is able to reproduce a broad range of morphologies depending on several factors: choice of initial conditions, relative tendencies of the various axis to orient with respect to each other, and cell proliferation rate. A clear strength of the computational model (I would not call it a theory) is that it provides a mechanistic link between morphology and cell-cell interaction as well as initial conditions. Thus, I find it a clear conceptual advance in the field which puts a fresh perspective on the role of generic interactions for the emergence of different tissue morphologies. The authors also list a few testable predictions. Here, I am not so convinced yet for very basic reasons: Testable predictions require either a clear link between a molecular and an emergent process, or a set of mutually exclusive alternatives that can be tested. Concerning the latter case (in the section of organoids) the authors analyse two possible scenarios (cell proliferation and a specific way to implement external forces). These are not two mutually exclusive alternatives but only two out of many possible conceivable scenarios. Concerning the former case, in the example on sea urchin gastrulation the authors seems to a large degree put in what they would like to get out. Finally, and this is more on the side of a personal opinion, I am not convinced that morphologies that strongly depend on initial conditions are a robust way for cellular development. For instance, how can an organism make sure that all these initial conditions are properly chosen?

In summary, I found the paper a well-written, creative, and thought-provoking contribution to the field of development. It puts the emphasis on local and rather generic interactions between cells' polarities. If the authors can convincingly reply to my above queries, I would be inclined to recommend it for publication.

https://doi.org/10.7554/eLife.38407.032

Author response

Reviewer #1:

This manuscript presents a cell-based model for the morphogenesis of 3D epithelial structures. Cells, described as point particles, are endowed with two-unit vectors representing apico-basal (AB) polarity and planar cell polarity (PCP). Cells interact through a pair potential that is modulated according to the orientation of the polarity vectors. The positions of the cells and the orientations of the polarity vectors follow a simple relaxational dynamics, with a stochastic term.

The introduction of a minimal model to study the collective dynamics of polarized cells, in the context of epithelial morphogenesis, is of potential interest. However, to be of value, such a model should present a physically consistent description of the phenomena of interest. I am unclear that this is the case of the present model. Although its energy function is invariant under simultaneous rotation of the cell positions and polarity vectors, the dynamics as written here (below equation 8) can give rise to a net torque and rotation of the system, which should not be if motion is driven by internal forces.

One of the implications of the coupling between position and polarity is that the rotation of the system can indeed result in absence of external forces on cell position, but if external force is applied to cell polarity (e.g morphogen gradients or change in cell shape, illustrated in Figure 2—figure supplement 2). In a sheet of cells, turning AB axis in one cell would result in rotation of its neighbors Figure 2—figure supplement 3 (and Figure 2—video 1). In case of two cells (Figure 2—figure supplement 3C–D and Figure 2—video 1C–D), a somewhat similar physical system would be two spheres connected by a rod. A finite rotation of one sphere (without change of its position) will rotate the rod and thus displace the other sphere.

We have now added this implication to our subsection ”Model implication”.

To rephrase my concern, if I understand the model correctly, the orientation-dependence of the energy means that the tissue as whole can rotate to orient in a specific direction, like a compass in a magnetic field (I think that's part of how flat sheets are made in one example). But for a freely-floating tissue there isn't an extrinsic force like the magnetic field to rotate it around, and it's not a motion that can occur through intercellular forces.

Author response image 1
Changing the polarized direction of a plane of cells does not rotate the plane as a whole but breaks it into smaller planes.

Time t = 0 shows a plane consisting of 500 cells with AB polarity pointing to the right. At time t = 0.1, the direction of the polarity is shifted by 45 degrees. Since the polarity is fixed in time, the planes break into smaller pieces, and at time 104 they have merged into two planes that are separated by several cell diameters.

https://doi.org/10.7554/eLife.38407.029

To parallel the example with a plane in a magnetic field, we have simulated a plane of 500 cells in which AB polarity is turned by 45 degrees and maintained at that new angle at all times. Here there is no external field on cell positions in our potential but a constant external constraint on the AB direction. Then the subsequent relaxation dynamics will move cells, disrupting the epithelial sheet in many tilted subplanes. Thus the sheet does not rotate as a connected structure.

This is particularly evident when the polarities are fixed: with a fixed AB polarity alone (as in Figure 4), a doublet of cells minimizing its energy would rotate to align perpendicular to their polarity vectors (cf. the coupling term in equation 6). As far as I understand it, the same effect is at play in the application of the model to sea urchin gastrulation (Figure 7): the imposition of a fixed PCP orientation, through a displacement of AB polarity, gives rise to an extrinsic torque that bends the tissue. If this interpretation of the model is correct, if there's an extrinsic force imposing the orientation of the surface at every point, the model is not a serious contender to explain how an urchin is made. It's maybe less critical for other examples but leaves me unclear how seriously the possibilities illustrated by the model can be taken. That being said, the model certainly stands as a stimulant for further thinking along the same lines.

The above intuition about cell movement/rotation is correct as is now cooperated by Figure 2—figure supplement 3B. Also it is correct that our previous gastrulation Figure 7 used PCP to twist AB polarity and thereby enforce invagination from the outside. As our whole description of this process was too complicated we have changed to a simple one step scenario where all the gastrulation occur due to one particular set of initial conditions, and where invagination directly is caused by an initial twist of AB polarity of a subset of cells. We cannot claim that this a realistic scenario, but we propose it may be an indirect way of simulating the effects of bottle cells, and keep this example to illustrate the potential of our approach. Also we maintain the neurulation figure in supplement (Figure 1—figure supplement 11) although this still relies on the speculative assumption that changed PCP could direct movement of AB polarity. In any case the main role of PCP in our sea urchin gastrulation model is to direct the convergent extension after an initial deformation.

Similarly, in the case of fold formation on a sphere (Figure 5), which the authors suggest can occur if cell proliferation is on a faster time scale than relaxation, the physical origin of a slow relaxation time scale, and the mechanism of the observed instability, would merit discussion. Although assuming relaxational dynamics is standard and justified when the local minima of the energy are of interest, it might be questioned when turning to dynamical phenomena. Do folds develop through buckling, because the sphere grows (in number of cells) faster than it can expand in the surrounding medium? In many contexts, mechanical equilibration is much faster (seconds or minutes) than proliferation (hours). What would be the origin of a longer time scale for relaxation? Resistance from the surrounding medium might be one, but the authors present fast proliferation as an alternative mechanism to pressure from growth in a stiff gel.

Author response image 2
At cell division, the daughter cell equilibrates by half a cell radius in one time unit, which is of order 1/1000 the generation time.

Here, we show how a new cell (in blue) reaches equilibrium (in red). Cell division happens at time 0. At time 5, the next consecutive division happens in the system (not shown). The systems consists of 200 cells placed on a hollow sphere which is the initial condition in our organoid simulations (Figure 5). A new cell is introduced half a cell radius in a random direction away from the mother cell. The red line is the mean distance from the mother cell to all it’s neighbor cells at t = 1000. The generation time is intermediate (Figure 5A).

https://doi.org/10.7554/eLife.38407.030

We appreciate the comment, which in part is associated to our poor description of the parameter K that should be divided by 10000 to become a growth rate. We have clarified this in the new Figure 5, and redefined the proliferation rate K as generation time tG. Buckling by growth in Figure 5 occurs when the doubling times of cells are between 140 and 4000 counted in our time units. From the potential in Figure 2 the typical displacement in one time unit is expected to be up to 0.5 cell radius (Figure 2, see also Author response image 2). Translated into a cell diameter of 10 micrometer and assuming that a cell generation of 1 day correspond to 1200 time units in the model, the cell movements correspond to up to 2 micrometer/minute.

A relatively fast local relaxation does not mean that the cells reach equilibrium at finite division rates, because the pressure from one division takes longer time to equilibrate on the longer distances that involves many cells. The entire system is driven out-of-equilibrium when the time between two consecutive divisions anywhere in the system becomes shorter than the relaxation time at the corresponding distances. Interestingly, by quantifying distribution of cell shapes it has been found that the growth of MDCK organoids is out-of-equilibrium and they implicated this to be a possible cause of organoid anisotropy (Cerruti et al., 2013).

I also have some issues of form with the manuscript (e.g. I do not find the connection that the authors draw between the model and evolution very substantive), but the above questions on the physical relevance of the model are the main issues that I would like to see clarified before the manuscript is given further consideration.

We agree, and have removed the evolution comment from the manuscript.

Reviewer #2:

1) This manuscript describes a mathematical model and associated simulations of large clusters of cells that can be assigned properties of apical basal polarity and planar cell polarity and may proliferate, to follow their 3D morphological evolution. The aim is to test the hypothesis that "stability, complexity, and diversity are emergent properties in populations of proliferating polarized cells." The authors conclude that "polarity enables complex shapes robust to noise but sensitive to changes in initial and boundary constraints." The work is of potential interest, but a number of substantial issues should be addressed. Most important, some essential definitions are not provided, the writing needs to be made more accessible to biologists who are not modelers, some key assumptions should be made explicit, and others.

We made numerous changes throughout the manuscript, as well as some additional supplementary figures that explain the way our potential acts (Figure 2—figure supplement 2 and Figure 2—figure supplement 3). See below for more details.

Essential revisions:

2) Definitions needed: What precisely is meant by "stability?" Are the authors referring to feature constancy over time? Do they mean that models will converge to a specified set of features? What features? For example, they suggest that the trajectory in Figure 3 is stable, but noise seems to qualitatively change the shape (Figure 2—figure supplement 1), though it still produces a complex folded shape. It's not clear what, exactly, the authors are claiming.

Thank you for pointing this imprecision out. By stability we mean both constancy in time and ability to converge to similar sizes and relative position of lumens and folds in presence of noise. We realized that the reader may have been confused in comparing Figure 3D with Figure 2—figure supplement 1E–F. The simulations in Figure 2—figure supplement 1 and Figure 3 were started from different initial conditions and should not be compared for shape similarity.

The intention was to separately illustrate stability in time (in Figure 3) and stability to noise (in Figure 2—figure supplement 1E–F, where we compare end-results of the simulations run at high, Figure 2—figure supplement 1F, and low, Figure 2—figure supplement 1E, noise). The shapes emerging under high and low noise are not identical, however a high level of similarity – relative position and sizes of the majority of channels and lumens – is preserved. Ideally, this larger scale qualitative comparison, sketched in Figure 2—figure supplement 1E–F would benefit from a more rigorous quantification, by e.g. mapping changes in network of channels and lumens, but it is non-trivial and we judged it to be outside of the scope of this article. (To compare the relative impact of noise vs. initial conditions, Figure 2—figure supplement 1, we used distributions of pairwise distances between cells that were in the same position at start. This simple method however does not work for comparing macroscopical features).

We have rewritten the subsection”The final shapes are robust to noise but sensitive to initial and boundary conditions” describing stability results. We also modified Figure 2—figure supplement 1E–G to visually point to macroscopic similarities.

While addressing this point, we realized that the observed changes due to high noise come from noise acting only at early time, during “expansion” stage. The same level of noise applied after the system reached metastable state does not cause any major macroscopic changes, further supporting the link between polarity and morphological stability. We have added a sentence and Figure 2—figure supplement 1G to cover this point.

3) I struggled with the idea of "boundary conditions." I eventually concluded that, here, the term appears to mean restrictions on polarity orientations, and maybe other parameters. However, to a biologist, "boundary conditions" means something very different.

We have now introduced biological context for what we mean by boundary conditions and then defined the term “boundary condition” in subsection ”Folding by pressure or rapid proliferation result in different fold-morphologies”.

4) Difficulty with this passage was compounded by the challenging sentence "At the same time, in these symmetric cases the differences in initial conditions (no boundary conditions) are not sufficient to generate different structures." This implies that differences in initial conditions means no boundary conditions. I think the authors mean to say, "differences in initial conditions but without imposed boundary conditions."

Yes, this was the original meaning. We have adopted the suggested formulation.

5) Finally, in the Discussion section the authors say, "boundary conditions, e.g. fixed polarities for a fraction of the cells in some time." A precise definition should be provided when the idea of boundary conditions is first introduced.

We now define boundary conditions and we have also changed this sentence to: “We find that the exact morphological details are defined by initial conditions, e.g. initial positions and orientation of polarities, and boundary conditions, e.g. polarities restricted to certain direction for a fraction of the cells” in subsection ”Phenomenological description bridges cell polarity to macroscopic morphologies”.

6) "Minima" is defined in the Materials and methods section, but the authors also use "maxima." By analogy, one can guess what this means, but please define.

The “pronounced local maxima” was used to describe the folded structures. To simplify, we have chosen to reformulate this sentence to “As cells divide faster, our simulations predict a transition from a smooth spherical shell to an increasingly folded structure with multiple pronounced folds” in subsection ”Folding by pressure or rapid proliferation result in different fold-morphologies”.

7) Improve writing to make accessible to biologists: I'm a biologist with ample experience collaborating with modelers. Still, I found the manuscript difficult to follow and needed several readings to understand. If the goal is to make this work accessible to a broad audience (and it should be if published in eLife), the authors should carefully revise the text to replace jargon with precise simple language.

We have replaced a number of technical terms with a simpler expressions (e.g. replaced “metastable”, “topology”). Defined the terms “boundary conditions”, “stability” and tubes, and tried to simplify the language throughout.

8) Describe implications of model structure: Several features of the model are not as explicit as they should be. One is the conflation of PCP with convergent extension behavior. PCP defines a vectorial molecular asymmetry, and PCP can have a number of outcomes. One is CE, but others include the asymmetric elaboration of cellular structures, or the asymmetric induction of gene expression that results in different cell fates, for example. Here, the model structure equates PCP to forces that drive CE.

We agree. To clarify model implications, we have added a the following paragraph to subsection”Model implications”: “The present formulation of PCP effects has several implications. First, we restrict the effects of PCP to directed (anisotropic) cell–cell adhesion and do not consider its other possible roles, in e.g. asymmetric cell differentiation, thus primarily focusing on its role in CE. Second, in our current formulation AB polarity and PCP influence each others orientation on equal footing. However, as PCP is typically established in the apical plane, and is thus secondary to AB polarity, the current view is that AB polarity can influence PCP orientation and not the other way around.”

9) Furthermore, the descriptions of this implementation are misleading. In the Introduction, the authors, in describing the implementation of AB and PCP polarity, say "In other words, the adhesion strength between neighbor cells is modulated by the orientation of their polarities." Are the authors describing how AB and PCP (=CE) are implemented in the model, or how they function biologically?

To clarify, we reformulated the sentence to: “In other words, in our approach the adhesion strength between neighbor cells is modulated by the orientation of their polarities.”

10) There is no evidence that molecular PCP functions by regulation of adhesion (involvement of Celsr proteins notwithstanding), and differential adhesion is only one part of a complex and only partially understood set of mechanisms that drive CE.

We would like to thank the reviewer for the healthy dose of skepticism, as while initially this was also our view – originally in subsection”Cell–cell adhesion depends on the orientation of polarity” we stated “However, unlike with tight junctions, there is no biological confirmation for preferred directional adhesion with PCP” – this comment motivated us to do a careful literature survey.

We found that both proteins regulating adherens junctions, e.g. Smash (Beati et al., 2018) as well as adherence proteins, proteins forming adherens junctions, e.g. Bazooka, E-cadherins (Simões et al., 2010; Tamada et al., 2017; Levayer and Lecuit, 2013; Warrington et al., 2013; Aigouy and Le Bivic, 2016) are planar polarized. These, we believe, support our assumption of anisotropic, planar polarized adhesion. In case, if the arguments below are notwithstanding, we would very much appreciate reviewer's opinion on where they fail.

We have now modified the paragraph in subsection”Cell–cell adhesion depends on the orientation of polarity”, below Equation 8, which now reads “Cells adhere to each other by membrane proteins assembled in adherens junctions just below the apical surface. Both proteins regulating adherens junctions, e.g. Smash (Beati et al., 2018) as well as adherence proteins forming adherens junctions can be planar polarized, e.g. Bazooka, E-cadherins (Simões et al., 2010; Tamada et al., 2017; Levayer and Lecuit, 2013; Warrington et al., 2013; Aigouy and Le Bivic, 2016). These indirectly support our assumption of anisotropic, planar polarized adhesion.“

We do agree that likely there are other mechanisms contributing to convergent extension, and have added this to the DDiscussion section which now reads “The presented formulations of model captures only some of the known events contributing to CE and does not include PCP-driven cell-shape changes mediated by e.g. apical constriction (Nishimura, Honda, and Takeichi, 2012) or contribution of external forces (Lye and Sanson, 2011).”

11) Another implicit part of the model as I understand it, and one that is probably very important to the outcomes, is the notion that altering PCP vectors can influence the AB vectors of cells. I know of no biological precedent to justify this behavior. This input assumption is inconsistent with my understanding of the biological reality. AB polarity is an initial condition that must be met for the subsequent establishment of PCP, which is then by definition orthogonal to the AB axis. There are examples in which PCP determines edges of groups of cells that then increase tension, thereby inducing folds, but this is a different mechanism from the one implied here.

We again appreciate the input, and have shown that constraining PCP to the apical plane (more details on this in our reply below) does not affect CE and formation of tubes neither when we start form spherical lumen (Figure 6—figure supplement 1) or simulate sea urchin gastrulation.

12) I think it would be interesting for the authors to explore how their model would operate if AB polarities orient PCP, but PCP cannot directly influence AB. I suspect that the examples of tube formation would fail. At a minimum, they should discuss the need for this assumption explicitly, noting that it is not currently supported by biological data (or citing the data if there is some that I'm not aware of).

Thank you for stressing this point. Indeed, in our formulation (described by Equation 7) PCP and AB polarity can influence each on equal footing, and we used this symmetry to propose an alternative hypothesis for invagination during sea urchin gastrulation (orientation of PCP set at the boundaries drive the orientation of AB polarity).

While we did explore the consequences of unidirectional, AB polarity to PCP, effects, we also found data supporting that PCP may be primary to, and direct AB polarity.

Consequences of constraining PCP: We have explored the consequences of constraining the hierarchy or polarities in two ways: First, in case of sea urchin gastrulation, without modifying our model we impose boundary conditions directly onto AB vectors (Figure 7; the entire section on sea urchin gastrulation was updated to reflect this change). Biologically, this is an implicit way to capture the effect of the bottle cells, where the change in cell shapes effectively reorients AB vectors in neighboring cells (see Figure 2—figure supplement 2). Second, we have directly tested how disabling PCP pull on AB polarity (setting λ2 to 0 when updating Equation 10) affects tube formation. We observe no qualitative difference compared to our original results in neither tube formation (Figure 6—figure supplement 1), nor sea urchin gastrulation (not included), thus suggesting that we could relax our assumption symmetry between PCP and AB influence.

Data supporting that PCP may influence AB polarity: The idea that AB polarity and PCP may act on equal footing is appealing both for its symmetry and simplicity, and we were excited to find a number of recent data supporting that PCP is not necessarily secondary to AB polarity and may regulate AB polarity and influence its orientation.

i) Cells can acquire PCP without AB polarity present. The two examples are the radial intercalation of the multiciliated cell progenitors in Xenopus epithelia (Ossipova et al., 2015) and transition from the solid cord to a tube in intestine morphogenesis of Zebrafish (Baer, Chanut-Delalande and Affolter, 2009; Zorn and Wells, 2009). In both of these examples, cells acquire AB polarity after convergent extension / intercalation events driven by PCP.

ii) Proteins required for AB polarity can become planar-polarized. These include: Baz/Par-3 protein (Warrington, Strutt and Strutt, 2013; Aigouy and Le Bivic, 2016) which in Drosophila development has a role in both AB polarity and PCP (Beati et al., 2018); PCP regulates Lgl, protein central for establishing AB polarity (Choi and Sokol, 2009; Dollar et al., 2005) and Lgl as well as other AB polarity proteins are planar polarized in Drosophila epidermis (Kaplan and Tolwinski, 2010).

iii) The changes in cell-shapes during invagination (e.g. sliding of adherens junctions and formation of bottle cells) are regulated by PCP. These changes in cell-shapes effectively reorient AB polarity (see Figure 2—figure supplement 2). Apical constriction leading to formation of bottle cells is under PCP control during neural tube closure (Ossipova et al., 2014; Nishimura, Honda and Takeichi, 2012; Kinoshita et al., 2008); gastrulation in C. elegans (Lee et al., 2006), sea urchin (Croce et al., 2006) and Xenopus (Choi and Sokol, 2009). Also, the sliding of adherens junctions leads to shorter AB axis and is mediated by a Baz (key protein in PCP and AB polarity) during Drosophila gastrulation (Wang et al., 2012). Thus, the data shows that PCP can drive changes of the cell-shape (and consequently orientation of AB polarity) and while does not directly explain how the orientations of PCP and AB polarity are coordinated, it supports our assumption that components of PCP pathway can influence the orientation of AB polarity.

These points are now included in the subsection ”Model implications”.

13) Some unlikely assumptions: While it may be interesting that complex forms such as the sea urchin-like gastrulation and the folded neural tube CAN be produced within the confines of the model, the conditions needed stretch the bounds of plausibility. The model simulating sea urchin gastrulation requires a precise choreography of PCP parameters that could happen but seem unlikely. Even more difficult to believe is the initial conditions invoked to simulate neural plate folding. Specifically, invoking "two rows where the PCP points parallel and antiparallel to AB polarity" invokes conditions for which I know of no precedent at all.

We should have been more explicit on how our assumptions for both cases relate to the current understanding of these processes. In case of sea urchin gastrulation, the process is believed to consist of two stages, primary invagination driven by swelling of the inner layer of extracellular matrix beneath invaginating cells (Lane et al., 1993) and formation of a ring of bottle cells due to apical constriction (Kimberly and Hardin, 1998) and secondary invagination were tube elongates due to CE and pull by mesenchymal cell (Lyons, Kaltenbach, and McClay, 2012). This was the motivation for two discrete stages in our original model formulation.

We have now converged to a simpler, less choreographed version, where all conditions can be set prior to invagination and PCP remains in apical plane. First, we set boundary condition such that planar cell polarities of the invaginating cells are oriented on a spiral, and are all in the apical plane. This constraint allows for convergent extension and is present at all times. While not documented, it is consistent with observed effects of WNT orienting PCP within apical plane (Humphries and Mlodzik, 2017). Second, instead of driving invagination by PCP orientation, we simulate the effect of bottle cells (see Figure 2—figure supplement 2) by reorienting AB polarities. This external force on AB polarity can be thought of as a combined effect of bottle cells and bending by swelling of neighboring extracellular matrix. Motivated by the observations in other species where the changes of cell shapes are sufficient for invagination (does not require swelling of neighboring extracellular matrix), we have tested the case where the AB polarities of bottle cell neighbors prefer to be tilted towards each other (as in Figure 2—figure supplement 2). This completely cell-intrinsic scenario produces similar invagination outcome.

These changes are now described in a new Figure 7 and updated subsection ”Phenomenological description bridges cell polarity to macroscopic morphologies”.

The most provoking and least likely assumptions in our initial formulation was the hypothesis that boundary conditions may set PCP to point out of the apical plane. We formulated it as direct effect of instructive cues by organizing morphogens (e.g. WNT), however this effect may be indirect and act through e.g. reported cell-shape changes where PCP drives apical constriction (see our reply to point 12), which tilts AB axis in neighboring cells, thus tilting PCP orientation.

In case of neurulation, the idea was to simulate the neural plate (cells in the middle, between the two rows with constrained PCP) surrounded by the epidermis (the rest of cells). The two rows of cells with PCP pointing out of epithelial plane would then correspond to the cells at the dorsolateral hinge points next to the neural plate (epidermis boundaries). In chick spinal neural tube can close with only these two hinge points (Nikolopoulou et al., 2017). The bending is driven by apical constriction and PCP is essential for bending, CE and closure (Nikolopoulou et al., 2017). We used very simple brute-force constrain on PCP (pointing out of the epithelial plane) at the boundaries, meant to represent hinge points, to show that this could drive the entire process. Here again, this may be an indirect effect of PCP → Apical constriction → tilt of AB polarity → PCP points out of the original plane of epithelium. Including an intermediate step linking PCP with apical constriction (by e.g. introducing potential favoring tilted AB polarities) and convergent extension is rather simple and would likely result in a more robust model and require less fine-tuning of boundary conditions but is beyond the scope of the current manuscript.

We have updated the Discussion section to include our reflections on PCP pointing out of plane.

14) Additional discussion of testable predictions needed: Re-emergence of folds in pancreatic organoids, pressure causes fewer and much deeper folds that form early during growth. This seems to be a qualitative measure. How would biological observations be assigned to one or the other category?

Thank you! We realized that we did not stress the quantitative aspect well enough and have edited the Discussion section.

In addition to the qualitative description, our results also suggest that the two scenarios can be distinguished quantitatively by the differences in depth and length of the folds (new Figure 5—figure supplement 1). We now refer to this figure also in the Discussion section.

15) Isn't the key difference that the pressure hypothesis predicts that the number of minima maxes out, whereas the proliferation model predicts that the number keeps growing with increase in size (Figure 5)?

Yes, it is indeed a key difference and is now clarified both in the Results section and in the Discussion section.

16) Re tubes formed by nonproliferating cells, "in PCP mutants with shorter and wider tubes one would expect higher fraction of "AB adhesion" relative to "PCP adhesion". However, as noted above, PCP is not directly related to any known adhesion events.

We have now elaborated on this prediction in the Discussion section. Please also see our reply to point 9 above. We also realized that our formulation “PCP adhesion “was probably confusing. We did not mean to claim this has to be a direct effect of PCP. Because we represent cells as point particles, the PCP in our model is a phenomenological description and is meant to represent the combined effect of anisotropic, both extra- and intracellular, factors (e.g. anisotropic apical constriction, polarized junctional actomyosin) that would direct cell movements and favor its position among other cells depending on orientation of these factors.

17) Re the conditions differentiating between tubes forming perpendicular (e.g. sea urchin gastrulation) or parallel (as in Drosophila gastrulation or neurulation) to the plane of epithelium, the measurements could in principle be made, but as noted above, the conditions seem to me to be highly implausible, and ask us to ignore the likely contribution of other signals.

We have now added an extensive discussion addressing both the implausibility of PCP pointing out of plane of epithelium as well as the contribution of other factors (Discussion section).

Reviewer #3:

1) This is a creative and innovative article on development. The authors introduce a simple theoretical model based on the idea that apical-basal (AB) polarity and planar cell polarity (PCP) are (the) main driving forces relevant for the formation of complex tissue morphologies. The model idealises cells as point particles with two independent axes (polarities), an AB and a PCP axis. Pairwise interactions between cells contain an attractive and a repulsive component where the magnitude of the attractive component depends on three additive factors: alignment of the AB axis, alignment between AB and PCP axis, and alignment of PCP axis. Noise comes into the model at various levels: Langevin noise in the equations of motion, noise at the level of choosing the axes orientation at cell proliferation, and noise at the level of initial conditions. In summary, this seems to be a computationally efficient model to explore the question: How can cell-cell adhesion that depends on the relative orientation of polarisation axes alone drive different morphologies during development?

Indeed, the authors simulations show that their model is able to reproduce a broad range of morphologies depending on several factors: choice of initial conditions, relative tendencies of the various axis to orient with respect to each other, and cell proliferation rate. A clear strength of the computational model (I would not call it a theory) is that it provides a mechanistic link between morphology and cell-cell interaction as well as initial conditions. Thus, I find it a clear conceptual advance in the field which puts a fresh perspective on the role of generic interactions for the emergence of different tissue morphologies.

We agree “theory” may overstate the nature of our work. The main motivation for calling it theory was to highlight that our approach is intentionally phenomenological and can work with few parameters. Computational model is a broad term and often associated with dozens of parameters. We have compromised the two by replacing “Theory” in the title with “Theoretical tool”.

2) The authors also list a few testable predictions. Here, I am not so convinced yet for very basic reasons: Testable predictions require either a clear link between a molecular and an emergent process, or a set of mutually exclusive alternatives that can be tested.

We agree that by its nature our phenomenological description is outside the realm of mechanistic predictions targeted at specific proteins. While the alternatives we consider are not mutually exclusive, one can in principle test how much each contributes to the final result e.g. varying the stiffness of the gel and generation time. We now added a more detailed discussion on this in the Discussion section.

The scenarios explored for pancreatic organoid folding are not limited to only two possibilities mentioned in the paper. (Other possibilities are the instabilities associated to secreted morphogens and the bias of cell differentiating into slow and fast growing cells depending on their position in the fold.) We see our model and predictions as a starting point and motivation for quantitatively testing organoid morphologies as our results suggest that different scenarios may result in quantitatively different shapes. The model will be rather easy to modify to either consider other alternatives or explore the combined effect of several, but we lack constraints by experimental data for this to be meaningful.

3) Concerning the latter case (in the section of organoids) the authors analyse two possible scenarios (cell proliferation and a specific way to implement external forces). These are not two mutually exclusive alternatives but only two out of many possible conceivable scenarios. Concerning the former case, in the example on sea urchin gastrulation the authors seems to a large degree put in what they would like to get out.

In regards to gastrulation, we have simplified the modeling of gastrulation greatly in the new version, reducing it to one step: one initial condition of twisted AB polarities and one boundary condition in form of an external constraint/forcing on PCP. This emphasizes the simple interplay of initial invagination and convergent extension, and makes the connection between input and output even more transparent. We very much agree that a more realistic approach should include forces that favour invagination from cell shapes (e.g. bottle cells) and various checkpoints and differentiation processes along the developmental pathway.

4) Finally, and this is more on the side of a personal opinion, I am not convinced that morphologies that strongly depend on initial conditions are a robust way for cellular development. For instance, how can an organism make sure that all these initial conditions are properly chosen?

Agree, but in Figure 2 we did not have any external constraints (boundary conditions). Overall, we try to convene that both initial conditions and external constraints are important. Initial conditions can be overruled by boundary conditions, i.e. persistent constraints on cell polarities by e.g. neighboring tissues through intra- and extracellular feedbacks. The robustness of the developmental process most likely stems from these constraints. All the finally obtained morphologies are subsequently robust to noise, growth and disruptions (see also the new Figure 3—figure supplement 2).

5) In summary, I found the paper a well-written, creative, and thought-provoking contribution to the field of development. It puts the emphasis on local and rather generic interactions between cells' polarities. If the authors can convincingly reply to my above queries, I would be inclined to recommend it for publication.

We are very happy for this overall positive view on our work, and hope that our response and the resulting changes, including in particular our simpler gastrulation treatment and more modest claims, sufficiently improved the manuscript.

https://doi.org/10.7554/eLife.38407.033

Article and author information

Author details

  1. Silas Boye Nissen

    1. Center for Models of Life, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
    2. StemPhys, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
    Contribution
    Conceptualization, Resources, Data curation, Software, Formal analysis, Supervision, Funding acquisition, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-9473-4755
  2. Steven Rønhild

    Center for Models of Life, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
    Contribution
    Software, Writing—review and editing
    Competing interests
    No competing interests declared
  3. Ala Trusina

    1. Center for Models of Life, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
    2. StemPhys, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
    Contribution
    Conceptualization, Resources, Data curation, Software, Formal analysis, Supervision, Funding acquisition, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    trusina@nbi.ku.dk
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-1945-454X
  4. Kim Sneppen

    Center for Models of Life, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
    Contribution
    Conceptualization, Resources, Data curation, Software, Formal analysis, Supervision, Funding acquisition, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    sneppen@nbi.dk
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-9820-3567

Funding

Danmarks Grundforskningsfond (DNRF116)

  • Silas Boye Nissen
  • Ala Trusina

Seventh Framework Programme (FP/2007/2013/ERC no. 740704)

  • Kim Sneppen

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

We thank Anne Grapin-Botton for insightful discussions on organoids and Sigurd Carlsen for discussion on tube formation. This research has received funding from the Danish National Research Foundation (grant number: DNRF116) and the European Research Council under the European Union's Seventh Framework Programme (FP/2007 2013)/ERC Grant Agreement number 740704.

Senior Editor

  1. Naama Barkai, Weizmann Institute of Science, Israel

Reviewing Editor

  1. Aleksandra M Walczak, École Normale Supérieure, France

Publication history

  1. Received: May 16, 2018
  2. Accepted: November 13, 2018
  3. Accepted Manuscript published: November 27, 2018 (version 1)
  4. Version of Record published: December 6, 2018 (version 2)

Copyright

© 2018, Nissen et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 1,080
    Page views
  • 197
    Downloads
  • 0
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading

    1. Computational and Systems Biology
    2. Developmental Biology
    Inna Averbukh et al.
    Research Article
    1. Developmental Biology
    Debadrita Bhattacharya et al.
    Research Article