1 Introduction

During morphogenesis, complex structures emerge starting from a single fertilized egg as the results of the collective organization of a large number of cells. Understanding principles that govern self-organization of cells into complex structures and organs is one of the major challenges of biology and biophysics. The collective behavior of cells relies on chemical signals [1,2,3,4], but also depends on cellular force generation and active mechanical processes as well as tissue mechanical properties [5, 6]. Morphogenesis, i.e., the generation of shape, is therefore a result of self-organized processes that couple chemical signalling with mechanical activity [7,8,9,10].

Change in tissue shape involves anisotropic active processes and cell rearrangements. The physics of tissue dynamics is based on a description of tissues as active viscoelastic and viscoplastic materials [11,12,13,14]. Depending on timescales, tissues can behave like solids, able to withstand external shear stresses, or like fluids, and can rearrange their cells and exhibit cell flows [12,13,14]. Such rearrangements permit the maintenance of mechanical integrity of a tissue while changing local connectivity and the overall shape. In tissues, rearrangements can result from cell divisions or extrusions, where new cells are added or removed from the tissue, which has been shown to permit tissue fluidization [15]. In addition, cells can also rearrange and change neighbors in so-called T1 transitions.

T1 transitions have been studied first in passive materials such as foams, where they occur in response to external shear forces that can drive material flow [16, 17]. Tissues, however, are active materials, that can deform spontaneously, driven by internally generated stresses, and can therefore also perform work on their environment. Such active deformations are for instance observed during convergence-extension, a widespread morphogenetic process driven by oriented T1 transitions that leads to anisotropic tissue deformation [18, 19]. In contrast to passive foams, where T1 transitions dissipate energy and relax elastic stresses resulting from external forcing, T1 transitions in tissues can be active and perform work, and therefore can build up stresses rather than relaxing them. The orientation of T1 transitions can be guided by tissue polarity cues, that are linked to chemical signals such as the planar polarity pathways [20, 21]. Thereby, tissues can extend along axes that are defined by chemical patterns. Such processes can be observed in developmental model systems. For example, during the germ-band extension of Drosophila embryo, experiments show that the tissue deforms anisotropically as a consequence of oriented T1 transitions, driven by active processes in the acto-myosin cytoskeleton. Two scenarios have been proposed: anisotropic accumulation of myosin II at cell–cell junctions [22,23,24,25], and anisotropic active stresses mediated by medial myosin pulses [26]. Similarly, data from the pupal wing blade of Drosophila reveal multiple roles of T1 transitions over time: they drive anisotropic cell and tissue elongation at early stages, while later they are responsible for a relaxation of cell shape elongation [8]. Finally, cells can actively propel themselves on a substrate. Such motion can also cause T1 transitions and tissue shape changes [27,28,29], and can trigger a solid-to-fluid transition in cell tissues model [30]. Here, we focus of anisotropic stresses generated in a tissue in the absence of active self-propulsion.

We use a two-dimensional vertex model to study how anisotropic tissue stresses can drive oriented cell rearrangements and anisotropic tissue shape changes. In particular, we discuss active T1 transitions that can perform mechanical work, in contrast to passive T1 transitions that relax elastic stresses. Vertex models provide simple models of tissue physics that can capture cell shape, packing geometry and cell rearrangements [8, 30,31,32,33,34,35,36]. The role of anisotropic internal stresses has been studied within the vertex model framework, for instance using anisotropy in cell bond tension [23, 25], or by introducing cell bond tension that depends on the identity of adjacent cells [24]. Following Ref. [37], we consider a cell network where a preferred axis is set by a nematic field that represents tissue polarity. We discuss two different physical realizations of anisotropic active stress: (i) anisotropic bond tension, where the contractility of bonds is increased along a preferred axis, and (ii) anisotropic cell stress, where the bulk of the cells exhibits an anisotropic stress that is contractile along a preferred axis. Surprisingly, we find that these two realizations, although involving an anisotropic active stress along the same direction, lead to cell rearrangements and cell elongation patterns which are very different. Our analysis suggests that the early stages of Drosophila pupal wing morphogenesis could be an example where anisotropic bond tension dominates, while both realizations of anisotropic active stress could contribute during germ-band extension in the Drosophila embryo. In the latter case, anisotropic cell stress could be a consequence of medial myosin II pulses that favors the opening of newly-formed cell bonds [26]. We complement our analysis and understanding by using a coarse-grained continuum description of the cell network. This description uses concepts from active matter theory [38,39,40] and has proven valuable to characterize the large-scale properties of tissues [41,42,43,44,45]. By considering the energetics of the cell network, we show that active stresses can induce T1 transitions that perform mechanical work. We call these active T1 transitions.

The paper is organized as follows. In Sect. 1, we present the vertex model and its modification to account for anisotropic bond tension and cell stress. We then introduce a linear anisotropic continuum model to capture tissue dynamics. In Sect. 2, we quantify the outcome of anisotropic vertex model simulations and highlight the differences between the two implementation of anisotropy. Fits of the continuum model to the simulation results provides us with a better understanding of the mechanisms at play. We finally discuss the energetics of the tissues, allowing us to provide a definition of active T1 transitions.

2 Mechanics of anisotropic cell networks

The apical junctions of an epithelial tissue can be described by a packing of convex polygons and its mechanics can be described by a vertex model, where cells are represented as polygons that are outlined by straight edges connecting vertices [46]. We consider a polygonal cell network consisting of \(N_{\mathrm{c}}\) cells. Each cell \(\alpha \) is characterized in terms of its area \(A^\alpha \), its perimeter \(L^\alpha \) and the lengths \({\mathcal {L}}_{mn}\) of the bonds that form the outline of the cell, where m and n label the vertices that they connect (see Fig. 1 for illustration).

We employ a quasistatic representation of epithelia where the cell network is at any instant in a mechanical equilibrium, while the parameters describing cell properties can slowly change with time. At each vertex m, the total force \({\mathbf {F}}_m = -\partial E_0/\partial {\mathbf {R}}_m\) vanishes, where \({\mathbf {R}}_m\) is the position of the vertex, and \(E_0\) is the vertex model work function and reads [46, 47]:

$$\begin{aligned} E_0= & {} \sum _\alpha \frac{1}{2} K^\alpha \left( A^\alpha - A^\alpha _0 \right) ^2 \nonumber \\&+ \sum _{\langle m,n \rangle } \Lambda _{mn} {\mathcal {L}}_{mn}+ \sum _\alpha \frac{1}{2} \Gamma ^\alpha (L^\alpha )^2. \end{aligned}$$
(1)

Note that for clarity, upper-case letters are used here and in the following for quantities related to the vertex model, while lower-case letters will be used for the continuum model. The first term describes an area elasticity contribution, with \(A^\alpha _0\) the preferred cell area and \(K^\alpha \) the area stiffness. The second term describes a contribution due to the tension of network bonds with length \({\mathcal {L}}_{mn}\) and line tension \(\Lambda _{mn}\). The third term describes an elasticity of the cell perimeter with stiffness \(\Gamma ^\alpha \).

Fig. 1
figure 1

Mechanics and dynamics of cellular networks. A Definition of the cell state variables. Left shows the cell area \(A^\alpha \) (blue patch), cell perimeter \(L^\alpha \) (green line) and bond length \({\mathcal {L}}_{mn}\) (red line) between the vertices with positions \({\mathbf {R}}_m\) and \({\mathbf {R}}_{n}\). Right shows the cell elongation tensor \({\mathbf {G}}\) which is constructed from the bond nematic tensors \({\mathbf {H}}_{mn}\) as defined in Eq. (6). B Cell dynamic processes can lead to tissue deformation as an effect of cell shape changes, T1 transitions, cell divisions or cell extrusions. C Large-scale tissue deformation can be driven by collective cell dynamics: cell shape changes (top), anisotropic T1 transitions (middle) and anisotropic cell divisions (bottom). The tissue may also deform as a result of changes in the mean cell shape of the cellular network

Full size image

Nonequilibrium dynamics of the vertex model are captured by a time-dependent line tension \(\Lambda _{mn}(t)\). The line tension dynamics of individual bonds in the network follows an Ornstein–Uhlenbeck process:

$$\begin{aligned} \frac{\mathrm{d}\Lambda _{mn}}{\mathrm{d}t} = -\frac{1}{\tau _\Lambda }(\Lambda _{mn}(t) - {\bar{\Lambda }}_{mn}) + \Delta \Lambda \sqrt{2/\tau _\Lambda } \, \Xi _{mn}(t) \, , \nonumber \\ \end{aligned}$$
(2)

where \(\Xi _{mn}(t)\) is a Gaussian white noise with zero mean \(\langle \Xi _{mn}(t) \rangle = 0\), and correlations \(\langle \Xi _{mn}(t)\Xi _{op}(t') \rangle = \delta _{\langle mn \rangle , \langle op \rangle } \delta (t-t')\) where \(\delta _{\langle mn \rangle , \langle op \rangle }=1\) if bonds \(\langle mn \rangle \) and \(\langle op \rangle \) are the same and 0 otherwise [46, 48]. The line tension of every bond relaxes toward its mean value \({\bar{\Lambda }}_{mn}\) with a characteristic time \(\tau _\Lambda \), which sets the timescale of the dynamics and is of the order of the acto-myosin cortex turn-over time.

As discussed for instance in Refs. [36, 37], the magnitude of bond tension fluctuations \(\Delta \Lambda \) has a crucial role in the rheological properties of cell networks. A low value of this fluctuation magnitude leads to a glassy dynamics and non-linearities dominate. In the following, we are interested in a regime where \(\Delta \Lambda \) is sufficiently large, such that the vertex model has linear viscoelastic properties.

A polygonal network described by the work function (1) has isotropic properties. In the following, we consider how this description can be extended to describe anisotropic cell networks.

2.1 Anisotropic cellular networks

Motivated by planar cell polarity in tissues [20, 21], we consider that the anisotropy of the network can be described by a unit nematic field \(\varvec{{\mathcal {P}}}\) assigned to each polygon and which gives locally a preferred axis. In two dimensions, the nematic field \(\varvec{{\mathcal {P}}}\) can be parameterized by a single angle \(\Psi \) which defines the direction of the anisotropy axis (see App. B4 of the Supplementary Material for details). For simplicity, we consider in the following that \(\Psi \) is constant and thus provides a global preferred axis in the tissue.

Anisotropic bond tension. To include the effect of such a nematic field in the dynamics of the vertex model, we first consider an anisotropic bond tension, such that its mean magnitude \({\bar{\Lambda }}_{mn}\) depends on the orientation of the bond with respect to the nematic field \({\mathcal {P}}\). We choose

$$\begin{aligned} {\bar{\Lambda }}_{mn} = {\bar{\Lambda }}_{mn}^0 \left( 1 + \beta \varvec{{\mathcal {P}}}: \hat{{\mathbf {H}}}_{mn} \right) \, , \end{aligned}$$
(3)

where \(\varvec{A}:\varvec{B}= \mathrm{Tr}(\varvec{A}\cdot \varvec{B})\) denotes the full tensor contraction, and where \(\beta \ge 0\) is dimensionless and sets the magnitude of the anisotropy. We have introduced the unit bond nematic \(\hat{{\mathbf {H}}}_{mn}={\mathbf {H}}_{mn}/{\mathcal {L}}_{mn}^2\) between vertices m and n. Here, \({\mathbf {H}}_{mn}\) is the nematic tensor constructed from the vector \(\varvec{{\mathcal {L}}}_{mn}\) pointing from vertex m to vertex n as \({\mathbf {H}}_{mn}=\varvec{{\mathcal {L}}}_{mn} \otimes \varvec{{\mathcal {L}}}_{mn} - {\mathcal {L}}_{mn}^2 \mathbb {1}/2\) with \(\mathbb {1}\) the unit tensor in two dimensions (see Fig. 1). With this definition, a bond making an angle \(\Phi \) with the local nematic field has a mean bond tension that reads:

$$\begin{aligned} {\bar{\Lambda }}_{mn} = {\bar{\Lambda }}_{mn}^0 \left( 1 + \beta \cos (2\Phi ) \right) \, . \end{aligned}$$
(4)

Consistently with our convention for the anisotropic cell stress, the definition given by Eq. (3) with \(\beta >0\) implies that cells have a higher bond tension along the axis set by \(\varvec{{\mathcal {P}}}\). As a consequence, they are more likely to undergo a T1 transition along this axis.

Anisotropic cell stress. An alternative description of anisotropic tissues can be obtained by considering an anisotropic cell stress \(\varvec{\Sigma }^{\mathrm {a}} = \Sigma ^{\mathrm {a}} \varvec{{\mathcal {P}}}\), where \(\Sigma ^{\mathrm {a}}\) is the magnitude of the active stress. The work performed by this anisotropic stress is added to the vertex model work function asFootnote 1:

$$\begin{aligned} E = E_0 + \sum _\alpha \frac{1}{2}A^\alpha \varvec{\Sigma }^{\mathrm {a}} : \varvec{G}^\alpha \, , \end{aligned}$$
(5)

where \(\varvec{G}^\alpha \) is the cell shape tensor of each cell \(\alpha \)

$$\begin{aligned} \varvec{G}^\alpha = \frac{1}{A^\alpha }\sum _{\langle m,n \rangle } {\mathbf {H}}_{mn} \, . \end{aligned}$$
(6)

The cell shape tensor \(\varvec{G}^\alpha \) quantifies the deviation of the cell shape from isotropic shapes, for which \(\varvec{G}^\alpha \) vanishes. See Fig. 1 for an illustration. With the definition of Eq. (5) and \(\Sigma ^{\mathrm{a}}>0\), the anisotropic cell stress implies a stronger contractility of the cells along the direction set by \(\varvec{{\mathcal {P}}}\), and cells therefore tend to elongate in the direction orthogonal to \(\varvec{{\mathcal {P}}}\).

Note that for simplicity, we use in the following the same constant values of the parameters \(K^\alpha \), \(A^\alpha _0\), \(\Gamma ^\alpha \) for all cells and the same value \({\bar{\Lambda }}_{mn}^0={\bar{\Lambda }}^0\) for all bonds. Note that the presence of bond tension fluctuations prevents crystalization of the cellular pattern that could be otherwise observed in a monodisperse system [49, 50]. In App. A of the Supplementary Material, we give details on the numerical implementation of the vertex model. Values of the (dimensionless) parameters used in the simulations are given in Table S1 of the Supplementary Material.

2.2 Dynamics of a polygonal cell network and shear decomposition

The deformation of a cellular network is quantified by its shear rate, which can be decomposed into cellular contributions. For flat polygonal networks, such a decomposition can be done exactly [8, 51]. Following Ref. [51], the large-scale shear-rate tensor \({\widetilde{V}}_{ij}\) of the cellular network can be decomposed as:

$$\begin{aligned} {\widetilde{V}}_{ij} = \frac{\mathrm{D} Q_{ij}}{\mathrm{D}t} + R_{ij} \, . \end{aligned}$$
(7)

Here and in the following, i and j corresponds to 2d Cartesian indices, \(Q_{ij}\) is the mean cell elongation tensor and \(\mathrm{D}/\mathrm{D}t\) is the corotational time derivative of a tensor (defined in Eq. (B3) of App. B of the Supplementary Material). The tensor \(R_{ij}\) accounts for shear rate due to topological rearrangements and is a sum of four contributions:

$$\begin{aligned} R_{ij} = T_{ij} + C_{ij} + E_{ij} + D_{ij} \, , \end{aligned}$$
(8)

where the tensors \(T_{ij}\), \(C_{ij}\) and \(E_{ij}\) account for shear rate due to T1 transitions, cell divisions and cell extrusions, respectively. The tensor \(D_{ij}\) is a shear rate associated with heterogeneities and fluctuations. If such fluctuations are correlated, they contribute to shear even if they vanish on average. In particular, the tensor \(D_{ij}\) includes shear stemming from correlations between triangle rotations and triangle elongation as well as correlations between triangle area changes and triangle elongation [8, 51]. Note finally that all the tensors introduced in Eqs. (7) and (8) are two-dimensional nematic tensors. It means that they are symmetric traceless tensors that define an orientation and a magnitude. They are fully characterized by two independent quantities: a norm and an angle with respect to the x-axis (see also App. B4 of the Supplementary Material).

Note that the trace of the velocity gradient tensor \(V_{kk}\) (summation over repeated indices is implied), which corresponds to isotropic tissue growth, can also be decomposed into cellular contributions [41, 51]. Here, we only focus on the anisotropic contributions. Finally, the tissue stress tensor \(\Sigma _{ij}\) in the simulations is symmetric and can be decomposed into an isotropic pressure and a symmetric traceless part, the shear stress \({\widetilde{\Sigma }}_{ij}\).

2.3 Hydrodynamic model for cellular networks under anisotropic active stress

The viscoelastic behavior of stochastic cellular networks can be captured by a continuum model of tissues [8, 37, 41, 52]. Such a coarse-grained description does not hold at a single-cell level but requires an averaging over many cells, as provided by the shear decomposition (7) of a triangulated network discussed above.

In the continuum description, we therefore introduce the anisotropic part of the deformation rate tensor \({\tilde{v}}_{ij}\), which can be decomposed into cellular contributions due to changes in the mean cell elongation tensor \(q_{ij}\) and shear \(r_{ij}\) caused by topological rearrangements. Note that we use lower-case letters for the continuum model description. We therefore have:

$$\begin{aligned} {\tilde{v}}_{ij} = \frac{\mathrm{D}q_{ij}}{\mathrm{D}t} + r_{ij} \, , \end{aligned}$$
(9a)

where \(\mathrm{D}/\mathrm{D}t\) denotes the corotational derivative defined in Eq. (B3) of the Supplementary Material.

We also include in our continuum description the fact that the axis of topological rearrangements is biased both by the axis of cell elongation and the axis of active anisotropic processes. This fact is captured by introducing linear relationships between the shear contribution from topological rearrangements \(r_{ij}\), the cell elongation \(q_{ij}\), and the anisotropic axis \(p_{ij}\). It reads [37, 41]:

$$\begin{aligned} r_{ij} = \frac{1}{\tau } \, q_{ij} + \lambda p_{ij} \, , \end{aligned}$$
(9b)

where \(\tau \) is the characteristic timescale of topological rearrangements and \(\lambda \) is the rate of anisotropic cell rearrangements.

We also introduce the tissue stress \(\sigma _{ij}\), which we decompose into an isotropic part and an anisotropic symmetric traceless part, the tissue shear stress \({\tilde{\sigma }}_{ij}\). We consider that the cellular network is an active elastic material, such that to linear order, the shear stress can be written as:

$$\begin{aligned} {\tilde{\sigma }}_{ij} = \mu \, q_{ij} + \zeta p_{ij} \, , \end{aligned}$$
(9c)

where \(\mu \) is the shear modulus of the tissue and \(\zeta \) is the anisotropic active stress magnitude.

3 Cell elongation and T1 transitions driven by active processes

We now discuss the role of anisotropy on the vertex model dynamics. For this purpose, we study the relaxation of the vertex model from an isotropic disordered steady state to an anisotropic steady state (see App. A of the Supplementary Material for details of the simulations). To understand the transient dynamics between these two steady states, we consider in the following a gradual activation of the anisotropic cell stress or of the anisotropic bond tension, given by:

$$\begin{aligned} \Sigma ^{\mathrm{a }}(t)&= \Sigma ^{\mathrm{a}}_0 \left( 1-\mathrm{e}^{-t/T_{\mathrm{a}}} \right) \Theta (t) \, , \nonumber \\ \beta (t)&=\beta _0 \left( 1-\mathrm{e}^{-t/T_{\mathrm{a}}} \right) \Theta (t) \, , \end{aligned}$$
(10)

where \(\Theta (t<0)=0\) and \( \Theta (t\ge 0)=1\). We have introduced an activation time \(T_{\mathrm{a}}\), and \(\Sigma ^{\mathrm{a}}_0\) and \(\beta _0\) are the steady-state anisotropic stress magnitude and bond tension magnitude, respectively. We have also considered an instantaneous activation where \(\Sigma ^{\mathrm{a }}(t)=\Sigma ^{\mathrm{a}}_0 \Theta (t)\) and \(\beta (t)=\beta _0 \Theta (t)\). This case is presented in App. C of the Supplementary Material.

In the following, we consider vertex model simulations with two types of boundary conditions: (i) fixed box boundary condition, for which the box size is fixed and the total tissue shear rate \(\widetilde{V}_{ij}\) vanishes; (ii) stress-free boundary condition, for which the total stress on the simulation box vanish and \(\Sigma _{ij}=0\), such that cells can rearrange and flow. Examples of realization of these two types of boundary conditions are shown in Movies 1 to 4.

3.1 T1 transitions driven by anisotropic bond tensions

We now focus on the case of anisotropic bond tension (see Eq. (3)) with a nematic tensor \(\varvec{{\mathcal {P}}}\) aligned with the vertical axis, such that bonds are more contractile along this direction. In the continuum model, we translate this choice by taking \(p_{xx}=-1\), \(p_{yy}=1\) and \(p_{xy}=p_{yx}=0\).

The left panel of Fig. 3 displays the outcome of vertex model simulations in the case of fixed box boundary conditions. See Movie 1 for an example of a vertex model simulation. In this case, the larger contractility of bonds along the y axis biases active T1 transitions along the same direction: bonds preferentially close along the y axis and new bonds are opened along the x axis. It results in a positive rate of rearrangements along the x axis: \(R_{xx}>0\) (red crosses). This is captured in the continuum model (solid lines) by the fact that the fitted active T1 rate \(\lambda \) in Eq. (9b) is positive (see App. D of the Supplementary Material for details of the fitting procedure), meaning that rearrangements are biased in the direction of the nematic tensor. Additionally, the fixed box imposes that the overall tissue is not sheared (\({\widetilde{V}}_{ij}=0\), blue crosses), and cells thus elongate in the direction orthogonal to that of the active T1 (green crosses). See also the lower left panel of Fig. 2 for a schematic explanation of the mechanism. Note that the tissue stress is along the y direction (\({\widetilde{\Sigma }}_{xx}<0\), grey curve). This is consistent with the fact that cells are elongated along the y direction, which induces an elastic stress along the same direction (\(\mu q_{xx}<0\) in Eq. (9c)). This elastic contribution adds up with the anisotropic one \(\zeta p_{xx}\) (with \(p_{xx}=-1\)), where \(\zeta \) is found to be positive from fits to the data, consistent with our definition of anisotropic bond tension which implies a larger stress along the elongation axis.

Fig. 2
figure 2

Relaxation dynamics after activation of anisotropic active stress under a fixed box boundary condition with anisotropic bond tension (left) or anisotropic cell stress (right). Top row. Total tissue shear (blue) decomposed into contributions of cell elongation change (green) and shear by topological rearrangements (red). The tissue stress is shown in grey. Only xx-components of the tensors are shown, xy-components are zero. Crosses are data from the vertex model averaged over 100 realizations. Error bars are smaller than the marker size. Solid lines are obtained by fits of the continuum model. Bottom row. Schematics of the cell rearrangement and elongation explaining the observed dynamics

Full size image

The case of stress-free boundary conditions is also illuminating, see left panel of Fig. 3 and Movie 2. The larger bond tension along the vertical axis leads to shorter vertical bonds, and to an active triggering of T1 transitions that close these short bonds and open horizontal bonds with lower tension (as sketched in the lower left panel of Fig. 3). These active rearrangements (red crosses) drive the shearing of the tissue along the x axis (blue crosses). In addition, cells are elongated along the x axis. Indeed, since \(\zeta \) is positive, and since \(\tilde{\sigma }_{xx}=0\) for stress-free boundary conditions, we deduce from Eq. (9c) that \(q_{xx}\) is positive. Note finally that a constant shear rate is obtained in the absence of external driving, which illustrates the active nature of the anisotropic bond tension.

Fig. 3
figure 3

Dynamics of tissue shear in a network with anisotropic bond tension (left) or anisotropic cell stress (right) under a stress-free boundary condition. Top row. Total tissue shear (blue) decomposed into contributions of cell elongation change (green) and shear by topological rearrangements (red). The tissue stress is shown in grey (and vanishes as stress-free boundary conditions are used here). Only xx-components of the tensors are shown, xy-components are zero. Crosses are data from the vertex model averaged over 100 realizations. Error bars are smaller than the marker size. Solid lines are obtained by fits of the continuum model. Bottom row. Schematics of the cell rearrangement and elongation explaining the observed dynamics

Full size image

3.2 T1 transitions driven by anisotropic cell stress

Interestingly, implementing anisotropy via an anisotropic cell stress as defined in Eq. (5) gives rise to a completely different behavior of the cellular network. As in the case of the anisotropic bond tension presented above, we consider a nematic tensor \(\varvec{{\mathcal {P}}}\) aligned with the vertical axis, such that cells elongated along the y axis experience a higher stress.

We first consider a fixed box boundary condition, see right panel of Fig. 2 and Movie 3. In this case, cells elongate in the direction orthogonal to the nematic axis since the anisotropic cell stress is higher along its direction (green crosses). This means, in the continuum model description, that the anisotropic stress magnitude \(\zeta \) is positive, as in the case of anisotropic bond tension. As a consequence, note that the tissue stress (grey crosses) changes sign during the simulation. At the beginning, tissue stress (9c) is dominated by the anisotropic cell stress \(\zeta p_{xx}<0\), which is higher in the y direction, leading to a \({\tilde{\sigma }}_{xx}<0\). As cells elongate in the x direction in response to this stress, \(q_{xx}\) grows and the elastic stress caused by this elongation starts overtaking the anisotropic one and the overall tissue stress changes sign. In a fixed box condition, cell elongation has to be compensated by T1 transitions in the opposite direction (\(R_{xx}<0\), red crosses), see also bottom right panel of Fig. 2 for a schematic explanation. Crucially, the orientation of these T1 transitions (that is, the direction in which new bonds are opened) is orthogonal to the orientation described in the previous section for the anisotropic bond tension (see also Table 1). This fact is reflected by the rate \(\lambda \) of anisotropic rearrangements in the continuum model, which is now found to be positive for anisotropic cell stress, whereas it was negative for anisotropic bond tension.

The consequences of an anisotropic cell stress can also be observed in the case of stress-free boundary conditions. Movie 4 shows an example of vertex model simulation in this case, and a quantification in terms of cumulative shear decomposition is displayed on the right panel of Fig. 3. Anisotropic cell stress drives cells to elongate in the direction orthogonal to the y axis, and we therefore have \(Q_{xx}>0\) (green crosses). As a consequence, cells have shorter bonds along their axis of elongation (the x axis) and T1 transitions close preferentially bonds along this axis and open new ones along the y axis (hence \(R_{xx}<0\), red crosses), see sketch in the lower right panel of Fig. 3. The tissue is therefore sheared along the vertical direction (blue crosses), which is opposite to the anisotropic bond tension case. Note also the change of sign of the tissue shear \(\widetilde{V}_{xx}\) at short times. At the beginning of the simulation, the anisotropic cell stress immediately drives the elongation of cells, implying \(\widetilde{V}_{xx}\simeq \mathrm{D}Q_{xx}/\mathrm{D}t>0\) at short time. With a delay, T1 transitions respond to this elongation and start contributing to the total tissue shear. They eventually dominate (at \(t\gtrsim 3\)), and account for the steady-state shear flow.

Table 1 Summary of the steady-state relative orientations. The tensor \(\varvec{{\mathcal {P}}}\) gives the direction of the tissue polarity. The tensor \(\varvec{R}\) indicates the direction of topological transitions (along which new bonds are opened), \(\dot{\varvec{Q}}=\mathrm{D}\varvec{Q}/\mathrm{D}t\) is the tensor for the rate of change of cell elongation and \(\widetilde{\varvec{V}}\) is the tissue shear rate, which indicates the direction along which the tissue elongates
Full size table

The different anisotropic stress realizations and boundary conditions lead to different relative steady-state orientations of cell elongation, T1 transitions and tissue shear with respect to the polarity axis. We summarize them in Table 1.

4 Energy balance in a tissue

In order to define the work performed by T1 transitions, we now discuss the energy balance of a cellular networks subject to active processes and external stresses. For simplicity and since this work is focused on shear, we limit our discussion to shape changes and shear deformations but do not include changes of tissue size. In the presence of an external shear stress \(\tilde{\varvec{\sigma }}^{\mathrm{ext}}\) applied to a tissue, the mechanical work per unit time \(\dot{w}_\mathrm{mech}\) performed on the tissue reads:

$$\begin{aligned} \dot{w}_{\mathrm{mech}} = \tilde{\varvec{\sigma }}^{\mathrm{ext}} : \tilde{\varvec{v}} \, . \end{aligned}$$
(11)

At force balance and for a homogeneous tissue, we have \(\tilde{\varvec{\sigma }}^{\mathrm{ext}}=\tilde{\varvec{\sigma }}\). Using the shear decomposition (9a) and the constitutive equations (9b)-(9c), the balance of elastic energy \(e=\mu \varvec{q}:\varvec{q}/2\) reads

$$\begin{aligned} \dot{e}&= \dot{q}_{\mathrm{heat}} + \dot{w}_{\mathrm{mech}} + \dot{w}_{\mathrm{chem}} \, . \end{aligned}$$
(12)

Here, we have defined the power \(\dot{q}_{\mathrm{heat}}\) supplied to the system in the form of heat, and the rate of chemical work by the environment on the system \(\dot{w}_{\mathrm{chem}}\). These quantities are given by

$$\begin{aligned} \dot{e}&= \mu \varvec{q} : \frac{\mathrm{D}\varvec{q}}{\mathrm{D}t} \, , \quad \dot{q}_{\mathrm{heat}} = -\eta \varvec{r}:\varvec{r} \, ,\nonumber \\ \dot{w}_{\mathrm{chem}}&= \eta {\tilde{\lambda }} \varvec{p}:\varvec{r} - \zeta \varvec{p}:\frac{\mathrm{D}\varvec{q}}{\mathrm{D}t} \, , \end{aligned}$$
(13)

where \(\eta =\mu \tau \) is the effective tissue viscosity, and we have defined the rate of active T1 transitions \({\tilde{\lambda }}=\lambda -\zeta /(\mu \tau )\). Note that the rate of heat production \(\dot{q}_{\mathrm{heat}}\) is always negative, indicating that the system releases heat to its surrounding.

The chemical power \(\dot{w}_{\mathrm{chem}}\) is an active contribution that would be vanishing for passive materials. There are two contributions to the chemical power stemming from different processes:

$$\begin{aligned} \dot{w}_{\mathrm{chem}} = \dot{w}_{\mathrm{T1}} + \dot{w}_{\mathrm{cell}} \, , \end{aligned}$$
(14)

where \(\dot{w}_{\mathrm{T1}}=\mu \tau {\tilde{\lambda }} \varvec{p}:\varvec{r}\) is the rate of work by T1 transitions, and \(\dot{w}_{\mathrm{cell}}=- \zeta \varvec{p}:\mathrm{D}\varvec{q}/\mathrm{D}t\) is the rate of work by cell deformations. Importantly, these contributions can be either positive or negative. A positive sign indicates that the process is performing work on the tissue, while a negative sign indicates that the process is typically dissipative, but it could also generate chemical free energy. Based on these considerations, we thus define active T1 transitions as T1 transitions for which \(\dot{w}_{\mathrm{T1}}>0\).

Determining the effective parameters from vertex model simulations reveal that the rate of active T1 transitions \({\tilde{\lambda }}\) is negative for the anisotropic bond tension but positive for anisotropic cell stress (see Table S2 of the Supplementary Material). However, the chemical work performed by T1 transitions \(\dot{w}_{\mathrm{T1}}\) is positive in both cases. This is because \(\varvec{p}:\varvec{r}\) is negative for anisotropic bond tension, while it is positive for anisotropic cell stress. We thus conclude that T1 transitions are active and perform chemical work on the tissue in both realizations of anisotropic active stress. The sign of \(\dot{w}_{\mathrm{T1}}\) could become negative if external stress induces shear along an axis perpendicular to the axis of spontaneous shear, i.e., by inducing shear (and rearrangements) along the y axis for the anisotropic bond tension case, or along the x axis in the anisotropic cell stress case. In this case, T1 transitions would be passive and the chemical energy of the active process would be dissipated.

Interestingly, the situation is slightly different for the chemical work performed by cells. Indeed, analysis of the vertex model simulations shows that \(\zeta \) is positive both for anisotropic bond tension and anisotropic cell stress (see Table S2 of the Supplementary Material). For fixed box boundary conditions, \(\varvec{p}:\mathrm{D}\varvec{q}/\mathrm{D}t\) is positive for anisotropic bond tension but it is negative for anisotropic cell stress (see Fig. 2). This reveals that the work performed by cell deformations \(\dot{w}_\mathrm{cell}\) is positive and active for anisotropic cell stress, while it is negative and passive for the anisotropic bond tension. In the case of anisotropic bond tension which drive active T1 transitions, cells elongate for fixed box boundary conditions along the y axis, thus increasing the length of bonds with large contractility, corresponding to \(\dot{w}_{\mathrm{cell}}\) being negative and typically dissipative. In contrast, in the case of the anisotropic cell stress, both T1 transitions and cell deformations are active and perform work. For stress-free boundary conditions, we find that for both anisotropic bond tension and anisotropic cell stress, \(\varvec{p}:\mathrm{D}\varvec{q}/\mathrm{D}t<0\) and therefore the work of cell deformations \(\dot{w}_{\mathrm{cell}}\) is always positive. Therefore, both T1 transitions and cell deformations are performing work on the tissue to shear it.

5 Discussion and conclusion

Using a vertex model with a preferred axis set by a prescribed nematic field, we have proposed two realizations of anisotropic active processes in tissues. The first one considers anisotropic bond tensions, for which cell bonds aligned with the nematic axis have higher contractility than those oriented perpendicularly. The second one involves an anisotropic cell stress aligned with the nematic axis. Importantly, although in both cases an active anisotropic stress exists that is contractile along the nematic axis, these two systems exhibit different orientations of T1 transitions and cell elongation (see Figs. 2 and 3).

In the case of anisotropic bond tension, cell bonds shorten and trigger T1 transitions. For fixed box boundary conditions, we therefore observe that cells elongate along an axis parallel to the nematic axis and orthogonal to the orientation of bonds opened by T1 transitions. For stress-free boundary conditions, both cell elongation and T1 transitions are perpendicular to the nematic axis. Anisotropic cell bond tension captures the behavior observed during the early stages of pupal wing development in Drosophila [8]. In that case, the increased contractility is oriented along the proximal-distal axis. The phenomenological parameters were measured as \(\zeta /\mu \simeq 0.33\), \(\tau \simeq 1.7\) h, and \(\lambda \simeq -0.11\) h\(^{-1}\), which corresponds to \({\tilde{\lambda }}\simeq -0.31\) h\(^{-1}\). This suggests that T1 transitions are active in these early stages and driven by anisotropic bond tension. Overall tissue shear was smaller than the cell shape change, corresponding to a case where the boundaries are slowly moving, not too far from the fixed box boundary conditions.

A different situation arises in the case of anisotropic cell stress, in which cells elongate and trigger T1 transitions. For fixed box boundary conditions, cells elongate perpendicular to the nematic axis, and T1 transitions are oriented parallel to the nematic axis. These orientations remain the same for stress-free boundary conditions (see Figs. 2 and 3). Anisotropic cell stress could contribute to the behavior observed during germ-band extension in the Drosophila embryo [22, 26]. The tissue extends along the anterior-posterior (AP) axis, which suggests that anisotropic active stress is orthogonal to this axis. Both in wild type and when tissue extension is obstructed by laser cauterization, cells elongate perpendicular to the AP axis. These observations are consistent with anisotropic cell stress both for fixed box boundary conditions (cauterization) and stress-free boundary conditions (rough approximation for wild type). In addition, it was reported that tissue elongation was driven by anisotropic medial myosin II pulses [26], which are expected to generate anisotropic cell stress.

We have considered these two realizations of anisotropic active processes separately. However, in biological tissues, both types of active stresses could coexist. This is likely the case during Drosophila germ-band extension, where an anisotropic accumulation of myosin II at cell junctions is observed [23,24,25] as well as anisotropic pulses of medial actin [26], indicating the existence of anisotropic cell stress. This suggests that both types of active anisotropic processes could be simultaneously relevant.

If both processes are at work at the same time with the polarity \(\varvec{{\mathcal {P}}}\) aligned to the same axis, they would be antagonistic and oppose each other. This case is similar to a tug-of-war situation where the resulting T1 transitions would be the net result of the two opposing anisotropic processes generating shear along orthogonal axes. One can speculate that the relative strength of these two opposing processes could be fine-tuned such that the resulting net rate of active T1 transitions would be vanishing, even though the system would still be chemically active and T1 transitions fluctuate strongly forward and backward. This could lead to a fluidization of the tissue or give rise to stress oscillations. Furthermore, at this balance point, a biological tissue could be capable of changing rapidly to one of the two steady states with orthogonal shear axis if the balance between the opposing activities is lost. Recent work has proposed more detailed descriptions of T1 transitions, including delay [53] or nonlinear dynamics [54]. While we expect the qualitative picture studied in this manuscript to be rather robust, it will be interesting to study the role of such nonlinear or delayed T1 transitions in the anisotropic active processes discussed here.

Using a linear continuum description that captures the anisotropic dynamics of the vertex model, we have shown that the difference between these two realizations of active stress is captured by a relative sign difference between the active stress magnitude \(\zeta \) (positive in both scenarios) and the active T1 rate \(\lambda \) (positive for anisotropic cell stress, negative for anisotropic bond tension). Despite these differences, an analysis of the energy balance in the system reveal that T1 transitions perform chemical work on the tissue in both cases, and can therefore be referred to as active. The determination of \(\zeta \) and \(\lambda \) experimentally is a challenge. However, the determination of the rate of active T1 transition \({\tilde{\lambda }}=\lambda -\zeta /\mu \tau \) may be accessible by state-of-the-art experimental techniques. Indeed, at steady state and for fixed box boundary conditions, the tissue shear stress reads \(\varvec{\tilde{\sigma }}=-\mu \tau {\tilde{\lambda }}\varvec{p}\) and could therefore provide a readout for the sign of this activity coefficient. This local tissue stress could for instance be measured by injecting liquid oil droplets, as recently shown in the zebrafish embryo [12, 14].