Wall mechanics and exocytosis define the shape of growth domains in fission yeast

The amazing structural variety of cells is matched only by their functional diversity, and reflects the complex interplay between biochemical and mechanical regulation. How both regulatory layers generate specifically shaped cellular domains is not fully understood. Here, we report how cell growth domains are shaped in fission yeast. Based on quantitative analysis of cell wall expansion and elasticity, we develop a model for how mechanics and cell wall assembly interact and use it to look for factors underpinning growth domain morphogenesis. Surprisingly, we find that neither the global cell shape regulators Cdc42-Scd1-Scd2 nor the major cell wall synthesis regulators Bgs1-Bgs4-Rgf1 are reliable predictors of growth domain geometry. Instead, their geometry can be defined by cell wall mechanics and the cortical localization pattern of the exocytic factors Sec6-Syb1-Exo70. Forceful re-directioning of exocytic vesicle fusion to broader cortical areas induces proportional shape changes to growth domains, demonstrating that both features are causally linked.


SUPPLEMENTARY NOTE
In this note, we present two models for the mechanics of wall expansion and cell morphogenesis. Wall expansion involves the deposition of new wall material and the mechanical deformation of the pre-existing wall fabric. While elastic deformations of 5% are common in many walled cells [8], most mechanical models of walled cell morphogenesis do not include their contribution explicitely. Such omission is understandable since small elastic strains, although they may be important in controlling growth and have been considered as such, can hardly alter the cell geometry directly. The clearest evidence for this conclusion comes from observing the shape of plasmolysed walled cells, which rarely differs significantly from the shape of their turgid counterparts. In contrast, plasmolysis of actively growing S. pombe cells reveals large elastic deformations often exceeding 25%. Moreover, the relaxation of these elastic strains by plasmolysis can have a drastic effect on cell shape ( Figure 2e in the main text). Therefore, we were compelled to include wall elasticity in our models of growth domain morphogenesis. Our simulation approach follows closely the approach adopted by Dumais et al. [9] and Rojas et al. [10] with one significant modification, the elastic loading of wall elements is also included explicitly. To implement the elastic loading of wall elements, we follow the approach of Goriely and co-workers [11] and decompose cell morphogenesis in two surfaces evolving in parallel ( Figure A). These are the observed turgid surface of the cell and the relaxed (plasmolysed) surface. Every material point S of the relaxed surface is associated with a point s(S) on the turgid surface. Under normal growth conditions, the cell maintains its turgor pressure throughout the cell cycle. Consequently, the relaxed geometry is never "experienced" directly by the cell; it is merely a ghost geometry that evolves in parallel with the turgid geometry. All physiological functions are performed within the confine of the deformed cell surface. Yet, a complete account of cell morphogenesis must include this virtual relaxed geometry since it is this surface that records faithfully the history of wall assembly.

Simulation of Cell Morphogenesis: Membrane Model
We model cell wall expansion as the product of two co-occurring processes. First, the wall is growing due to the action of enzymes which can break load-bearing bonds in the glucan network and insert new glucan chains. This growth process affects directly the relaxed geometry of the cell since new material becomes integrated in the wall. Second, the wall can respond elastically to any change in the forces it experiences. This elastic effect is encapsulated by the mapping between the relaxed and turgid cell geometries. This decomposition can be applied directly to the meridional velocity of material points as visualized by the displacement of Qdots. We have: v(s(S)) = ds(S)/dt = (ds(S)/dS)(dS/dt) = λ s (S)V (S), where V (S) captures the flow associated with the incorporation of new wall material in the relaxed geometry and λ s is the elastic stretching of wall element when the relaxed geometry is deformed by turgor pressure.
Substituting v(s(S)) = λ s (S)V (S) in the equations for the strain rates in an axisym- metric shell, we get the following decomposition of the strain rates: where the stretch ratios λ s = ds(S)/dS and λ θ = r(s(S))/R(S) were used.
Simulation algorithm -We have implemented the model described above with a sequences of seven steps applied recursively.
i) We first input the observed curvature (i.e. the growth domain geometry) and fluorescence profile from a particular experiment. This step sets the initial conditions for the simulations and provides a benchmark to evaluate the ability of a given cortical marker to reproduce the morphogenesis of the OE and NE.
ii) We compute the elastic strains and stresses for the geometry. The membrane stresses are given by the force balance between the cell's turgor pressure (P ) and the tensions in the wall surrounding the cell. Given that the wall is thin compared to the typical radius of the cell, the meridional and circumferential stresses can be expressed directly in terms of the turgor pressure and the local geometry of the cell surface [9]: where δ = 0.2µm is the wall thickness, and κ s and κ θ are the principal curvatures of the surface.
To compute the elastic strains, we used the stress relations above and the material properties (ν = 0.3 and E/P = 40) which, within a broad domain of morphogenetically compatible material properties, were closest to the material properties measured in the plasmolysis experiments ( Figure 2d in the main text). The equations for the elastic strains are: For simplicity, these equations assume the cell wall to be isotropic, homogeneous and linearly elastic. It is likely that these assumptions are not perfectly satisfied although any more complex model would add free parameters which so far have not been measured.
iii) We compute the strain rates associated with growth. The growth strain rates are: where γ(s) is the fluorescence intensity profile of the marker of interest and α is a factor relating the fluorescence intensity with the activity of the marker in question. We note that this factor has no effect of the shape of the growing cell end but sets the rate of elongation.
iv) We compute the meridional velocity of material points. The meridional velocity is simply the spatial integral of the meridional strain rate˙ g s computed above.
v) We compute the rate of elastic loading due to the growth process. The rate of elastic loading is the second contribution to the deformation of wall elements. It is given by the relations:˙ e s = 1 vi) We compute the Lagrangian velocity field from the total strain rates. The total strain rates are the sum the contribution of strain rates arising from growth of the undeformed cell geometry (˙ g s and˙ g θ ) and elastic loading of wall elements (˙ e s and˙ e θ ). For an axisymmetric cell, the two total strain rates can be integrated to give the normal and tangential speed of displacement of material points [9]. The integration gives: We take a small step forward to find the new geometry. Given the velocity of every point on the meridian, the new, deformed, geometry is determined by finding the position of each point after a small time interval. This process maps the material point trajectories as if the were Qdots. During a growth interval, some material points are displaced from the growing dome to the non-growing cylinder. These points are eliminated and replaced by new points within the growth region such that the length of the growth region and its spatial resolution are preserved. This is done by fitting the deformed meridian with a cubic spline and remeshing with a uniform spacing between points.
Steps ii to vii are repeated until a steady state is reached.
Validation -Since we are using the growth simulations to test the compatibility between the distribution of molecular markers and the geometry of the fission yeast cell, it is imperative that the simulation protocol be exact at the kinematic level; that is, all the geometrical relations embedded in the implementation of the model must be exact to any accuracy desired. We first confirmed the convergence of our algorithm when run with the analytical solution for the isotropic growth of a spherical cap ( Figure B). According to this solution, the strain rates must vary as the cosinus of the angle of the surface normal. Second, we verified that the kinematics is "closed"; that is, it is possible to reproduce the geometry of a given cell exactly when the strain rates extracted for this cell are used as input to our simulation protocol ( Figure B).

Simulation of Cell Morphogenesis: Bending Model
The membrane assumption, which states that turgor pressure is supported only by inplane tensions in the wall, is widely employed to model the growth of walled cells such as fungal hyphae and pollen tubes. Its main attributes are its robustness and simplicity. The aim of the bending model is describing not only the growth of the cell end but the whole 'morphogenetic' cell cycle including the deformation of the septum into the NE following cell division. Because of the large change in curvature at the septum, the description of division must include both bending and transverse shear terms in the shell equation. Our bending model of yeast cell growth is based on the model of Su and Taber [12]. The equations were solved for an axisymmetric shell approximating the yeast cell geometry. However, as yeast cells grow without twist, we discarded all the terms in Su and Taber's equations which correspond to twist. The boundary conditions are that the two curvatures (κ 11 , κ 22 ) and two deformations (λ 11 , λ 22 ) are equal on the two points laying on the axis of rotational symmetry, while the transverse shear Q 2 is 0 at those points. We solved the governing equations according to the method outlined by Kempski et al. [13] and making use of the function fsolve in Matlab (see below). While testing our model on an inflating sphere, we observed that the set of equations proposed by Su and Taber does not maintain the spherical geometry even at quite low deformation. We found that this asymmetry originates from an oversimplification proposed by Reissner [14] using r 0 rather than r in the derivative of the ODE. After removing the simplification, the inflating sphere remained spherical. The system we ultimately used corresponds to the one proposed originally by Reissner (before the simplification) including the transverse shear term of Su and Taber.
The system of equations behind the bending model is: where the variables are as follow: r, z are the radial and axial coordinates of the axisymmetric contour. ii=11, 22 indicates the circumferential and meridional directions, respectively. λ ii is the deformation along direction ii. γ 2 corresponds to the shear deformation. κ ii is the difference between rest state curvature and deformed state curvature along direction ii. Other parameters are used for numerical purposes (H 2b , V 2b , M 22b , α). They are defined by the following relations: The following linear elastic energy including both transverse shear and bending was sufficient to describe plasmolysis experiments. 4(1+ν) , and C = Eh 3 24(1−ν 2 ) . The system of ODE was closed by calculating the following variables λ 22 , κ 22 , γ 2 : Simulation algorithm -For matter of clarity, the preceding system of equations is rewritten in the following form: f stands for the differential equations (eqns. 11-16). g stands for the algebraic equations (eqns. 22-24). l stands for the total arclength of the plasmolysed contour. Two functionals of u are introduced: Solving the boundary condition problem is equivalent to finding a function u which is a zero of F and g on the whole interval ]0, l[ and a zero of BC at both poles.
To solve the equations, we first defined the initial solution U 0 using the plasmolysed, stress-free contour (p 0 = 0). Then the pressure is increased by a small step δ : p n = p n−1 + δ. The problem is discretised in the following way. Two new functionals G and H are defined: The problem H(U, p n ) = 0 is solved using fsolve of Matlab among the cubic splines whose nodes are linearly spaced between 0 and the arclength l of the initial contour. As f diverges at both poles, the function G inside the function H is evaluated at the s i which are evenly distributed between ( , l − )( = 10 −5 ). The s i are distributed on the whole interval but are denser close to the pole. The guess function used to initiate the function fsolve is the cubic spline U n−1 . This numerical code was validated by calculating the solutions for geometries where the boundary condition problem admits analytical solutions for the whole contour (the sphere) or for part of the contour (a very elongated cylinder with two hemispherical cell ends). For these geometries both analytical and numerical solutions match very well.
In order to model the transition between NE and OE, a growth model including bending terms was implemented. In this second model the strains induced by the turgor are directly calculated from the plasmolysed state whereas they were evaluated from the deformed configuration in the first model. The set of equations used to describe growth conserves the ODE system (eqns. 11-16) but each of the strains and the bending deformations are decomposed in two subvariables: The "residual" part is determined at the precedent growth step. The "step" part is calculated at the current step. The system of ODE (eqns. 11-16) is now closed by the seven following algebraic equations: The following steps were followed.
i) The initial conditions used for the simulations were the geometry of the plasmolysed cell with the material properties (Poisson's ratio ν = 0.033 and a wall thickness h = 0.2µm). The initial geometry is assumed to be relaxed (free of bending and membrane stresses): ii) The initial geometry subsequently swelled up to a ratio E/P of 58. The "step" part of the strains and bending strains equal the strains and bending strains induced by this swelling.
iii) The new curvilinear coordinate s n is calculated: iv) The residual strain and bending strain are defined: These prefactors correspond physically to an increase of the rest length proportional to the deformation (λ ii,n−1 − 1) times the fluorescence profile K of the cortical marker we want to study. The proportionality factor δ quantifies the amplitude of each step of the growth.
v) The grown shape of cell is obtained by solving eqn. 32 with fsolve where the function g inside the formula (eqn. 31) is replaced by the condition (eqns. 35-41). The guess function used to initiate fsolve is: U guess n = (r n−1 , z n−1 , ϕ n−1 , V 2b,n−1 , H 22b,n−1 , λ step 11 = 1, κ step 11 = 1, κ step The last three steps are reiterated until an asymptotic cell end shape is reached. Provided that δ is small enough (< 0.3), its values does not affect the asymptotic shape. We in fact found that the simulations converge to the same asymptotic cell end geometry irrespective of the initial conditions ( Figure C).
The way of including the bending in the growth of shell is still a matter of debate. For this reason we adopted the simplest approach with no term of bending in the growth but a term of bending in the energy. As a control of robustness of the hypothesis, we carried out two supplementaries sets of simulations one with no term of bending in the energy and another one with the most naive way of including the growth in the bending strain (eqn. 48) by multiplying the rest radius of curvature by the same factor as the rest length in (eqn. 47). The first method did not affect the attractor at all for both the asymptotic curvature and the asymptotic diameter but it changes the dynamic of the curvature because of the initial bending momentum. The second method gives a slightly The simulations converge to the same steady-state geometry irrespective of the initial cell geometry. The contours on the left represent the relaxed (red) and loaded (blue) geometries used as initial conditions for the simulations.
smaller asymptotic diameter as the main method (the decrease lays between 5% and 10 % depending on the cortical marker) and increases sharply the curvature at the apex giving an aspect pointy which is not physiological. The conclusion of these two controls is that the asymptotic diameter is quite robust whereas the cell end curvature is a litlle less robust and depends on the way the bending term in included in the growth equations.
Simulating cell division -To simulate cell division, an axisymmetrical shell whose profile is an idealized S. pombe shape is swelled using the material properties estimated from the simulation of the plasmolysis experiments (see above). An initially straight septum is added at the middle of the swelled cell. The rest length of the septum can be tuned by setting λ 11 = λ 22 = λ * for the septum in equation (42). The size and the curvature Γ trans of the transition between the horizontal septum and the vertical cell wall can also be tuned. The value of the rest length in this zone of transition varies linearly between the value of the swelled cell wall at the insertion of the septum and the value in the septum. The shape of the septum once equilibrated with turgor is obtained by solving for the shell equations listed above.
Two independent parameters have to be estimated for this model: the curvature of the transition Γ trans and the rest length of the septum λ * . Γ trans was estimated by fit-ting circles to the transition zone between the septum and the sides of the cell ( Figure  D, panel b). The average curvature is: 2.27 ± 0.11µm −1 (n = 39). λ * was estimated indirectly. Simulations were run for the measured value of Γ trans and for a wide range of λ * . For each simulation, different measurements of the deformed septum geometry were recorded: D s the diameter after division, D 0 the diameter before division, H the height of the bulged septum, Γ ap the curvature at the apex. The same parameters were measured on cell images after and before division. First the contour was taken on the brightfield images (the selected points were located at the middle of the black thick line at the contour of the cell). As it was difficult to determine where were the septum extremities after division, the diameter after and before division was calculated by a geometrical method (see Figure D, panel a). The contour was first symmetrized by looking for the best fit by a symmetrical contour. Then the symmetrized contour was oriented vertically (the coordinate x is along the radius, y is along the long axis of the cell). Contour before and after division were superimposed. The daughter cell was divided in two half parts: an old part far from the division plane and a division zone. Daughter cells, whose old part did not superimpose well with the mother cell, were first excluded. The diameter after division was automatically calculated as the maximum of the x coordinates of the divided contour in the division zone. The coordinate of this maximum after division are x max and y max .The diameter before division was chosen as the x coordinate of the contour before division situated at y max .The average ratio between the two diameters was 1.029 (n = 27). Then to include more cells in the statistics, a dilatation factor which best superimposed the two non-dividing zones was calculated. This dilatation factor was then used to dilatate the whole contour. The measurements were then repeated. The average ratio only changed slightly (1.034) (n = 39). H and Γ ap were estimated on the same contours. The same rest length λ * was sufficient to predict the different geometrical parameters of the bulging septum ( Figure D).
Comparison of the results for the membrane and bending models Finally, we compared the membrane and bending models to ascertain that our results are robust. Both the asymptotic radius and asymptotic curvature of the cell end are in good agreement between the two models (Table and Figure E).