Implementation of cellular bulk stresses in vertex models of biological tissues

Abstract

abstract

Vertex models describe biological tissues as tilings of polygons. In standard vertex models, the tissue dynamics result from a balance between isotropic stresses, which are associated with the bulk of the cells, and tensions associated with cell–cell interfaces. However, in this framework it is less obvious how to describe anisotropic stresses arising from the bulk of cells. In epithelia, such bulk anisotropic stresses could arise for instance through medial myosin fluctuations. Two recent publications—Tlili et al. (Proc Natl Acad Sci USA 116(51):25430–25439, 2019) and Comelles et al. (eLife 10:e57730, 2021)—have proposed different schemes to implement bulk anisotropic stresses in vertex models. Here we show that while both schemes transform in the same way under affine deformations, they lead to significantly different tissue dynamics. Our results are consistent with the interpretation that the Tilli et al. scheme describes bulk stresses that are uniform within each cell, while the Comelles et al. scheme corresponds to non-uniform bulk stresses. Finally, we wondered whether a standard vertex model can be fully expressed in terms of bulk cellular stresses alone. We find that, in general, neither scheme can mimic the vertex forces created by cell–cell interface tensions

Graphic abstract

Appendices

Appendix A: Cell and tissue aspect ratio change

In Figs. 2 and 3, we have quantified the relative cell aspect ratio change as $\varepsilon _{\mathrm {cell}} = \mathrm{AR}_\mathrm{final} / \mathrm{AR}_\mathrm{initial} - 1$, where $\mathrm AR_{initial}$ ($\mathrm AR_{final}$) is the aspect ratio of the cell in the initial (final) state. $\mathrm AR$ is defined as $\mathrm{AR} = \sqrt{{I_1}/{I_2}}$, where $I_1$ and $I_2$ ($I_1 > I_2$) are the two eigenvalues of the matrix $\varvec{M}$ defined in Eq. (11).

Similarly, for free boundary conditions, we analogously define the relative tissue aspect ratio change $\varepsilon _{\mathrm {tissue}}$. In this case however, the tensor $\varvec{M}$ is defined via Eq. (11) using all margin vertices of the tissue.

Table 1 List of parameter values used in vertex model simulation

Appendix B: Consistency checks

1.1 Appendix B.1: Zero net force and torque across each cell

Zero net force We check that the net force on all vertices of a cell satisfies the relation

$$\begin{aligned} \varvec{F}_G&= \sum \limits _{i\in cell}{\varvec{F}_{i}^{\left( T,C \right) }} = \varvec{0}, \end{aligned}$$

(B.1)

in both schemes. Indeed, substituting Tlili et al. forces Eq. (4) into Eq. (B.1), we find that

$$\begin{aligned} \varvec{F}_G =\frac{1}{2}\left[ \varvec{\hat{z}}\times \sum \limits _{i\in cell}{\left( {{\varvec{r}}_{i+1}}-{{\varvec{r}}_{i-1}} \right) } \right] \cdot \varvec{\sigma }^{(\mathrm {b})}=\varvec{0}, \end{aligned}$$

(B.2)

Likewise, substitution of Comelles et al. definitions Eq. (13) into Eq. (B.1) leads to

$$\begin{aligned} \sum \limits _{i\in cell}{\varvec{F}_{i}^{\left( C \right) }}&= -A\left( \sum \limits _{i\in cell}{{{\varvec{\rho }}_{i}}} \right) \cdot {{\varvec{M}}^{-1}}\cdot \varvec{\sigma }^{(\mathrm {b})}=\varvec{0}. \end{aligned}$$

(B.3)

Zero net torque We also check that the overall torque acting on a given cell satisfies the relation

$$\begin{aligned} \varvec{T} = \sum \limits _{i\in cell}{{{\varvec{r}}_{i}}\times \varvec{F}_{i}^{\left( T,C \right) }} = \mathbf{0} \end{aligned}$$

(B.4)

in both schemes. Substituting Tlili et al. scheme definition Eq. (4) into Eq. (B.4), and adopting Einstein notation with Greek dimension indices, we obtain:

$$\begin{aligned} T_\alpha&=-\frac{1}{2}\sum \limits _{i\in cell}{ \varepsilon _{\alpha \beta \gamma }{r}_{i,\beta }\varepsilon _{\eta z\mu } ({r}_{i+1,\eta }-{r}_{i-1,\eta })\sigma _{\mu \gamma }^{(b)}}. \end{aligned}$$

(B.5)

To further simplify this expression, we note that:

$$\begin{aligned} \sum \limits _{i\in cell}{ {r}_{i,\beta }\varepsilon _{\eta z\mu } ({r}_{i+1,\eta }-{r}_{i-1,\eta }) } = 2A\delta _{\beta \mu }. \end{aligned}$$

(B.6)

This relation can be verified by testing each dimension index combination for the pair $(\beta ,\mu )$. Using Eq. (B.6) in Eq. (B.5), we find that

$$\begin{aligned} T_\alpha&=-A\varepsilon _{\alpha \beta \gamma }\sigma _{\beta \gamma }^{(b)}. \end{aligned}$$

(B.7)

For symmetric bulk stress tensors $\varvec{\sigma }^{(\mathrm {b})}$, we thus obtain $\varvec{T}=\varvec{0}$.

Substitution of Comelles et al. force definition Eq. (13) into Eq. (B.4), adopting index notation, yields

$$\begin{aligned} T_\alpha&=-A\varepsilon _{\alpha \beta \gamma }M^{-1}_{\eta \mu }\sigma ^{(b)}_{\mu \gamma } \sum \limits _{i\in cell}{ r_{i,\beta }\rho _{i,\eta } }, \nonumber \\&=-A\varepsilon _{\alpha \beta \gamma }M^{-1}_{\eta \mu }\sigma ^{(b)}_{\mu \gamma } \sum \limits _{i\in cell}{ \rho _{i,\beta }\rho _{i,\eta } }, \nonumber \\&=-A\varepsilon _{\alpha \beta \gamma }M_{\beta \eta }M^{-1}_{\eta \mu }\sigma ^{(b)}_{\mu \gamma }, \nonumber \\&=-A\varepsilon _{\alpha \beta \gamma }\sigma ^{(b)}_{\beta \gamma }. \end{aligned}$$

(B.8)

Hence, again, for symmetric $\varvec{\sigma }^{(\mathrm {b})}$, the overall torque is zero: $\varvec{T}=\varvec{0}$.

1.2 Appendix B.2: Self-consistency: the input bulk stress equals the output Batchelor stress

Here we show that in both schemes, Batchelor stress Eq. (15) is identical to the input bulk cell-stress $\varvec{\sigma }^{(\mathrm {b})}$, i.e. that: $\varvec{\sigma }_B\left[ \varvec{F}_{i}^{\left( T,C \right) } \right] =\varvec{\sigma }^{(\mathrm {b})}$.

Substituting Tlili et al. force expression Eq.  (4) into Eq. (15) we immediately obtain:

$$\begin{aligned} \varvec{\sigma }_B\left[ \varvec{F}_{i}^{\left( T \right) } \right]&=\frac{1}{2A}\sum \limits _{i\in cell}{{{\varvec{r}}_{i}}\otimes \left[ \left( {{\varvec{r}}_{i+1}}-{{\varvec{r}}_{i-1}} \right) \times \varvec{\hat{z}} \right] \cdot \varvec{\sigma }^{(\mathrm {b})}} \nonumber \\&= \varvec{\sigma }^{(\mathrm {b})}, \end{aligned}$$

(B.9)

where in the second step, we have applied Eq. (B.6).

Similarly, substitution of Comelles et al. force definition Eq. (13) into Eq. (15) leads to

$$\begin{aligned} \varvec{\sigma }_B\left[ \varvec{F}_{i}^{\left( C \right) } \right]&= -\frac{1}{A}\sum \limits _{i\in cell}{{{\varvec{r}}_{i}}\otimes \left[ -A{{\varvec{\rho }}_{i}}\cdot {{\varvec{M}}^{-1}}\cdot \varvec{\sigma }^{(\mathrm {b})}\right] } \nonumber \\&= \left( \sum \limits _{i\in cell}{{{\varvec{r}}_{i}}\otimes {{\varvec{\rho }}_{i}}} \right) \cdot {{\varvec{M}}^{-1}}\cdot \varvec{\sigma }^{(\mathrm {b})}, \nonumber \\&= \left( \sum \limits _{i\in cell}{{{\varvec{\rho }}_{i}}\otimes {{\varvec{\rho }}_{i}}} \right) \cdot {{\varvec{M}}^{-1}}\cdot \varvec{\sigma }^{(\mathrm {b})}\nonumber \\&= \varvec{\sigma }^{(\mathrm {b})}. \end{aligned}$$

(B.10)

with the matrix $\varvec{M}$ defined in Eq. (11).

Appendix C: Capacity of both schemes to reproduce standard vertex model forces

1.1 Appendix C.1: Identity of forces for regular polygonal cells

Here we consider a regular polygonal cell composed of N vertices $\{ \varvec{r}_i \}_{i = 1, \dots , N}$ with ${{\varvec{r}}_{i}}=r{{\varvec{e}}_{i}}$, with $r >0$, ${{\varvec{e}}_{i}}=\left( \cos {{\theta }_{i}},\sin {{\theta }_{i}} \right)$, and the vertex angles ${\theta }_{i}=2\left( i-1 \right) \pi /N$. Area and perimeter of such a cell are:

$$\begin{aligned} A&=\frac{1}{2}N{{r}^{2}}\sin \left( \dfrac{\text {2 }\!\!\pi \!\!\text { }}{N}\right) , \end{aligned}$$

(C.11)

$$\begin{aligned} P&= 2Nr \sin \left( \dfrac{\pi }{N}\right) . \end{aligned}$$

(C.12)

To simplify standard vertex model forces Eq. (35), we use the relations

$$\begin{aligned} \left( \varvec{r}_{i+1}-\varvec{r}_{i-1} \right) \times \varvec{\hat{z}}&= 2r\sin \left( \frac{2\text { }\!\!\pi \!\!\text { }}{N}\right) {{\varvec{e}}_{i}} \end{aligned}$$

(C.13)

$$\begin{aligned} {{\varvec{t}}_{i,i+1}}-{{\varvec{t}}_{i-1,i}}&= -2\sin \left( \frac{\pi }{N}\right) \varvec{e}_i \end{aligned}$$

(C.14)

and obtain:

$$\begin{aligned} \varvec{F}^{(\mathrm {svm})}_i =&-K_A\left( A-{{A}_{0}} \right) r\sin \left( \frac{2\text { }\!\!\pi \!\!\text { }}{N}\right) {{\varvec{e}}_{i}} \nonumber \\&-2K_P\left( P-P_0\right) \sin \left( \frac{\text { }\!\!\pi \!\!\text { }}{N}\right) {{\varvec{e}}_{i}}. \end{aligned}$$

(C.15)

Fig. 4

Differences between the forces according to both schemes and the standard vertex model forces for irregular cell shapes. Vertex positions correspond to a randomly perturbed hexagonal cell according to Eq. (C.17). a The Tlili et al. forces $\varvec{F}_i^{(\mathrm {T})}$ (green arrows), the Comelles et al. forces $\varvec{F}_i^{(\mathrm {C})}$ (blue arrows), and the standard vertex model forces $\varvec{F}^{(\mathrm {svm})}_i$ (red arrows) at each cell vertex for an example configuration. b Dependence of the relative force deviations $\zeta _\mathrm{T}$ and $\zeta _\mathrm{C}$ defined in Eq. (C.18) on the amplitude $\eta$ of the deviation of the cell shape from a regular hexagonal shape. Vertex model parameters are provided in Table 1

The Tlili et al. vertex forces in Eq. (38) simplify since we consider regular hexagons and thus $\varvec{S} = \varvec{I} / 2$. Inserting Eqs. (C.11), (C.12) and (C.13) in Eq. (36), we recover indeed an expression identical to the one in Eq.  (C.15).

To compute the Comelles et al. forces, we fist evaluate each component of the $\varvec{M}$ tensor defined in Eq. (11). We find for $N\ge 3$:

$$\begin{aligned} {{M}_{xx}}&={{r}^{2}}\sum \limits _{i=1}^{N}{{{\cos }^{2}}{{\theta }_{i}}}=\frac{1}{2}{{r}^{2}}\sum \limits _{i=1}^{N}{\left( 1+\cos 2{{\theta }_{i}} \right) }=\frac{N}{2}{{r}^{2}}, \\ {{M}_{yy}}&={{r}^{2}}\sum \limits _{i=1}^{N}{{{\sin }^{2}}{{\theta }_{i}}}=\frac{1}{2}{{r}^{2}}\sum \limits _{i=1}^{N}{\left( 1-\cos 2{{\theta }_{i}} \right) }=\frac{N}{2}{{r}^{2}}, \end{aligned}$$

while

$$\begin{aligned} {{M}_{xy}}&={{r}^{2}}\sum \limits _{i=1}^{N}{\cos {{\theta }_{i}}\sin {{\theta }_{i}}} =\frac{1}{2}{{r}^{2}}\sum \limits _{i=1}^{N}{\sin 2{{\theta }_{i}}}=0, \end{aligned}$$

which leads to $\varvec{M}=N{{r}^{2}}/2 \ \varvec{I}$. Therefore, the vertex forces given by the Comelles et al. scheme are

$$\begin{aligned} \varvec{F}_{i}^{\left( \text {C} \right) } =&-r\sin \left( \frac{\text {2 }\!\!\pi \!\!\text { }}{N}\right) {{\varvec{e}}_{i}}\cdot \varvec{\sigma }_B\left[ \varvec{F}_{i}^{{(\mathrm {svm})}}\right] . \end{aligned}$$

(C.16)

Injecting stress expression Eq. (36) into Eq. (C.16), we find that the Comelles et al. expression is identical to Eq. (C.15); thus $\varvec{F}_i^\mathrm{(T)} = \varvec{F}_i^\mathrm{(C)} = \varvec{F}^{(\mathrm {svm})}_i$ for regular polygonal cells.

1.2 Appendix C.2: Numerical comparison of forces for irregular hexagons

We numerically compare the Tlili et al. and Comelles et al. forces with the standard vertex model forces for non-regular hexagonal cells. We define the corners of the hexagonal cell as $\varvec{r}_i=(x_i,y_i)$ with $i=1,2,\cdots ,6$ and:

$$\begin{aligned} {{x}_{i}}=\cos \frac{i\text { }\!\!\pi \!\!\text { }}{3}+\eta {{\vartheta }_{ix}}, \qquad {{y}_{i}}=\sin \frac{i\text { }\!\!\pi \!\!\text { }}{3}+\eta {{\vartheta }_{iy}}. \end{aligned}$$

(C.17)

Here, $\vartheta _{ix}$ and $\vartheta _{iy}$ are independent, zero-mean, unit-variance Gaussian random variables, and $\eta$ is a parameter tuning the deviation of the cell shape from that of a regular hexagon.

Discrepancies between the two schemes and the standard vertex model forces are quantified in terms of

$$\begin{aligned} {{\zeta }_{X}}=\left\langle \left| \varvec{F}_{i}^{\left( X \right) }-{{\varvec{F}}^{(\mathrm {svm})}_{i}} \right| /\left| {{\varvec{F}}^{(\mathrm {svm})}_{i}} \right| \right\rangle , \end{aligned}$$

(C.18)

where $X \in \lbrace T,C\rbrace$, indicating Tlili and Comelles et al. schemes, respectively. The average $\langle \cdot \rangle$ is over all cell vertices and over 1, 000 realizations for the set of the Gaussian random variables $\vartheta _{ix}$ and $\vartheta _{iy}$.

We find that the deviations of both Tlili et al. and Comelles et al. forces from the standard vertex model forces increase with increasing cell shape deviation from a regular hexagonal shape, see Fig. 4b. However, the Tlili et al. forces deviated less from the standard vertex model forces than the Comelles et al. forces did.