Mechanical role of the spatial patterns of contractile cells in invagination of growing epithelial tissue

Authors: Yasuhiro Inoue, Tadashi Watanabe, Satoru Okuda, Taiji Adachi

Abstract

Epithelial invagination is one of the fundamental deformation modes during morphogenesis, and is essential for deriving the three‐dimensional shapes of organs from a flat epithelial sheet. Invagination occurs in an orderly manner according to the spatial pattern of the contractile cells; however, it remains elusive how tissue deformation can be caused by cellular activity in the patterned region. In this study, we investigated the mechanical role of the spatial patterns of the contractile cells in invagination of growing tissue using multicellular dynamics simulations. We found that cell proliferation and apical constriction were responsible for expanding the degree of tissue deformation and determining the location of the deformation, respectively. The direction of invagination depended on the spatial pattern of the contractile cells. Further, comparing the simulation results of surface and line contractions as possible modes of apical constriction, we found that the direction of invagination differed between these two modes even if the spatial pattern was the same. These results indicate that the buckling of the epithelial cell sheet caused by cell proliferation causes the invagination, with the direction and location determined by the configuration of the wedge‐shaped cells given by the spatial pattern of the contractile cells.

Introduction

In morphogenesis, deformation of epithelial tissue is essential for the formation of three‐dimensional (3D) epithelial organs. Invagination is a key deformation mode required to fold a planar epithelial sheet into a 3D shape. Epithelial tissue invagination is passively induced by ambient mechanical influences and actively induced by cellular activities such as cell proliferation and apical constriction. Cell proliferation changes the cell volume during cell growth and increases the number of cells by cell division (Mao et al. 2011, 2013; Amar & Jia 2013). Elastic energy can be stored in the tissue during cell proliferation, inducing subsequent tissue deformations. Apical constriction affects the deformation of epithelial tissue by causing the cells to actively change from a columnar shape into a wedge‐like shape. In addition, it is known that there are variations in the cellular mechanisms of apical constriction (Martin & Goldstein 2014). Surface contraction occurs when cells accumulate F‐actin on their apical side to form a thick F‐actin band that squeezes their apical surface (Lee & Harland 2007; Sawyer et al. 2010; Martin & Goldstein 2014). Line contraction occurs when apical constriction occurs locally on an edge between adjacent cells. Examples of line contraction include supra‐cellular myosin cables (Nishimura et al. 2007; Kondo & Hayashi 2015) and planar‐cell‐polarity based machinery (Mao et al. 2011; Nishimura et al. 2012). It remains to be elucidated how a particular apical constriction mode affects tissue deformation.

In the process of epithelial tissue invagination, cellular activities are known to occur in an orderly manner according to known spatial patterns, such as band‐shaped (Sweeton et al. 1991; Smith & Schoenwolf 1997), circular (Simões et al. 2006), or circumferential patterns (Nishimura et al. 2007). These patterns have been studied by regulating expression of genes controlling cellular activities (Solnica‐Krezel 2005; Fuccillo et al. 2006; Suzuki et al. 2012; Jidigam et al. 2015). In addition to uncovering the molecular mechanism for these patterned cellular activities, there should be a mechanical basis for understanding morphogenesis (Heisenberg & Bellaïche 2013). For instance, mechanics can be used to determine how the spatial pattern controls the force distribution in the tissue, and what kind of tissue deformation results.

In this study, we investigate the mechanical role of the spatial pattern of contractile cells in the epithelial invagination of growing tissue using multicellular dynamics simulations. Further, we examine how a particular apical constriction mode affects epithelial invagination. A 3D vertex model was developed to perform multicellular dynamics simulations (Honda et al. 2004; Okuda et al. 2013a), and tissue deformations were successfully reproduced based on a force balance at the cellular level (Honda et al. 2008a,b; Okuda et al. 2013b; Inoue et al. 2016). This model can be used to examine mechanical effects of the spatial pattern of cell activity under various conditions on tissue deformation.

Mathematical model

3D vertex model expressing multicellular dynamics

The 3D vertex model represented the shape of a cell as a polyhedron consisting of vertices and edges (Honda et al. 2004; Okuda et al. 2013a). The tissue was represented as an aggregate of multiple cells (Fig. 1A,B), where the vertices and edges of each cell are shared with neighboring cells. The vertices and edges composed a network that represents the entire shape of the tissue (Fig. 1C,D).

The equation of motion of vertex i with position vector xi at time t
(1) ηdridt−Vi=−∇iU.

Here, η is the friction coefficient and Vi is the mean velocity vector of vertex i. In the 3D vertex model, vertex i was directly connected to four adjacent vertices by edges. Indexing the directly connected vertices as j(i), the mean velocity vector can be defined
(2) Vi=15dridt+∑j(i)drj(i)dt.

The right hand side of Equation (1) represents a force exerted on vertex i derived from the total energy function U, which represents the mechanical properties and morphogenetic events of the cells
(3) U=∑jccellujcv+ujcs+ujch+ujcac, where ∑jccell indicates a summation over all of the cells. The energy function U includes the cell volume elastic energy ujcv, cell surface elastic energy ujcs, cell height elastic energy ujch, and apical constriction energy ujcac. The mathematical expressions for these energy functions are defined
(4) ujcvvjc=12kvvjcvjceq−12,

(5) ujcssjc=12kssjcsjceq−12,

(6) ujchhjc=12khhjchjceq−12,

(7) ujcacpjc=12kacpjc2(surfacecontraction),

(8) ujcacljce=12kac(e)ljce2(linecontraction).

These mathematical expressions were derived from a 3D vertex model of epithelial cells (Okuda et al. 2013b; Inoue et al. 2016). The jc‐th cell's volume vjc, surface area sjc, height hjc, apical circumferential length pjc, and edge length ljce between the jc‐th cell and its adjacent cells are represented as variables. The superscript eq in several variables in Equations (4), (5), (6) indicates the value in the stress‐free state. The constants kv, ks, kh, and kac (kac(e)) are the volume elasticity, surface elasticity, height elasticity, and apical constriction constant of the cell, respectively.

Equation (7) indicates that the entire perimeter of the apical surface contracts, while Equation (8) indicates that only the boundary line with a specific cell contracts. These represent the two apical constriction modes. Equation (7) captures the surface contraction mode, and Equation (8) captures the line contraction mode (Fig. 1E,F). To examine how a particular apical constriction mode affects epithelial invagination, Equation (7) was employed in Mechanical role in apical constriction and cell proliferation and Role of the spatial pattern of contractile cells in invagination , and Equation (8) was used in Mode of apical constriction in invagination .

Equations (7) and (8) are quadratic functions of variables pjc and ljce, respectively, and do not include the linear terms that would express a line tension (Fletcher et al. 2014). However, since the linear terms only change the net natural length that gives the minimum energy value with respect to apical constriction, for simplicity we do not consider them here.

All model constants are listed in Table 1. Here, to focus on invagination caused by the spatial pattern of the contractile cells, we ignored the asymmetry of physical properties that depend on epithelial polarity except for the contractile cells, and assumed a flat homogeneous epithelial monolayer sheet. To solve Equation (1), the parameter values were normalized by unit length (veq)1/3, unit time η(veq)2/3/kBT, and unit energy kBT. Here, kBT is formally introduced to normalize the energy parameters. To express cell proliferation, we employed a polyhedron‐division model and its parameters for the 3D vertex model (Okuda et al. 2013c), where we adopted a cell cycle time with mean τcc=100 and standard deviation 0.1τ cc.

Simulation conditions

The initial shape of the tissue was a flat epithelial monolayer sheet as shown in Figure 2. There were 400 hexagonal columnar cells in the calculation area. The polarity of the epithelial tissue in the initial state was assumed to coincide with the z‐axis. Periodic boundaries were employed for boundary conditions.

Tissue deformation

The magnitude of the tissue deformation at time t, md(t), was defined as the absolute value of the difference between the maximum and minimum values of the z‐component of the centroid position vector, gc, of each cell in the tissue. (9) md(t)=maxZc(t)−minZc(t).

Here, Zc(t) is the set of z‐components of the centroid position of each cell at time t.

To quantify the central position of the deformation, we defined the deformation region as the set Pd(t)
(10) Pd(t)=gc(t)|gc(t)·ez−zavegc(t)>zsdgc(t).

Then, the central position of the deformation, pd¯(t), at time t was defined as the arithmetic mean of the element gc(t) of the set Pd(t). Here, the set Pd(t) is a set of ones in which the deviation of the z‐component at the centroid position of the cell exceeds the standard deviation, zsdgc(t). The vector ez is a unit vector in the z‐direction, and zavegc(t) is obtained by averaging the z‐component of the centroid position in all the cells.

The depth of the invagination was quantified by Bessel functions that represent solutions of the Helmholtz equation in cylindrical coordinates (i.e., the radial component of the out‐of‐plane deformation of a circular elastic sheet). We measured the displacement in the z‐direction with respect to the radial axis, r, and expressed displacements u(t,r) using a superposition of Bessel functions by the least squares method (Arfken et al. 2012). Then, we obtained the displacement, uo(t):=u(t,0), at the origin as the depth of the invagination at time t.

Results

Mechanical role in apical constriction and cell proliferation

We examined tissue deformation as a function of (i) apical constriction, (ii) cell proliferation, and (iii) both apical constriction and cell proliferation. A circular pattern was used for the spatial pattern of the contractile cells under conditions (i) and (iii). Here, a regular hexagonal lattice represented the circular pattern, with each side initially consisting of nine cells. Cell proliferation occurred in whole for conditions (ii) and (iii). Simulation results are shown in Figure 3. Conditions (i), (ii), and (iii) resulted in a shallow basally‐convex invagination of the pattern (Fig. 3A), a large undulation of the tissue (Fig. 3B), and a large basally‐convex invagination of the pattern (Fig. 3C), respectively. Figure 3B shows that the tissue is shaped like a wave. This is due to buckling as a result of cell proliferation under periodic boundary conditions, because tissue deformation was limited to conform to the periodic boundaries. As will be shown below, such deformation occurs in random places.

The magnitude and location of the deformation along the z‐direction was quantified (Fig. 4A) by the method described in Tissue deformation . The magnitude of the deformation, md(t), for conditions (i)–(iii) is shown in Figure 4B. Conditions (ii) and (iii), in which cell proliferation occurred, resulted in a larger deformation compared with condition (i), in which cell proliferation did not occur. Comparing conditions (ii) and (iii), apical constriction resulted in an earlier deformation. Because apical constriction slightly bent the sheet to be basally convex as shown in Figure 3A, tissue deformation originated in this region and expanded using elastic energy stored in the tissue from cell proliferation. Consequently, the buckling of the epithelial cell sheet due to cell proliferation switched to invagination. Therefore, the mechanical roles of cell proliferation and apical constriction were the amplification of the tissue deformation and determination of the location of the tissue deformation, respectively.

We further examined the central position of the deformation, pd¯, for 150 samples with different random number seeds. Figure 4(C,D) shows the central position of the deformation with respect to time for each sample, which typically coincided with the centroid of the spatial pattern with apical constriction, while it was random in the absence of apical constriction, confirming that the spatial pattern of the contractile cells determined the location of the invagination.

Role of the spatial pattern of contractile cells in invagination

We next determined how differences in the spatial pattern of the contractile cells participating in apical constriction affected invagination. We considered (i) circular and (ii) circumferential spatial patterns (Fig. 5A,B). To use the same cell proliferation condition in both patterns, the area in which the cell proliferation occurs was the same in both patterns. Therefore, the only difference between the two conditions was whether the pattern topologically contained a hole or not. This topology difference caused a large difference in invagination; therefore it has important implications in the deformation of epithelial tissue. Figure 5(A–C) shows that the directions of invagination were opposite. A basally convex shape emerged from the circular pattern, while an apically convex shape emerged from the circumferential pattern. Because the epithelium sheet has an apico‐basal polarity, the direction of invagination is very important for establishing organs. This result indicated that the tissue could be deformed to either direction solely as a function of the spatial pattern of the contractile cells. What remains unclear is exactly how the spatial pattern determines the invagination direction.

To investigate the differences in cell shapes between the two patterns, the shape of the contractile cells was quantified by the ratio of apical to basal area, λab, and is shown in Figure 5D. Under both conditions, λab decreased with time, indicating that the cells changed to an apically‐narrow wedge shape. Therefore, on the cellular scale, almost the same shape was obtained from either pattern. However, since the arrangement of such wedge‐shaped cells is different for the two patterns, we considered the relationship between the arrangement of the cells and the direction of the invagination. These cell shapes restrict how the surface of the tissue on the apical side connects. Because wedge‐like cells are located next to each other in the circular pattern, the tissue surface area on the apical side is smaller than on the basal side, causing the tissue shape to become basally convex (Fig. 6A). However, wedge‐like cells are not located next to each other along the radial axis in the circumferential pattern, while non‐contractile cells are located next to each other. The normal direction of the contact surface between the contractile and non‐contractile cells inside the pattern was inclined in the apical direction so as to conform to the shape of the contractile cell. Because the reaction force due to the cell proliferation occurred along this normal direction, the tissue inside the pattern showed an apically convex invagination (Fig. 6B).

Mode of apical constriction in invagination

In the previous section, we used surface contraction by apical constriction, in which the apical ring squeezed the apical surface of the cell. It is known, however, that there is another mode of apical constriction termed line contraction, which may occur locally on an edge between adjacent cells. In this section, we examine how differences in the apical constriction mode affect tissue invagination even when the topology of the spatial pattern of the contractile cells is equivalent. For the simulations, the circumferential pattern was used, and a contractile force was generated at the boundary line in contact with the cell inside the circumferential pattern. Cell proliferation occurred in the region enclosed by the circumferential pattern.

Figure 7(A,B) shows that the contractile line along the circumferential pattern produced a basally convex invagination, while surface contraction in the same pattern produced an apically convex invagination as shown in Figure 5(B). Therefore, the directions of invagination are opposite.

The shape of the contractile cells was quantified by the ratio of apical to basal area, λab, and is shown in Figure 7(C). λab increased with time in the case of the line contraction, indicating that the cells became basally narrow wedge‐like shapes; this is because centripetal motion is forced on arbitrary points on the circumference only on the apical side when the circumference shrinks. The normal direction of the bounding contact surface between the contractile and inner cells across the circumference was inclined in the basal direction (Fig. 8), which is opposite of the surface contraction mode. Quantitatively, we estimated the inclination angle of the contact surface between the contractile cell and the inner cell by measuring the cell height and apical and basal widths when the tissue invagination seemed to start (at t=1.0 cell cycle). The contact surface was inclined toward the apical side by about 12° in the case of surface contraction, while it was inclined toward the basal side by about 5° in the case of line contraction. Therefore, the reaction force due to the cell proliferation occurred along this normal direction, resulting in an apically (Fig. 6B) and a basally convex (Fig. 8) invagination in the case of surface or line contraction, respectively. Thus, even if the spatial pattern is the same, there is a possibility that the direction of invagination may change because of a different mode of apical constriction.

Discussion

In this study, we examined the mechanical roles of cell proliferation and apical constriction, and determined that they are responsible for expanding the degree of tissue deformation and determining the location of the deformation, respectively. We showed that the direction of invagination depended on the spatial pattern of the contractile cells. In addition, the direction of invagination differed between surface and line contraction modes of apical constriction, even when the spatial pattern was the same. The invagination direction can be understood from the alignment of wedge‐shaped cells.

From a physical point of view, surface contraction by apical constriction can also be expressed by surface energy proportional to apical area. By performing simulations of invagination using a surface energy function, we confirmed that the choice of energy function expressing the surface contraction did not affect the directionality of the invagination (Appendix I). This is because the surface energy function can be converted mathematically to the form of Equation (7) under the assumption that the apical surface is a regular hexagon. From a mechanical point of view, because buckling of tissue depends on tissue size, we examined whether the system size affected the directionality, and qualitatively confirmed that the result did not change (Appendix II).

Basal convex invagination was observed in the Drosophila spiracle (Simões et al. 2006), in which an enrichment of Myosin II and activated Rho 1 at the apical side occurred in a circular pattern. Basal convex invagination was also observed in the Drosophila tracheal placode, in which a circumferential pattern with line contraction was seen via accumulation of Myosin at the boundary of cells with high and low EGFR activity (Nishimura et al. 2007). During the optic cup formation, apical constriction in the hinge area along the circumference of the tissue and apically convex invagination has been observed (Eiraku et al. 2011). Our simulation results are consistent with these experimental results.

In our simulations, the cell division axis was assumed to follow the long axis of the cell (Hertwig's rule); therefore the tissue growth was isotropic. However, the mitotic axis may be controlled by the morphogen concentration gradient, for example, in the Drosophila wing imaginal disc, the mitotic spindle tends to be oriented along the proximal‐distal axis (Mao et al. 2011). The determination of the direction of invagination will not be affected by the orientation of the cell division axis if the cell division generates no bending moment in the tissue. Because the dividing nucleus is known to translocate apically (Ragkousi & Gibson 2014), the biased location of the nucleus may cause the generation of a bending moment in the tissue, possibly affecting the invagination direction.

From the mechanical point of view, because invagination due to cell proliferation can be interpreted as a type of sheet buckling due to the increase of in‐plane elastic stress, the direction of invagination is indeterminate. To determine the direction of the invagination, mechanical asymmetry can be required along the apical‐basal polarity. One possible source of asymmetry is an apically narrow wedge‐like cell shape induced by apical constriction as shown in this study. Another possible source of asymmetry is a difference in the mechanical properties between the apical and the basal sides of the epithelium. In addition, the surrounding mechanical properties are important for both the invagination direction and the tissue shape after invagination. In experiments on 3D tissue formation from Madin‐Darby canine kidney (MDCK) cells using an in vitro system, it was observed that the height of the protrusion structure (tulip‐hat‐like shape) decreased as the viscosity of the substrate increased (Imai et al. 2015). Investigation of the cell division orientation axis and the surrounding mechanical properties are future research directions.

A basal constriction induced the morphogenesis of a tulip‐hat‐like shape in the MDCK cell sheet (Imai et al. 2015). Although we did not incorporate basal constriction in this study, a comparison can be made between the in vitro experiments and our simulation results from a mechanical point of view. With active force generation, the tulip‐hat‐like shape protruded to the opposite side of the contractile side. Additionally, phosphorylated myosin regulatory light chains were observed on the outermost edge of the MDCK cell sheet, indicating a circumferential pattern. According to our simulation results, the invagination that satisfies these geometric conditions is the circumferential pattern with line contraction, which coincides with the experimental results.

Other mechanical effects can be attributed to the surroundings. When an invagination of bottle cells occurs during gastrulation in Xenopus embryos, their apical surfaces shrink mainly along the animal‐vegetal direction (Keller et al. 2003). However, bottle cells in isolated explants without adjacent tissues show an isotropic contraction and become circular apices. This experiment suggests that even if the intrinsic contraction is isotropic, an anisotropic deformation may occur because of mechanical influences from surrounding tissues. Throughout our study we assumed that all cells have the same physical properties and used a periodic boundary condition to enforce isotropy, therefore anisotropic deformation did not occur. To reproduce the actual in vivo phenomena, it is necessary to incorporate differences in boundary conditions and cell physical properties from tissue to tissue. How anisotropic deformation emerges from an isotropic spatial pattern can then by analyzed.

Finally, we would like to point out that actively contracting the cell surface through apical or basal constrictions is not the only way to trigger invagination. If the cell sheet is adhered to the matrix or another tissue and the elastic moduli of the surroundings are high enough, it is expected that cell height elongation would result in an invagination to the matrix side (Clausi & Brodland 1993; Inoue et al. 2016). For example, if an incompressible basement membrane is formed, cell elongation causes preferential reduction in the apical area of the cell to maintain the cell volume, and thus, the cell becomes an apically narrow wedge‐like shape, causing invagination. Furthermore, from the physical point of view, it is suggested that when the cell sheet proliferates, a folding pattern is generated because of the difference in growth rate between the cell sheet and the matrix (Hannezo et al. 2011; Amar & Jia 2013; Tallinen et al. 2016). Such a folding pattern is of interest for understanding intestinal villi and the wrinkle structure of the brain.

Funding

Y. I. acknowledges support from MEXT KAKENHI (15H05861, 15H05856). Y. I. and T. A. acknowledge support from the Platform Project for Supporting in Drug Discovery and Life Science Research (Platform for Dynamic Approaches to Living System) from the Japan Agency for Medical Research and Development (AMED).