Modelling Cell Orientation Under Stretch: The Effect of Substrate Elasticity

When cells are seeded on a cyclically deformed substrate like silicon, they tend to reorient their major axis in two ways: either perpendicular to the main stretching direction, or forming an oblique angle with it. However, when the substrate is very soft such as a collagen gel, the oblique orientation is no longer observed, and the cells align either along the stretching direction, or perpendicularly to it. To explain this switch, we propose a simplified model of the cell, consisting of two elastic elements representing the stress fiber/focal adhesion complexes in the main and transverse directions. These elements are connected by a torsional spring that mimics the effect of crosslinking molecules among the stress fibers, which resist shear forces. Our model, consistent with experimental observations, predicts that there is a switch in the asymptotic behaviour of the orientation of the cell determined by the stiffness of the substratum, related to a change from a supercritical bifurcation scenario, whereby the oblique configuration is stable for a sufficiently large stiffness, to a subcritical bifurcation scenario at a lower stiffness. Furthermore, we investigate the effect of cell elongation and find that the region of the parameter space leading to an oblique orientation decreases as the cell becomes more elongated. This implies that elongated cells, such as fibroblasts and smooth muscle cells, are more likely to maintain an oblique orientation with respect to the main stretching direction. Conversely, rounder cells, such as those of epithelial or endothelial origin, are more likely to switch to a perpendicular or parallel orientation on soft substrates.


Introduction and Biological Background
In the 1980s the study of cardiovascular diseases gave rise to the need of understanding the behaviour of cells undergoing periodic stretch. In fact, heart beats induce a periodic stretch not only of its own cells, but also of the arterial walls due to the induced periodic pressure and then inflation. It was then observed that smooth muscle cells in the intima of the aorta are oriented along specific directions (Buck 1979;White et al. 1983). Specifically, they either align along the vessel axis (Zhao et al. 1995) or obliquely forming helical-like structures with an angle of 20 • -40 • with respect to the vascular axial direction (Driessen et al. 2004;Gasser et al. 2006;Osborne-Pellegrin 1978). Starting from this need, Buck (1980) examined the behaviour of fibroblasts seeded in vitro on a cyclically stretched rubber plastic substrate and observed that they tended to reorient more perpendicularly with respect to the stretching direction than along it. Similar experiments were then performed by many other researchers with many cells types (see Giverso et al. 2023, for a review).
In general, it is observed that in most cell lines, such as fibroblasts, myofibroblasts, cardiomyocytes, and endothelial cells, periodic stretching triggers a dynamic response of the cytoskeleton that results in the formation of a preferred orientation of both the cytoskeleton and the cell shape itself, at either an oblique or a perpendicular angle with the main stretching direction. Actually, the aspect ratio of the cell can even reach values close to 10 with the nucleus that achieves an ellipsoidal shape with the long axis along the same direction.
In this process the leading role is played by the cytoskeleton, that is characterised by the presence of strong filamentary structures, named as stress fibres (SFs), which are bundles of aligned actin filaments cross-linked by several molecules, such as fascin, fibrin, and actinin. Then, cell focal adhesions (FAs) cluster at the ends of assembled SFs, anchoring the cell to the substrate.
Generally speaking, the action of mechanical stress induces an internal remodelling process involving the disruption, formation, and reorganization of both SFs and FAs. In this way the cell is able to adapt to the external mechanical cue, fostering changes in shape and orientation of the whole cell. Specifically, the remodelling process slowly continues untill the cell reaches a stable configuration characterised by a clearly visible mean angle between its major axis (i.e., the one defined by the main assembly of SFs) and the direction of largest stretching. It is useful to remark that cell internal orientation can occur even in the absence of effective spatial migration, i.e., it does not necessarily lead to a head-and-tail differentiated cell configuration. For this reason, we will use the word orientation, rather than polarization.
After an initial attempt to model such behaviour as a strain avoidance mechanism, based on the minimization of the sensed strain (De et al. 2008;Faust et al. 2011;Morioka et al. 2011;Safran and De 2009;Wang 2000;Wang et al. 1995), Livne et al. (2014) showed that a better fit was obtained by minimizing the elastic energy. This approach, firstly developed in a linear elasticity framework, was recently generalised to large deformations in Ciambella et al. (2022), Lucci and Preziosi (2021) using a very general nonlinear orthotropic elastic model.
Then, in Giverso et al. (2021) a linear viscoelastic model was proposed to explain the experimental evidence that cells reorient only if the frequency of the deformation applied to the substrate is sufficiently high, say, greater than 0.1 Hz for fibroblasts (Hsu et al. 2009;Jungbauer et al. 2008;Liu et al. 2008). In fact, this experimental fact requires the introduction of viscoelastic effects with a characteristic response time that need to be compared with the period of oscillations.
A different approach was used by Xu et al. (2016), as well as in Mao et al. (2021) and Xu et al. (2018). They schematize cells as formed by two parallel elastic struts along the orientation axis and two other struts perpendicular to them. A string connecting the ends of the struts is also considered, but then neglected. The energy of the system is then evaluated, but unfortunately it presents some flaws in the computation of the work done by the applied force. In addition, there is some confusion with some signs (see, for instance, Eqs. (13) and (18) in Xu et al. 2016). Fortunately, apart from some signs, that are a posteriori adjusted, the energy discussed is consistent with a special case of the linear elastic energy. This allows them to recover the same ordinary differential equation for the temporal evolution of cell orientation obtained by Livne et al. (2014), and thus the same equilibrium states for time independent stretches.
In the present article we use a similar approach, schematizing the cell as formed by two elastic elements that respectively model stress fiber assemblies in its main and transversal orientations (Sect. 2). The two elements are hinged at their center through a torsional elastic spring. Specifically, the elastic elements are intended to model the behaviour of both stress fibers and focal adhesions, while the torsional spring mimics the action of crosslinking molecules connecting stress fibers and resisting to shear. The mathematical model proposed is then based on the assumption that cell reorientation results from processes occurring at different time scales, i.e., the cell readily deforms to comply to the deformation of the substrate, and actively reorients slowly remodelling its internal cytoskeletal structure to minimise its internal energy. In Sect. 3 the analytical solution of the proposed model is given and its qualitative behaviour for long times is discussed. It is shown that, in the limit of very stiff substrates, our mathematical model reduces to the models mentioned above, e.g. Livne et al. (2014), Lucci and Preziosi (2021), being then able to reproduce the related experimental results. However, the main focus of this article is to discuss the influence of the stiffness of the substrate on cell orientation. In this respect, the first fact that is observed is that the substrate need to be stiff enough to induce cell reorientation. We show that, when the stiffness of the substrate becomes much smaller than that of the stress fibers, the time to reach the equilibrium configuration increases considerably, so that reorientation is not observed in the time of the experiment, being moreover hidden by random behaviours as explained in Loy and Preziosi (2023).
The model is also able to explain the empirical observation that when the substrate is very soft, cells elongate along the main stretching direction or perpendicularly to it, rather than achieving an oblique orientation, as shown in the experiments reported in Terracio et al. (1988), Thodeti et al. (2009). In Sect. 4 we discuss in detail how the model is able to explain this switch in the asymptotic trend and we relate it to a change of bifurcation scenario from a supercritical one, with a stable oblique equilibrium, to a subcritical one, in which the oblique configuration becomes unstable.
Finally, in Sect. 5 we investigate the impact of cell elongation demonstrating that an increase in cell elongation leads to a larger region in the parameter space leading to either oblique or perpendicular orientations, at the expenses of the oblique orientation.
This indicates that the stability switch is more likely to occur in rounder cells, such as epithelial or endothelial cells, compared to elongated cells, such as fibroblasts and smooth muscle cells.

Mathematical Model
In order to model the cell re-orientation process described in the Introduction, it is important to observe that the phenomenon is characterised by different time scales, which allow to identify the following phases: 1. As soon as the substrate is deformed, the first reaction of the cell is to immediately adapt to the deformation of the substrate. This results in a stretching of the cell and occurs at the same time scale as the period of the imposed external deformation, which is typically of the order of one second or less. 2. Under the action of stretching, the cell then re-adjusts its cytoskeletal orientation in order to minimize the stored internal energy. This process occurs on a much longer time scale, which is in the range of several minutes. As a result of this internal remodelling process, the cell may acquire an orientation that is different from the initial one.
With this in mind, let us define a bidimensional reference frame (O, i, j) with the origin at the center of the mass of our representative cell and the axes aligned with the main stretching directions. If the reaction of the cell to the externally applied deformation were negligible, the deformation gradient in the plane that is applied to the substrate would be homogeneous and given by where r is a fixed constant called biaxiality stain ratio. In Eq. (1), the strain amplitude ε(t) is assumed to be small and continuous in time. Actually, it was experimentally observed that strains larger than 40% can cause cell death Boccafoschi et al. (2007) and usually a strain of 10% is applied Giverso et al. (2023). The range of r usually tested in experiments is r ∈ [0, 1], where r = 0 corresponds to a substrate clamped along the sides parallel to the main stretching direction, while r = 1 implements an isochoric planar deformation (i.e., ε yy = −ε x x ). Very few results are reported for negative values of r , which corresponds to stretching in both directions (see, for instance, Li et al. 2016). However, the fact that the x-axis is assumed parallel to the maximum stretching direction implies that r ≥ −1. At the other extremum, r > 1 would correspond to an untested strong compression of the substrate along the y-axis.
For this reason, we hereafter focus on the range r ∈ [−1, 1]. As shown in Fig. 1, the cell cytoskeleton is schematised by two elastic elements, named a and b, with elastic moduli k a and k b and end points A, A and B, B , respectively. They are in charge of describing the response of the cell to the external mechanical input: the element a along the main orientation of SFs, hereafter characterised through the time-dependent angle , and the element b along a transversal (i.e., rotated of an angle θ ) direction. These two elements are further assumed to be: Fig. 1 Schematization of the cell as two elastic elements that model stress fibers-focal adhesions assemblies in its main and transversal orientations. The main orientation of the cell is determined by the angle (a) that slowly evolves due to the remodelling of the actin cytoskeleton. Due to the cell reaction to stretching, during the deformation the anchoring points X A and X B go to x A = x S A and x B = x S B (b), where, for instance, x S A is the location of the point that would be occupied by X A as a result of the imposed deformation in absence of cell reaction (Color figure online) • linked to the substrate at their end points, which indeed model focal adhesion complexes. Specifically, we denote as S Q the point of the substrate that is connected with the end point Q of one of the two cytoskeletal elements of the cell (being Q ∈ {A, A , B, B }); • hinged at their centers (which are superimposed to the domain origin O); • interconnected by a torsional spring, with modulus equal to k θ , that is at rest when the two elements are orthogonal, i.e., when θ = π 2 ; • the two cytoskeletal elements have relaxed lengths equal to 2L a and 2L b , respectively, which can be interpreted as mean dimensions of the cell in the absence of stresses. They are assumed to be constant.
Let us now model the above-introduced first phase of the process, i.e., the cell instantaneous response to the imposed substrate deformation. To do this, we hereafter take advantage of the symmetry of the system and focus on the upper half of the domain. Before stretching, the cell is in a rest configuration. This, as seen, implies orthogonal elastic elements with semi-lengths L a and L b , respectively, with main SF orientation defined by the angle . The two end points of interest occupy domain locations X A and X B that coincide with the locations of the substrate points to which they are linked to, i.e., X A = X S A and X B = X S B (see Fig. 1a).
As a consequence of the imposed deformation, in absence of cell reaction, the two substrate points of our interest would go to x S A = FX S A and x S B = FX S B (see Fig. 1b). However, the cell is not passively dragged by such a substrate displacement, but it rather resists to the underlying stretch. So, the end points A and B shift to x A = x S A and x B = x S B . In this respect, the points x A and x B are identified by the angles ( − θ a ) and ( + π/2 + θ b ) and by the distances from the origin L a = L a + ξ a and L b = L b + ξ b , respectively. In particular, θ a and θ b give the variation of the relative orientation of the two cytoskeletal elastic elements, whereas ξ a and ξ b give the variation of their lengths. The effective values of these quantities can be obtained by the minimization of the total energy stored by the system as a consequence of the imposed substrate deformation where the last term is an elastic contribution that accounts for the above-explained mechanical resistance of the cell to the dragging due to the substrate deformation, being k s the corresponding elastic modulus. As already mentioned, the Hooke constants of the elements along a and b are respectively k a and k b , while k θ is the torsion coefficient of the torsional spring in O. In the limit of small deformations, i.e., for (see the "Appendix" for more details on the technical calculations), which is minimised by a cell configuration characterised by So, for a given , under the externally imposed stretch, we assume that the cell readily achieves a cytoskeletal arrangement defined, in terms of SF lengths and orientations, by the equilibrium quantities in Eq. (4). In this respect, being an equilibrium, the viscosity of the cytoskeleton does not play any role. The influence of the cell membrane is also neglected. However, if it is schematised as an elastic spring connecting the end points of the two elements it would lead to a contribution that can be absorbed in the coefficients of Eq. (4).
Let us now focus on the slower remodeling process. In order to do that, we assume that the temporal evolution of is driven by the virtual work of the forces exerted on the cell cytoskeleton due to the externally imposed deformation of the substrate. In the elastic case such a virtual work can be written in terms of the derivative of the internal potential energy In addition, there is a dissipation force that resists to remodelling and eventually results in the application of a viscous-like moment proportional to the angular velocity. So, neglecting inertia, the evolution of cytoskeletal reorientation can be described by the following equation where h is related to the characteristic time for stress fiber and focal adhesion renewal dynamics.
Recalling Eq. (4), we have that and Eq.(5) therefore rewrites as We firstly notice that the case r = −1, corresponding to the equi-biaxial deformation with ε yy = ε x x , will lead to no evolution of . In fact, from the mathematical viewpoint the r.h.s. of Eq. (7) vanishes and from the mechanical viewpoint it corresponds to an isotropic deformation with equal strain in all directions.
We also observe that in the limit k s k a , k b , (7) is identical to the evolution equation considered in Livne et al. (2014) for the case k θ = 0, and it can be also obtained from the model proposed in Giverso et al. (2021) by neglecting viscoelastic effects and second order mixed energy terms in the quadratic form of the elastic energy.
In addition, softening the substrate leads to a decrease of the forcing term in Eq. (7). Therefore, the dynamics slow down with the characteristic time that increases like k −2 s as k s decreases. Lastly, we remark that we are only focusing on orientational changes and not on shape change, that is, the rest length of the axes of the cells do not remodel.
For the following discussion, it is useful to scale the mechanical coefficients in Eq. (7) with k a and introduce the dimensionless quantitieŝ and the aspect ratio λ = L 2 b /L 2 a . Then, by scaling times with τ = h k a L 2 a , we can rewrite Eq. (7) as where Interestingly, α is independent from the biaxiality ratio r , while the mechanical parameters, i.e.,k bkθ ,k s , and the aspect ratio λ influence both α and β.

Analytic Solution and Asymptotic Behaviour
The analytical solution of Eq. (9) can be written in an implicit form as where C is the integration constant depending on the initial value (0) and Though it is not immediate to visualize the behaviour of the solution, it can be noticed that, for any ε(t), for large times the r.h.s. of Eq. (11) will eventually vanish. So, the Fig. 2 a Summary of the asymptotic trends in the (α, β) plane. In the panels on the right, the bifurcation diagrams for a time indepedent strain are reported for b α > 0 and c α < 0. In the former case the oblique orientation is stable when it exists and bifurcations are supercritical. In the latter it is unstable and bifurcations are subcritical. In particular, only the parallel and/or perpendicular orientations are stable (Color figure online) asymptotic behaviour of (t) is governed by the sign of the exponents in the factors of the l.h.s. of Eq. (11), i.e. it is substantially dictated by the sign of α, β, and γ . In fact, by recalling that |ε x x | > |ε yy | implies r > −1 and defining eq obl such that one has the following Proposition Proposition For r > −1 • If α > β > 0, then (t) → eq obl ; • If α < β and β > 0, then (t) → eq = 0; • If α < β < 0, then either (t) → eq = 0, or (t) → eq ⊥ = π 2 ; • If α > β and β < 0, then (t) → eq ⊥ = π 2 . Proof Let us start observing that (t) = π 2 can nullify the l.h.s. of Eq. (11) only if β < 0 (pale blue and magenta regions in Fig. 2a).
Regarding the oblique orientation eq obl , Eq. (13) implies that it exists only if 0 ≤ β α ≤ 1 (i.e., the pale blue and yellow regions in Fig. 2a). So, the fact that the ratio β/α must be positive implies that the l.h.s. of Eq. (11) can be nullified when = eq obl only for γ > 0, i.e. α > β. This rules out the possibility of negative values of α and β, and leaves only the possibility α > β > 0, that is the yellow region in Fig. 2a.
We observe that Eq. (13) is identical to the oblique orientation found in Giverso et al. (2021), Livne et al. (2014), Lucci and Preziosi (2021), where, however, the effect of substratum deformability was not taken into account. We also notice that the asymptotic values just identified remain the same if one looks at the limit case of a time independent strain. This allows to discuss the possible scenarios using classical stability tools and terminologies, as done in previous papers (see, for instance, Giverso et al. 2021;Lucci and Preziosi 2021).
Before proceeding with it, we highlight that, taking separately the three terms on the r.h.s. of Eq. (9), for a time independent stretch, the first two terms tend to destabilize the equilibrium orientations eq and eq ⊥ , while the last one tends to stabilize them. The opposite holds true for eq obl . As we will see in Sect. 4, this will lead to important switches in the asymptotic behaviour also seen in some experiments (see, for instance, Terracio et al. 1988;Thodeti et al. 2009). Specifically, referring to Fig. 2b, we can observe that for α > 0 • eq is stable for β ≥ α; • eq ⊥ is stable for β ≤ 0; • eq obl is stable whenever it exists, i.e. when 0 ≤ β ≤ α. On the other hand, referring to Fig. 2c, for α < 0 • eq is stable for β ≥ α; • eq ⊥ is stable for β < 0; • eq obl is always unstable. Therefore, in the case of a time independent strain the system is characterised by a switch from a supercritical to a subcritical bifurcation scenario when α changes sign from positive to negative. In the following, we will thus extend this classification, calling supercritical (resp. subcritical) scenario the one with α > 0 (resp. α < 0).

Effect of Substratum Stiffness
As α and β both depend onk s , the effect of the elasticity of the substratum is not immediately clear, because it influences both dimensionless parameters at the same time. On the other hand, the aim of the article is to highlight the effect of the substratum on the evolution of cell orientation. Therefore, hereafter we will translate the discussion done in the previous section in terms of substrate stiffness to prove that it can determine a switch in the asymptotic behaviour. In particular, we want to show that for decreasing values of the substrate stiffness, the bifurcation landscape can change from a supercritical one as in Giverso et al. (2021), Lucci and Preziosi (2021) to a subcritical one with an unstable oblique orientation and stable parallel and/or perpendicular orientations, as observed in Terracio et al. (1988), Thodeti et al. (2009).
Before doing that we make some comments on the ammissible physiological range of the parameters, so that we can limit the discussion to the biologically relevant dynamics. First of all, from phenomenological observations, we can certainly restrict our attention to λ ∈ (0, 1] (i.e., 0 < L b ≤ L a ) because experiments show that cells are more elongated along the main direction, i.e., along the element a rather than along b. In addition, from the mechanical point of view, it can also be observed that since the response of cells to deformations along the main orientation is stronger Fig. 3 Surfaces γ (k b ,k θ ,k s , r ) = 0 in blue and β(k b ,k θ ,k s , r ) = 0 in brown that delimit the different scenarios in the case of round cells, i.e., with λ = 1, for different values of r , positive in the top row, negative below than along the perpendicular one (even much stronger), we can takek b ∈ (0, 1) (i.e., 0 < k b < k a ). Similarly, we expectk θ ∈ 0, 1 2 . This can be realised comparing, for the same deformation involving the opening of the angle in O for a round cell, the energetic terms related to the torsional spring (proportional to k θ ) with the energy that would be stored by a spring with rigidity k a connecting the extrema of the elements a and b.
So, having in mind the above biologically sound constraints, we will focus the following discussion on the subspace of parameters defined by though it could be easily extended outside this range. In Fig. 3and then in Fig. 6we plot, respectively for round cells λ = 1 and elongated cells with λ = 0.5, the surfaces corresponding to γ (k b ,k θ ,k s , r ) = 0 in blue and to β(k b ,k θ ,k s , r ) = 0 in brown for several fixed values of r . These surfaces subdivide the space of parameters in regions characterised by the different asymptotic behaviours, that are identifiable by the same colours as in Fig. 2a.
Focusing on the coordinate planes, we have that in the limitk θ = 0, if r > 0 then β and γ are both positive, and so, the solution can only tend to eq obl . If, instead, r < 0, that corresponds to imposing a tension in both directions, then the solution will tend 79 Page 12 of 20 A. Colombi et al.
to eq ⊥ above this interval and to eq below it.
On the other hand, in the limitk b = 0, if r < 0, then β is always negative. So, eq ⊥ is always attractive and eq might be for softer values ofk s below the blue surface γ = 0, i.e, fork s < k low where = λ + 1 2λ , if k low s is positive. We explicitly remark that in this case a trend toward eq obl is ruled out. In addition, if r < −1 + 2 2 , then k low s becomes negative for k θ > 1 + r 2 . In this case also the trend toward eq is ruled out and the solution can only tend to eq ⊥ . Instead, in order to describe what happens in the limitk b = 0 when r > 0, corresponding to the usually tested cases, we start noticing that since γ = β + 1 − r 1 + r > β (we recall that r < 1), the brown surface β = 0 is always above the blue one γ = 0. Hence, the solution will tend toward eq obl above the brown surface, i.e., for and to eq ⊥ for values ofk s below k up s , with eq that becomes also attractive for suitable initial conditions when the condition (15) is satisfied. In order to describe what happens fork b ,k θ = 0, for decreasing values ofk s , we use as a representative example the case with λ = 1 and r = 0.5, corresponding to an isocoric deformation, and we fixk θ = 0.25 (see Fig. 4).
According to where the parameter values fall, one has a trend toward eq obl (yellow region), eq ⊥ (magenta region), eq (green region) or either eq ⊥ or eq (cyan region).
In particular, fixingk b = 0.5 one has that for stiffer substrates (namely,k s > 0.46) the solution will tend toward eq obl . For softer substrates with 0.41 <k s < 0.46 the solution will tend toward eq and for even softer substrates withk s < 0.41 the solution may also tend toward eq ⊥ according to the initial orientation. For very low values ofk b (e.g.,k b = 0.15) one has that for a very high stiffness, i.e.,k s > 1, again the solution will tend toward eq obl . For softer substrates (i.e., Fig. 6 Surfaces γ (k b ,k θ ,k s , r ) = 0 in blue and β(k b ,k θ ,k s , r ) = 0 in brown that delimit the different scenarios in the case of elongated cells with λ = 0.5, for different values of r , positive in the top row, negative below 0.05 <k s < 1) the solution will tend toward eq ⊥ , with the parallel orientation that adds up when 0.05 <k s < 0.58. Then, for very soft substrates (i.e.,k s < 0.05) the solution can only tend toward eq .
From the central row of Fig. 5, it can be observed that variations in the value ofk θ = 0.25 (that is the ratio between cell resistence to shear and the cell elastic modulus in its main direction) lead to similar trends, obviously with different values of the parameters leading to transitions. In particular, generally speaking, lower values ofk θ lead to larger regions characterised by a trend toward eq obl and smaller regions characterised by a trend toward eq and/or eq ⊥ . This means that resistence to shear favours the alignment of the representative cell in the parallel and perpendicular orientations. Instead, in the first row of Fig. 5 it can be for instance observed that the regions corresponding to admissible parallel asymptotic orientation (i.e., the cyan and green ones) get larger decreasingk b , that is the ratio between the elastic moduli of the two elements a and b. This suggests that the presence of fibers not all aligned in the same direction favours oblique asymptotic orientations. The fact that softer substrates, i.e., with lowerk s , are characterised by larger regions with parallel and perpendicular orientations (see the third row in Fig. 5) has been already discussed above.

Effect of Cell Elongation
In this section we will focus on the effect of the aspect ratio λ of the representative cell on its asymptotic behaviour. Generally speaking, referring to Fig. 6 it can be observed that the cyan region characterised by the trend either to eq or eq ⊥ (with eq obl unstable) becomes much smaller. In particular, the magenta region where the asymptotic trend of cell orientation is toward eq ⊥ enlarges considerably. The yellow region related to eq obl slightly changes w.r.t. to the case with λ = 1 when r is larger than 0.5. Conversely, for smaller values of r , the yellow region decreases at the benefit of magenta region related to eq ⊥ . This is even more evident for negative values of r . In fact, the region even disappears as r gets closer to −1. Further details on the evolution and on the regions are given in Figs. 7and 8.
This means that, on soft substrates, round-like cells more likely orient themselves in a parallel or perpendicular fashion with respect to the main stretching directions, than elongated cells, that prefer the perpendicular orientation also with respect to the oblique one. So, cells like fibroblasts and smooth muscle cells, that are typically characterised by a large aspect ratio, are more prone to keep an oblique or perpendicular orientation even on soft substrates with respect to endothelial or epithelial cells, that usually keep an almost round shape, even under stretch.

Conclusions
We have modelled the behaviour of cells seeded on an elastic substrate undergoing a periodic deformation using a discrete model for the cell interacting with the continuum substrate, focusing on the effect of the mechanical properties of the cells and above all of the substrate. Specifically, we schematised the cell as formed by two elastic elements, modelling its stress fiber/focal adhesion complexes along its main and transversal orientations, hinged at the center through a torsional spring.
The main focus of the present work was to describe the possible asymptotic scenarios as a function of the mechanical characteristics of the substrate. Consistently with the experiments, it is found that the supercritical scenario characterised by the presence of a stable oblique configuration is favoured by sufficiently large substrate stiffness, while the subcritical scenario characterised by the stability of parallel and perpendicular orientations is favoured by softer substrates. In addition, the region characterised by subcriticality (i) increases ask b (measuring the ratio of the cell response along the direction perpendicular to the main cell orientation axis versus the latter) decreases and (ii) increases ask θ (measuring the ratio of the cell response to shear versus the one along the main cell orientation axis) increases.
Finally, we studied the effect of cell elongation, finding that the region of the mechanical parameters characterised by a subcritical scenario decreases with cell elongation. This means that oblique orientations with respect to the main stretching direction are more easily mantained by elongated cells such as fibroblasts and smooth muscle cells, than by rounder cells such as epithelial or endothelial cells.
We also showed that the r.h.s. of the ordinary differential Eq. (9) governing cell reorientation presents a factork 2 s (1 +k s ) 2 that leads to a considerable increase of the reorganization time as the substrate softens. Specifically, this means that if the stiffness decreases by one order of magnitude (e.g., from about 1MPa for silicon elastomers and polydimetylsiloxane (PDMS), to about 100 KPa of collagen gels, or even less), we can expect the reorientation time to increase by two orders of magnitude (i.e., from a couple of hours to more than a week). This practically implies that over very soft substrates reorientation is not observed in the time of the experiment and is easily hidden by random behaviours, as explained in Loy and Preziosi (2023) As the discrete cell model is based on elastic-like elements, it is not able to describe the dependence of cell behaviour on the frequency of the imposed deformation, i.e., the experimental fact that, roughly speaking, oscillations need to be fast enough (say, higher than 0.1 Hz) to trigger cell reorientation. In order to do that, following what done in Giverso et al. (2021), one should include viscoelastic effects relative to the remodelling of focal ahesions and to the viscoelastic behaviour of the cytoskeleton. This might also take into account of the behaviour of slip bonds and of the development of catch bonds in response to the imposed traction. However, the modelling approach proposed here can be properly generalised to include such effects by adding viscous elements to the struts and to the torsional spring, possibly time dependent and stress dependent, so that they become Maxwell-like elements. Actually, with respect to Giverso et al. (2021), the modelling structure proposed here has the advantage to allow to discriminate between the reaction of the cytoskeleton from the one due to focal adhesions reorganization, here treated as a single structure. In a similar way, it would be interesting to include the active response of the acto-myosin machinery to the imposed deformation, once this effect is clarified in experiment.
Equations (17)- (18) are obtained recalling that, before the imposed external stress, Due to the elastic characteristics of the cell cytoskeleton, the end points of the elements instead go to x A = l a cos( − θ a ) i + l a sin( − θ a )j, In the limit of small deformations, Eqs. (21)-(22) simplify to x A ≈ [L a cos + ξ a cos + L a sin θ a ] i + [L a sin + ξ a sin − L a cos θ a ] j, being l a = L a + ξ a and l b = L b + ξ b . From Eqs.