Contractile network models for adherent cells
Abstract
Cells sense the geometry and stiffness of their adhesive environment by active contractility. For strong adhesion to flat substrates, two-dimensional contractile network models can be used to understand how force is distributed throughout the cell. Here we compare the shape and force distribution for different variants of such network models. In contrast to Hookean networks, cable networks reflect the asymmetric response of biopolymers to tension versus compression. For passive networks, contractility is modeled by a reduced resting length of the mechanical links. In actively contracting networks, a constant force couple is introduced into each link in order to model contraction by molecular motors. If combined with fixed adhesion sites, all network models lead to invaginated cell shapes, but only actively contracting cable networks lead to the circular arc morphology typical for strongly adhering cells. In this case, shape and force distribution are determined by local rather than global determinants and thus are suited to endow the cell with a robust sense of its environment. We also discuss non-linear and adaptive linker mechanics as well as the relation to tissue shape.
pacs:
87.16.Ln, 87.10.-e, 87.17.RtI Introduction
During the last decade, it has been increasingly realized that adherent cells actively explore and respond to the geometry and stiffness of their adhesive environment, with dramatic consequences for fundamental cellular processes such as survival, proliferation, differentiation and migration Janmey et al. (2009); Geiger et al. (2009). Using chemical inhibitors for both myosin motors and actin filaments, it has been shown that active contractility of the actin cytoskeleton is an indispensable requirement for the observed cellular response to the physical properties of the environment. Although more and more details are revealed regarding the molecular mechanisms underlying these sensing processes Zaidel-Bar et al. (2007); Zaidel-Bar and Geiger (2010); Hoffman et al. (2011), what is still missing is a systems-level understanding of how cells integrate the way force is generated, distributed, and sensed over the whole cell. Therefore theoretical models are required which describe how force is propagated through the actin cytoskeleton as a function of environmental geometry and stiffness. Because the actin cytoskeleton is an integral part of the sensing capabilities of cells, one expects that its mechanical properties have evolved to support these important cellular functions.
Although it is a standard procedure in experiments to block force generation and propagation with chemical inhibitors or by RNA interference, it is very difficult to measure how force is distributed inside cells and between cells and their environment. Different experimental approaches have been developed to meet this challenge. Traction force microscopy allows to measure the forces transmitted to the substrate Dembo and Wang (1999); Butler et al. (2002); Schwarz et al. (2002); Tan et al. (2003); Sabass et al. (2008). Recently this technique has been extended in such a way that also cell-cell forces can be estimated from cell-matrix forces Liu et al. (2010); Maruthamuthu et al. (2011). However, it is important to note that many forces balance inside the cell and are not transmitted to the substrate, so the forces existing inside cells might be much higher than appreciated from traction force microscopy Bischofs et al. (2009). Laser cutting allows to estimate forces from the mechanical relaxation after cutting load-carrying elements like microtubules in the mitotic spindle Grill et al. (2001, 2003) or stress fibers in the actin cytoskeleton Kumar et al. (2006); Tanner et al. (2010); Colombelli et al. (2009); Smith et al. (2010). Laser ablation can be used for subcellular analysis of cortical tension Mayer et al. (2010). However, these experiments only probe local relaxation events of prominent cytoskeletal structures and therefore might miss the global effects of more distributed and less visible structures. Micromanipulation can be used to distort the mechanical balance of the cell globally Riveline et al. (2001); Paul et al. (2008); Heil and Spatz (2010); Kirchenbüchler et al. (2010), but the resulting changes in force distribution can only be estimated indirectly from its effects, e.g. growth of focal adhesions. To achieve a more systematic understanding, these experimental approaches have to be complemented by theoretical approaches.
Because the mechanical properties of cells, extracellular matrix and tissue are strongly determined by filamentous networks of proteins like actin, tubulin, lamin, spectrin or collagen, mechanical networks are widely used theoretical models for cell and tissue mechanics Howard (2001); Boal (2002). One of the best studied cases is the red blood cell, whose shape and mechanics has been studied with network approaches in very large detail Hansen et al. (1996, 1997); Boey et al. (1998); Discher et al. (1998); Lim H. W. et al. (2002); Noguchi and Gompper (2005); Li et al. (2005, 2007); Fedosov et al. (2010). Modern computer power permits to simulate each of the roughly spectrin links separately and with molecular detail, for example using the appropriately parametrized force-extension curve of a semiflexible polymer Li et al. (2005). In the limit of small extensions, these models usually reduce to Hookean networks.
For adherent cells, the main structural determinant is the actin cytoskeleton, whose mechanics differs in several important aspects from the one of the spectrin network of red blood cells. In general, the molecular structure of the actin network is much less defined. Its most prominent feature in adhesion is strong contractility due to activity of myosin II motors. This observation implies that the mechanical links between the nodes of the network cannot be simple actin filaments, but have to be bundles of actin filaments crosslinked and tensed by myosin II motors. The simplest model for prestress in a mechanical network is the introduction of a finite resting length which is smaller than the typical extension of each link. Indeed, one- and two-dimensional spring networks with prestress are widely used for modeling cell migration DiMilla et al. (1991); Gracheva and Othmer (2004); Bottino et al. (2002); Dokukina and Gracheva (2010).
Network models are conceptually very appealing due to their multi-scale nature: by changing the microscopic rules for the mechanics of the links, one can explore how the macroscopic behavior of the whole network changes. In particular, important biological effects like viscoelasticity of the links or coupling to diffusion fields can be incorporated Li et al. (2005, 2007); Fedosov et al. (2010); DiMilla et al. (1991); Gracheva and Othmer (2004); Bottino et al. (2002); Dokukina and Gracheva (2010). Spring networks offer the additional advantage that homogenization techniques can be used to arrive at continuum models. Recently, the interplay between force generation and the geometrical and adhesive properties of the environment have been addressed using the powerful framework of finite element models (FEM) Deshpande et al. (2006); Legant et al. (2009a), which can be considered as the continuum limit of appropriate network models. Most FEM-models use constitutive equations which correspond to Hookean networks.
Although conceptually very appealing, modeling cell mechanics with Hookean networks does not reflect the fact that the actin cytoskeleton does not provide much resistance to compression. This is especially true for two-dimensional networks for cell adhesion and migration, because in this case the network might contract laterally, while the cytosol flows into the third dimension. In this situation, the network links do not behave as springs, but rather as cables, which are characterized by an asymmetric force-extension relation. There are several microscopic reasons for this effective behavior: not only do thin actin bundles easily buckle under load, they also tend to telescope in due to filament sliding and even to depolymerize once tensile stress is relieved. Cable networks have been successfully used to model the prestress-dependent mechanical response of adherent cells to local mechanical perturbations Coughlin and Stamenović (2003). The same model has also been used to describe how mechanical stress is propagated from the nuclear region through the cytoskeleton towards focal adhesions, where changes in load lead to changes in adhesion size Paul et al. (2008).
One striking feature of strongly adhering cells is the fact that retracted contours often take the shape of circular arcs Zand and Albrecht-Buehler (1989); Bar-Ziv et al. (1999); Lehnert et al. (2004); Bischofs et al. (2008). Although for cable networks the resulting shapes are strongly invaginated, it has been shown that the circular invaginations observed for cells pinned at discrete sites of adhesions can only be explained if an additional contractile force is introduced for each mechanical link Bischofs et al. (2008). This additional force in an actively contracting cable network does not vanish at the resting length and represents the fact that contractility arises mainly from myosin II motors, which in steady state operate close to a non-vanishing stall force.
In this paper, we systematically compare the different network models introduced before for adherent cells (spring models, cable networks, and actively contracting networks) in regard to the predicted shapes and force distributions. Our main conclusion is that actively contracting cable networks share many crucial features with adherent cells. Due to their linear nature, actively contracting spring networks are equivalent to passive spring networks with a reduced resting length. Passive networks (both from springs and cables) have a well-defined reference state even in the absence of adhesion constraints and in general give similar results regarding shape and force distribution, which is determined mainly by global inputs like the spatial distribution of the adhesion points. In contrast, actively contracting cable networks do not have a well-defined reference state because without adhesion constraints, they contract onto a point. In this case, we find that shape and force distributions are determined mainly by the local distribution of adhesion sites. The internal force distribution is constant in the bulk and strongly localizes to the contour, where forces jump by orders of magnitude. This motivates a detailed study of two contour models, which allow us to derive analytical predictions which we then compare with the results from the computer simulations. We also discuss how actively contracting cable networks can be extended to model non-linear or adaptive linker mechanics, and comment on the relation of our network models to tissue mechanics.
Ii Network Models
ii.1 Link Mechanics
Tissue cells adhering to discrete adhesion points on a flat substrate usually become very flat and therefore effectively two-dimensional. Only the nucleus, which rises in the middle, makes them fully three-dimensional, see sketches in Fig. 1a,b. However, here we focus on lateral contraction and contour effects and thus the nucleus is expected to play a minor role. In the following we therefore restrict ourselves to two dimensions and model the cytoskeletal network as a two-dimensional mesh of mechanical links joined at discrete nodes. Nodes are labeled with indices, e.g. and . Microscopically links and nodes may represent filament bundles and local accumulation of cross-linkers, respectively, but in a more general sense, these mechanical elements are simple representatives of an unknown network architecture which we model in a statistical sense. The network is subject to internal forces originating from molecular motor activity, , and the mechanical resistance of filaments to strain, , compare Fig. 1c.
We introduce three fundamentally different kinds of mechanical models for the network links. The simplest case is a Hookean spring network (HSN) composed of links with resting length , which represent linear springs with spring constant , where is the Young’s modulus of the link and its cross-section. The restoring force acting on a node due to elastic strain in the link then reads:
(1) |
where is the dimensionless unit vector along the link . Here, and specify the node positions, is the length of the link and is the strain in the link. The linear force-extension curve of a single link in a HSN is shown in Fig. 1d as a dashed line with short dashes.
The mechanical properties of the cytoskeleton are attributed mainly to the actin part. Actin is a semi-flexible filament prone to buckling under compression and thus behaves like a cable, which can be stretched but not compressed. The Hookean assumption of a symmetric elastic response is therefore not valid. The mechanical properties of actin networks on a coarse-grained scale are more accurately described by assuming a finite resistance of filaments to tension, but no resistance to compression. The mechanical restoring forces originating from a link connecting two nodes and in the passive cable network (PCN) are therefore given by:
(2) | |||||
(3) |
We show the asymmetric force-extension relation of the PCN links in Fig. 1d as a dashed line with long dashes.
Let us assume that molecular motors are homogeneously distributed in the network. Because they are arranged in a parallel fashion, their individual forces add up. We therefore assume that a link contracts with a force proportional to its initial length, where is force density per length. For computational simplicity, here we assume that this force does not change as the filament contracts, although in practise, it might well be that the line density rather than the total number of active motors is constant. For a link we therefore have:
(4) |
where is the tensile force per initial length applied by the motors. The finite force at zero length is unphysical and we avoid it by introducing an additional rule such that force is diminished if two neighboring nodes come closer to each other than some small distance :
(5) |
Combining PCN and active contraction, we obtain what we call the active cable network (ACN), compare the solid line in Fig. 1d:
(6) | |||||
(7) | |||||
(8) |
Neither the detailed choice of nor the assumption of a linear force reduction below are crucial for our results.
Active contraction can also be combined with the HSN. However, this simply shifts the straight dashed line in Fig. 1d, i.e. this reduces the resting length to
(9) |
as long as . Therefore active contraction does not change the basic definition of the HSN and therefore in the following the actively contracting HSN is not discussed further.
ii.2 Mechanical Equilibrium
For nodes within the network, the total force exerted on a node is the sum of all forces applied by neighboring nodes :
(10) |
In mechanical equilibrium, the force on each non-adherent node has to vanish, . The existence of adhesion sites is modeled by fixing the positions of the adherent nodes. Thus the adhesion site geometry will enter through the boundary conditions.
To reduce the number of parameters, we scale all lengths with respect to , e.g. . All forces are scaled as . We define the ratio of active to elastic forces as
(11) |
For an ACN we therefore may rewrite the forces acting on a node in the non-dimensionalized form as:
(12) | |||||
(13) | |||||
(14) |
In the computer simulations we use .
Mechanical equilibrium requires the forces on each non-adherent node to vanish
(15) |
with the summation over all neighbors. For a two-dimensional network of non-adherent nodes, the system of equations (15) consists of coupled non-linear equations. If the lhs of system (15) consisted of arbitrary functions of the , the method of choice to solve it would be the Newton-Raphson method Gerald and Wheatley (1989). However, since the lhs of system (15) is a force which has a potential, it is also a -dimensional gradient vector. Therefore we solve the minimization problem for the potential with the conjugated gradient method Fletcher and Reeves (1964). We stop iterating as soon as the force on every node (except the periphery nodes) is smaller by at least two orders of magnitude than the smallest link force.
ii.3 Parametrization
Due to its dynamic and multiscale organization and the known limitations of microscopy, a detailed model of the actin cytoskeleton is currently out of bounds. In the face of these uncertainties, our model is not meant to represent the details of the organization of the actin cytoskeleton. Nevertheless, for practical purposes it is helpful to parametrize our model using some benchmark values for the actin cytoskeleton.
The elastic modulus of an actin filament, which has cross-section area , was experimentally found to be Gittes et al. (1993), while typical values for stress fibers are a radius around (corresponding to ) and an effective modulus of Deguchi et al. (2006). Hence, the Young’s modulus of stress fibers is three orders of magnitude smaller than the one of single actin filaments. This suggests that cross-linkers like -actinin and myosin II are the main contributors to elasticity and not the actin filaments themselves. However, the values for the one-dimensional modulus, and , are effectively very similar, so the one-dimensional modulus is expected to be of the order of .
The mesh size of the cytoskeleton is expected to be typically around . This is an intermediate value introduced in Coughlin and Stamenović (2003) based on experimental observations of the actin cytoskeleton in adherent endothelial Satcher Jr and Dewey Jr (1996) and fibroblast cells Luby-Phelps et al. (1986). For the active force, we estimate that around 1.000 myosin II motors are active in one effective link. With a stall force of 2 pN per motor head Molloy et al. (1995), we have for the motor force per length. Using Eq. (11) and , we estimate for the active tension in the network.
Iii Results
iii.1 Equilibrium Shapes
In the following we analyze how the different types of networks act under tension for given adhesion constraints. We first discuss passive networks. In order to understand the role of geometry, as tension-free reference states we use both a square and a triangle as shown in Fig. 2. If the corners are fixed and the network is set under tension by reducing rest length, invaginated shapes appear as shown in Fig. 3 for Hookean spring networks (HSN) and passive cable networks (PCN). Contraction is quantified by . In Fig. 3, the results for HSN and PCN are the same for the square shape, Fig. 3a,c, because these two kinds of networks behave identical as long as all links are tensed. However, the results are different for the triangle shape, Fig. 3b,d. In this case, the PCN gives a significantly flatter contour due to the missing response to compression in the thin extensions leading to the adhesion points. The most prominent examples for compressed links in Fig. 3b,d are indicated by arrows.
In order to understand our numerical results in more detail, we first note that the HSN with triangular network topology has a well-defined continuum limit, in which it corresponds to a two-dimensional sheet with isotropic linear elastictiy Boal (2002); Hansen et al. (1996). The two corresponding elastic constant are a Young’s modulus of and a Poisson’s ratio of , where is the spring constant of the links. The HSN with simple cubic topology does not have such a rigorous limit, but in our context it works in a similar way as the triangular lattice. Therefore similar results as obtained here for the HSN are also obtained with continuum elasticity theory applied to two-dimensional cell shapes Deshpande et al. (2006). Without any adhesion constraint, the HSN contracts isotropically to a finite size, i.e. the network is uniformly scaled and has a new side length . This shape we call the unconstrained reference shape and it is key to understand the results for HSN. The same shape as shown in Fig. 3a results if the network contracts away from its initial state under adhesion constraints or if the network starts from its unconstrained reference state and its corners are dragged to the desired adhesion points. This explains the main feature of the force distribution shown by the color coding in Fig. 3a, namely the strong localization of stresses and strains to the regions around the adhesion points. With the amount of tension used here, the network can attain its unconstrained reference state away from the adhesion points and therefore its contour is essentially flat in the middle parts.
In Fig. 4a we directly compare the calculated network shapes for the HSN for the square geometry to the unconstrained reference shape. In addition, we demonstrate the role of the link length . As explained above, in this case the PCN gives the same results. For the small value of tension, (lower set), the contracting network can reach the unconstrained reference shape over a large region where it is therefore essentially flat. The smaller , the faster this contour is reached due to an increased force density along the contour. For the large value of tension, (upper set), the unconstrained reference shape is not reached by the contracted network and it stays non-flat along its whole contour even for rather small values of .
Next we turn to actively contracting networks. Figs. 5a-d show the equilibrium shapes of an active cable network (ACN) with reference state from Fig. 2a and increasing tension . As tension increases, the shape becomes more and more invaginated, until it collapses onto the zero-area network in Fig. 5d. This network basically consists of a centrally contracted region which is connected to the adhesion points by long extensions. We therefore call this network the center tree (CT). Note that this network still retains aspects of the two-dimensional network, because an effectively one-dimensional structure would collapse onto the so-called Steiner tree of minimal length, which for a square shape is not four-fold symmetric Hildebrandt and Tromba (1996).
In contrast to the passive networks, where tension ceases as the unconstrained reference state is reached, for the ACN no such unconstrained reference shape exists and without adhesion constraints the shape would collapse onto a single point. This explains why flat parts are not observed in the contours of the networks shown in Fig. 5a-d. This is also demonstrated in Fig. 4b, where we show contour shapes for (lower set) and (upper set). Like for the HSN in Fig. 4a, increasing lattice constant increases the force density along the contour and therefore leads to stronger invagination. However, in marked contrast to the HSN, no flat parts appear in the contour as no unconstrained reference state exists.
Fig. 5a-d also shows that the formation of inward directed arcs now corresponds to a much more inhomogeneous density distribution of filaments: in the bulk of the network the distance between nodes and thus the filament density remains unchanged, while at the edges filaments start to bundle strongly along the edge. The color code in Fig. 5 shows that stress is strongly localized at the periphery. In the interior, the only forces acting are the motor forces which balance each other at every node. At the periphery, the force jumps up from to much higher values , compare Fig. 5a. The forces are largest close to the adhesion points and decrease towards the center of the boundary. For large tension , tubes are formed near the adhesion points and the stress distribution along the contour becomes more inhomogeneous.
Fig. 5e shows the effect of changed discretization for the same tension value as Fig. 5a for the complete network. Fig. 5f demonstrates that for ACN, shape and force values do not depend significantly on network topology. As an instructive example here we use a disordered network topology obtained by a Voronoi construction. Even the presence of relatively large elements in the discretization does not change the invaginated shape feature of the contracted network. We conclude from Fig. 5e and Fig. 5f that ACN are surprisingly robust in regard to the details of the network topology. In general, similar results as obtained in Fig. 5 for the square shape are also obtained for other initial cell shapes.
In Fig. 6, we investigate another striking property of ACN, namely its robustness in regard to addition of new adhesion points. In Fig. 6a,b we show the equilibrium shapes of the HSN from Fig. 3a and the ACN from Fig. 5b with one adhesion point added in the middle of the bottom line. In the case of the ACN model, this change at the bottom of the network has little influence on the positions of nodes not directly connected to the bottom line. In contrast, for the HSN the additional adhesion point affects the shape of the opposite arc, becoming more curved in the center. In Fig. 6c we plot by which distance the nodes in the vertical middle line are pulled down in the negative y-direction upon addition of the new adhesion point. The plot of versus is not smooth for numerical reasons, but clearly shows that the effect decays much more rapidly for the ACN versus the HSN. In addition to square and triangular network topologies, here we also show results for rotated square (diamond) and hexagonal networks. Intriguingly, stress in the contour behaves very differently, compare the color coding of Fig. 6a,b. While in the HSN stress in the bottom line stays approximately the same, in the ACN it decreases to half its value, indicating a strong effect on contour forces.
In summary, ACNs behave very differently from HSNs (and therefore also from the mostly equivalent PCNs). Roughly speaking, they respond more locally than globally. They are more robust in regard to network topology and adhesion geometry and show strong localization of the stress to the periphery. Addition of new adhesion points leads to little global change, but to a strong change in local stress distribution. Because for ACN shape and contour stress seem to be mainly determined by the local adhesion geometry, they will now be analyzed in more detail.
iii.2 Contour shape and tension-elasticity model
In contrast to HSNs, the contour of ACNs appears to be more circular. Indeed, a circular arc morphology has been noted before for the shapes of cells adhering to micropatterned substrates and therefore this shape feature is an important motivation to study ACNs Bischofs et al. (2008). We now investigate this important aspect in more detail. In Fig. 7a and b we show contours of the HSN and ACN from Fig. 3 and Fig. 5, respectively, together with fits to straight lines or circular arcs. As we increase tension, we observe more invaginated shapes, compare Fig. 4. For the HSN in Fig. 7a, small tension allows the network to reach the unconstrained reference shape and therefore the best fit to the middle part of the contour is a straight line. For larger tension, the unconstrained reference shape cannot be reached anymore and circular shapes become better fits. This crossover is in marked contrast to the ACN from Fig. 7b, where circular arcs fit very well for all values of . For large , the overall contour starts to deviate from the perfect arc shape because the networks starts to collapse into tubes near the adhesion points. However, locally (in between the tubes) the contour stays circular. Another difference between the two network types lies in the observation that for HSN, network shape strongly depends on lattice constant , while for ACN, the equilibrium contour is relatively independent of .
It has been argued before that the circular arc shape feature of the ACNs can be explained by an analytical theory, the tension-elasticity model (TEM) Bischofs et al. (2008). For clarity, here we repeat this analysis and compare it in detail with our network simulations. Because ACNs do not propagate compression and the motor forces represent a constant pull in the network, in the TEM the bulk contractility is modeled by a structure-less surface tension . However, elasticity is crucial to understand how the contour reacts to the internal pull. Therefore the elastic nature of the mechanical network is represented by an elastic line tension , which prevents the contour from collapsing under the inward pull of the bulk network. This line tension is written as
(16) |
where is the contour length and is its resting length. Note that we use the same value like for the single links because the elastic line tension will be dominated by the contribution from the most peripheral line of links. We further assume , where is the initial (spanning) distance between two neighboring adhesion points (which in the simulations above has been chosen to be ) and is a dimensionless resting length parameter (compare Fig. 8a for a schematics). In the following we restrict ourselves to , that is we assume that a completely relaxed contour is straight, but without internal tension. This implies that we neglect the contribution of the active contractility in the periphery to the line tension.
The relation between surface tension and network tension , which depends on network topology and discretization, can be obtained numerically. For this purpose, we simulate the pulling of a rectangular sheet of network. Surface tension then follows as total force on the pulling boundary divided by its width. Due to this normalization, the result does not depend much on link length , thus we use for the simulations. For all considered network geometries, we find a linear relation:
(17) | |||||
(18) | |||||
(19) | |||||
(20) |
The constant for the square lattice is close to because here all links pull essentially perpendicular to the boundary.
Given the forces assumed by the TEM, one can derive the shape of the contour from the force balance. While the surface tension acts in the direction of the normal , resulting in a pulling force , the line tension acts in the tangential direction (compare Fig. 8a). Because the elastic line tension is a global quantity, it does not vary with the contour length and therefore the contour tension is . Then the force balance reads
(21) |
where for the second part we have used the geometrical relation with being the radius of curvature. We thus conclude that the TEM predicts circular arcs with a radius
(22) |
Although this results looks like a simple Laplace law in two dimensions, it is more complicated, because the arc radius will depend on global properties like spanning distance through the elastic line tension from Eq. (16).
In order to arrive at an expression for arc radius as a function of adhesion geometry and network tension, we use the trigonometric relation
(23) |
to replace contour length by spanning distance , compare Fig. 8a. In combination with Eq. (16) (in dimensionless form) and Eq. (22), this gives
(24) |
Since this equation cannot be solved analytically for , it has to be solved numerically for given values of and . For geometrical reasons, must always be larger than . Therefore, a critical exists above which Eq. (24) cannot be solved anymore. For small values of , the invagination is small and one can expand the geometrical relation in . This leads to the analytical result
(25) |
In Fig. 9 we compare the results from computer simulations for over a large range of network tension to the results of the TEM with the corresponding range of surface tension , both for the numerical solution of Eq. (24) and the analytical solution Eq. (25). We note that the predicted power law behavior applies over a very large range of tensions, and only breaks down at very large tension , compare Fig. 9. The inverse relation between and represents a modified Laplace law for ACNs and thus demonstrates that the concept of an isotropic surface tension works well to explain cell shape. With the linear relation between and , this implies that .
We also find excellent agreement between computer simulations and TEM upon variation of spanning distance, compare Fig. 9. Thus the elastic effects mediated by the spanning distance are captured well by the concept of an elastic line tension. In summary, the analytical TEM results in a surprisingly good description of the contour shape of ACNs. As we will discuss in the next section, however, agreement is less good regarding contour forces.
iii.3 Contour forces and elastic catenary model
We now discuss the forces resulting from our computer simulations. In the simulated networks of Fig. 3 for HSN and Fig. 5 for ACN, the stresses in the network are coded by colour. In Fig. 10a we plot these stresses along a line horizontally crossing the networks with square shapes slightly above the middle line. In the HSN, stress gradually decays into the sample, while for ACN, it jumps up at the periphery. Fig. 10a also shows that forces in HSNs are much larger than those in ACNs with a comparable equilibrium shape. The stress distribution in the boundary of an ACN depends on the lattice constant of the network, as shown in Fig. 10b. The contour forces are minimal in the middle and increase towards the sides. At an adhesion site, all network forces add up to the overall force being transmitted to the substrate. For both passive and active networks, the force which is exerted on an adhesion site follows a power law as () is increased, see Fig. 10c. For , adhesion force is much smaller in the HSN than in the ACN. For active and passive networks of a comparable shape, however, e.g. and , we observe the opposite behavior.
Although on an absolute scale the variation is not very strong, Fig. 5c and Fig. 10b both demonstrate that for ACNs under large network tension, peripheral force varies along the contour. In contrast, the tension-elasticity model (TEM), which is very successful in explaining shape, predicts homogeneous force along the boundary. Fig. 5c suggests one reason which could explain this discrepancy. For ACN, the links essentially telescope in under contraction and therefore their density along the contour varies for strong curvature along the contour. This suggests that in order to explain the spatially varying force in the contour, one has to revisit the assumption of an isotropic surface tension creating a homogeneous force density along the contour.
As an alternative to the TEM, we now investigate another analytical model, which incorporates the effect of varying link density, namely the elastic catenary Lockwood (2007). In the elastic catenary, the pulling force on the elastic contour is not along the normal, but along the direction perpendicular to the original contour, similar to the situation in networks with square shape and square topology, compare Fig. 8b. Due to the linear elasticity in the contour, the line density of links along the contour varies in proportion to the contour tension. In dimensionless units, a length element is expanded to length in the presence of tangential contour tension . With an initial inward force per unit length , the effective force density is therefore . In contrast to the TEM, we now assume not only a heterogeneous force distribution, but also a vertical pulling direction. Therefore the pulling force is . The force balance again reads
(26) |
but now the solution is more difficult than for the TEM. Because the tangent is normalized, it can be written as , where is the tangential angle along the contour Whewell (1849), compare Fig. 8b. Different from the TEM, we now have to solve two equations:
(27) | |||||
(28) |
Eq. (27) can directly be integrated, leading to , while Eq. (28) can then be solved via the substitution . This gives
(29) | |||||
(30) |
The integration constants and are determined by the positions of the adhesion sites.
While in the TEM we assume the line tension to be constant along the whole boundary line, for the elastic catenary only the -component of tension, , is constant. Its value can be determined numerically from the above equations for given and . One can show that for small surface tension , the contour becomes parabolic with a radius of curvature . For the relative deviation to the prediction of the TEM is only , that is in this regime, the elastic catenary model leads essentially to the same result as the TEM with circular arcs. However, in contrast to the TEM, this model predicts a spatially varying boundary tension of
(31) |
The curve without symbols in Fig. 10b shows that this model qualitatively predicts the observed minimum in the stress distribution. The force acting on an adhesion dot is predicted to be
(32) |
Fig. 10d shows that this prediction is quite accurate.
iii.4 Strain Stiffening
So far, we have treated the network links as cables or springs, in which force increases linearly with elongation. This implies that the elastic modulus of the network is independent of strain. However, it is well known that the cytoskeleton of cells shows a strong increase of elastic modulus with strain (strain stiffening) Janmey (2011). This property of cytoskeletal networks has been reconstituted in vitro with crosslinked actin filaments of physiological length, which is of the order of their persistence length (around 1 ) Gardel et al. (2006). Under these conditions, strain stiffening results mainly from the mechanical properties of the crosslinkers. For larger filament lengths and stiff crosslinkers, there exists another mechanism for strain stiffening, namely the non-linear force-extension curve of semiflexible polymers as described by the worm-like-chain (WLC) model Doi and Edwards (1986). The WLC-model has been used before to model semiflexible biopolymers like DNA Marko and Siggia (1995), actin Kroy and Frey (1996); Heussinger and Frey (2006) and spectrin Li et al. (2005, 2007). Because our network model is especially suited to study the effect of link mechanics, we now study strain stiffening based on the WLC-model, although for the physiological relevant case, a model for crosslinker mechanics might be more relevant.
While the WLC proper has vanishing resting length, here we combine it with a finite resting length to also include the effect of compression. Thus we use the non-linear WLC model to describe the mechanical links as they are tensed away from their reference state, while the compressed state is modeled as above (linear response for springs and no response for cables). Complementing Eqs. (2) and (3) for the PCN, we get in dimensionless form
(33) | |||||
(34) |
Force is now given as multiples of , while length is again scaled with . Strain is now defined as , in which gives the dimensionless difference between maximal extension and reference length. Effectively there is only one difference to the original model, namely the additional term which diverges if the strain approaches . Without strain, this term vanishes. For ACNs, the original model is extended in the same way.
We choose , that is the maximal extension is . For small tension, , the non-linearity does not affect the shape of the PCN much and we observe the same invagination as in Fig. 3. Fig. 11a,b shows for the PCN a comparison between linear and WLC-networks for a large value of tension, . Obviously the strain-stiffened network shows a much larger resistance to invagination. We also note that forces are two orders of magnitude larger in the non-linear model. ACNs are affected less by the non-linearity, as shown in Fig. 11c. Here we use and again . Comparison with Fig. 5d reveals that the ACN collapses to a lesser degree than without strain stiffening. This can be quantified by the arc height , defined as the maximum distance between initial and current edge in the equilibrium shape. is significantly reduced by the non-linearity, Fig. 11d. Thus much higher motor forces are needed to reach the collapsed state.
In Fig. 12a, the contour of the strain stiffening PCN is analyzed in more detail. Circular and linear fits are shown as lines. The bottom two are given by straight lines. For the contour of the PCN-WLC is qualitatively the same as that of the linear one, shown in Fig. 7b. At the contour cannot be fitted well by circle or line. If is increased beyond , the network does not contract any further, but again expands outward. This surprising effect does not occur for ACNs. For , the arcs appear to be circular. The ACN-WLC contour, Fig. 12b, only differs little from the linear ACN contour, Fig. 7b. Arcs are always circular (except at the regions where tubes form). With increasing they continuously move inward. Comparison with Fig. 7a reveals that radii typically are larger in the non-linear case.
iii.5 Link Adaption
Adherent cells are known to strongly adapt their cytoskeleton to the physical properties of their environment. During recent years, it has become clear that the actin cytoskeleton tends to reinforce under load. In addition, mechanical loading of adhesion contacts leads to regulatory signals which increase myosin motor activity inside the cell. Network models are especially suited to study these biologically important effects in a theoretical framework. In the following, we will investigate which changes occur in the network if the elastic constant and the force density resulting from myosin II activity are increasing with load.
For simplicity, we assume that both and first increase with force in a linear fashion and then saturate at constant values, which is the simplest assumption for a process based on enzymatic regulation:
(35) | |||||
(36) |
where the force scale determines when half the maximal increase has been reached. We again use dimensionless parameters. Forces , , are measured in units , i.e. we define (the same for and ), (the same for , ) and (the same for ).
Fig. 13 demonstrates that the effect of adaptation is fundamentally different for passive versus active networks. Here stiffer links are represented by thicker lines. Fig. 13a shows the result for a HSN (as shown by Fig. 3, for the square shape a PCN gives the same results). In this case, only Eq. (35) must be considered. While the peripheral links show the largest values of , the rigidity decays smoothly from there into the bulk. This is strikingly different for the ACN shown in Fig. 13b. Here Eqs. (35) and (36) have been used. We find that the stiffness is strongly localized to the periphery, as found before for the internal stress. Thus for passive networks the adaption response is spatially continuous, while for actively contracting networks, it is strongly localized to the rim. This nicely agrees with experimental observations that strong peripheral actin bundles typically line the cell contour Bischofs et al. (2008). In particular, our model suggests that this effect is strongly determined by the mechanical properties of the underlying networks.
iii.6 Relation to Tissue Shape
Tissue contraction with discrete pinning sites is very similar to cell contraction since adherent tissues also show invaginated arcs Bischofs et al. (2008); Legant et al. (2009b). However, in this case the spatial dimensions of the pinning objects tend to be relatively large. In order to include this effect, we have simulated a circular network contracting around four circular dots of finite size, see Fig. 14a. All nodes on the circles are fixed in space. Motor force is increased stepwise and nodes coming closer to each other than are glued together and in the following act as one Paul et al. (2008). Relative dimensions are taken from Legant et al. (2009b), where a microtissue tethered to four cylindrical posts is analyzed. With the ACN we are able to reproduce the typical arc morphology of the contracted microtissues, see Fig. 14b. We note that this is not possible with a FEM approach, as this leads to flat contours as for HSN Legant et al. (2009b). Thus ACN are a useful model both for cells and tissues.
Iv Conclusions and Outlook
Motivated by the network nature of the actin cytoskeleton and its effectively 2D organization in mature adhesion to flat substrates, we have modeled adherent tissue cell contraction by 2D network models. The main aim of this work is to achieve a detailed comparison of the shapes and force patterns for different network types, namely Hookean spring networks (HSNs), passive cable networks (PCNs), and active cable networks (ACNs).
The shape of a HSN can be understood best by considering the shape of its unconstrained reference shape. If tension is not too large, the network contour follows the unconstrained reference shape at regions sufficiently far away from the adhesion sites. Closer to the adhesion sites, the network deforms and stress and strain accumulate. In contrast, the ACN does not have an unconstrained reference shape and without adhesion constraints would contract into a point. Therefore no signature of the unconstrained reference shape (like flat parts for a square-shaped lattice) appear in the contour. Because it does not resist compression, stress and strain are not propagated much into the network and are strongly localized to the contour.
One of the most striking difference between the different network types revealed by our analysis is the fact that in passive networks, local changes to the adhesion geometry changes the network globally. This is in marked contrast to the active network, where the addition of local adhesions has only a local effect on the boundary. However, in this case the change in spanning distance has a large effect on the stress in the contour, as predicted by the tension-elasticity model (TEM). The TEM is especially suited to quantitatively predict the shape of an ACN, namely the circular arcs observed between neighboring adhesion points and the scaling of their radius with spanning distance and surface tension.
Despite this success, the TEM does not capture all aspects of the network model. While the TEM assumes constant contour tension, the computer simulations reveal that tension varies along the contour. An elastic catenary theory qualitatively predicts that tension decreases towards the middle of an invagination due to local changes in link density along the contour. However, it does neither predict the quantitative details of the contour stress nor the circular arc morphology.
In the case of very large tension, both network types develop different features. For the passive networks, the invaginations tend to become more round, as the unconstrained reference shape becomes so small that the contour cannot reach it anymore. In contrast, the active network develops straight features, because the network collapses into tubes at the adhesion points. Indeed the formation of tubes has been observed experimentally and eventually leads to pearling through a Rayleigh-Plateau instability Bar-Ziv et al. (1999). The region between the tubes always stays circular for ACN.
Network models are ideally suited for multi-scale modeling because physical properties can be easily added on the level of single link and lead to non-trivial effects on the level of cell shape and forces. In order to demonstrate this important aspect, we have studied two important additional features of the cytoskeleton. First non-linear links were introduced via the worm-like-chain model (WLC). In the WLC case the passive square network first contracts and then expands again as tension is increased. The ACN requires much larger values of tension to contract compared to the linear case. Otherwise the arc morphology is essentially the same as in the linear case.
As a biologically very relevant aspect of the cellular cytoskeleton, we also have studied the adaptation response of network links. For certain parameter values, a saturation response for both elasticity and tension leads to a strong difference between and of boundary links, which are strongly increased, and and of internal links, which are increased much less. Similar aspect have been addressed before in the framework of Finite Element Modeling (FEM) Deshpande et al. (2006). In this case, the biochemical regulation has been modeled with more detail. Both the resulting cell shapes and the formation of stress fibers inside the cell demonstrate that the FEM model strongly resembles the HSN studied here. Therefore it would be interesting to combine the detailed biochemical model with the actively contracting cable network studied here. While it appears to be very challenging to develop a homogenization strategy for the ACN, it is interesting to consider if similar features as resulting from the ACN could be obtained in a FEM-framework.
Active cable networks have also been shown to describe the circular arc morphology of tissues pinned at discrete sites Bischofs et al. (2008). Because here arc radius also scales with spanning distance as for the arc radius of strongly adhering cells, the tension-elasticity model seems to capture all essential element of this situation. In the tissue case, the cable network represents the fibrous nature of the collagen matrix and the active contractility corresponds to cell contraction. Because in addition water can flow out of the contracting cell-matrix composite, volume is not conserved and compression is not propagated. Therefore the standard models of elasticity are strictly speaking not appropriate. Indeed they do not predict the circular arc morphology, but rather show flat contours corresponding to the unconstrained reference shape of the elastic model Legant et al. (2009b).
In summary, HSN, PCN and ACN are simple model systems which however show surprisingly rich responses to internal contractility and therefore lead to interesting conclusions about the physical elements required to endow cells with a sense of geometry. ACN seem to be very appropriate to model strongly adherent tissue cells as they not only implement some of the most important fundamental features of the cytoskeleton (asymmetry under tension and compression, contraction by molecular motors), but also lead to functions which are very reminiscent of real cells (robustness under structural re-arrangements and adaptation to local adhesion constraints).
Acknowledgements.
USS is a member of the Heidelberg cluster of excellence CellNetworks and was supported by the MechanoSys-grant from the Federal Ministry of Education and Research (BMBF) of Germany. IBB was supported by the German Science Foundation (DFG) through the Emmy Noether program (BI1213-3/1).References
- P. A. Janmey, J. P. Winer, M. E. Murray, and Q. Wen, Cell Motil. Cytoskel. 66, 597 (2009).
- B. Geiger, J. P. Spatz, and A. D. Bershadsky, Nat. Rev. Mol. Cell Biol. 10, 21 (2009).
- R. Zaidel-Bar, S. Itzkovitz, A. Ma’ayan, R. Iyengar, and B. Geiger, Nat. Cell Biol. 9, 858 (2007).
- R. Zaidel-Bar and B. Geiger, J. Cell Sci. 123, 1385 (2010).
- B. D. Hoffman, C. Grashoff, and M. A. Schwartz, Nature 475, 316 (2011).
- M. Dembo and Y.-L. Wang, Biophys. J. 76, 2307 (1999).
- J. P. Butler, I. M. Tolić-Nørrelykke, B. Fabry, and J. J. Fredberg, Am. J. Physiol. Cell Physiol. 282, C595 (2002).
- U. S. Schwarz, N. Q. Balaban, D. Riveline, A. Bershadsky, B. Geiger, and S. A. Safran, Biophys. J. 83, 1380 (2002).
- J. L. Tan, J. Tien, D. M. Pirone, D. S. Gray, K. Bhadriraju, and C. S. Chen, Proc. Natl. Acad. Sci. USA 100, 1484 (2003).
- B. Sabass, M. L. Gardel, C. M. Waterman, and U. S. Schwarz, Biophys. J. 94, 207 (2008).
- Z. Liu, J. L. Tan, D. M. Cohen, M. T. Yang, N. J. Sniadecki, S. A. Ruiz, C. M. Nelson, and C. S. Chen, Proc. Natl. Acad. Sci. USA 107, 9944 (2010).
- V. Maruthamuthu, B. Sabass, U. S. Schwarz, and M. L. Gardel, Proc. Natl. Acad. Sci. USA 108, 4708 (2011).
- I. B. Bischofs, S. S. Schmidt, and U. S. Schwarz, Phys. Rev. Lett. 103, 048101 (2009).
- S. W. Grill, P. Gönczy, E. H. K. Stelzer, and A. A. Hyman, Nature 409, 630 (2001).
- S. W. Grill, J. Howard, E. Schäffer, E. H. K. Stelzer, and A. A. Hyman, Science 301, 518 (2003).
- S. Kumar, I. Z. Maxwell, A. Heisterkamp, T. R. Polte, T. P. Lele, M. Salanga, E. Mazur, and D. E. Ingber, Biophys. J. 90, 3762 (2006).
- K. Tanner, A. Boudreau, M. J. Bissell, and S. Kumar, Biophys. J. 99, 2775 (2010).
- J. Colombelli, A. Besser, H. Kress, E. G. Reynaud, P. Girard, E. Caussinus, U. Haselmann, J. V. Small, U. S. Schwarz, and E. H. K. Stelzer, J. Cell. Sci. 122, 1665 (2009).
- M. A. Smith, E. Blankman, M. L. Gardel, L. Luettjohann, C. M. Waterman, and M. C. Beckerle, Dev. Cell 19, 365 (2010).
- M. Mayer, M. Depken, J. S. Bois, F. Jülicher, and S. W. Grill, Nature 467, 617 (2010).
- D. Riveline, E. Zamir, N. Q. Balaban, U. S. Schwarz, T. Ishizaki, S. Narumiya, Z. Kam, B. Geiger, and A. D. Bershadsky, J. Cell. Biol. 153, 1175 (2001).
- R. Paul, P. Heil, J. P. Spatz, and U. S. Schwarz, Biophys. J. 94, 1470 (2008).
- P. Heil and J. P. Spatz, J. Phys. Condens. Matter 22, 194108 (2010).
- D. Kirchenbüchler, S. Born, N. Kirchgeßner, S. Houben, B. Hoffmann, and R. Merkel, J. Phys. Condens. Matter 22, 194109 (2010).
- J. Howard, Mechanics of Motor Proteins and the Cytoskeleton (Sinauer Associates, Inc., Sunderland, Massachusetts, 2001).
- D. Boal, Mechanics of the Cell (Cambridge University Press, Cambridge, 2002).
- J. C. Hansen, R. Skalak, S. Chien, and A. Hoger, Biophys. J. 70, 146 (1996).
- J. C. Hansen, R. Skalak, S. Chien, and A. Hoger, Biophys. J. 72, 2369 (1997).
- S. K. Boey, D. H. Boal, and D. E. Discher, Biophys. J. 75, 1573 (1998).
- D. E. Discher, D. H. Boal, and S. K. Boey, Biophys. J. 75, 1584 (1998).
- G. Lim H. W., M. Wortis, and R. Mukhopadhyay, Proc. Natl. Acad. Sci. USA 99, 16766 (2002).
- H. Noguchi and G. Gompper, Proc. Natl. Acad. Sci. USA 102, 14159 (2005).
- J. Li, M. Dao, C. T. Lim, and S. Suresh, Biophys. J. 88, 3707 (2005).
- J. Li, G. Lykotrafitis, M. Dao, and S. Suresh, Proc. Natl. Acad. Sci. USA 104, 4937 (2007).
- D. A. Fedosov, B. Caswell, and G. E. Karniadakis, Biophys. J. 98, 2215 (2010).
- P. A. DiMilla, K. Barbee, and D. A. Lauffenburger, Biophys. J. 60, 15 (1991).
- M. E. Gracheva and H. G. Othmer, Bull. Math. Biol. 66, 1470 (2004).
- D. Bottino, A. Mogilner, T. Roberts, M. Stewart, and G. Oster, J. Cell Sci. 115, 367 (2002).
- I. V. Dokukina and M. E. Gracheva, Biophys. J. 98, 2794 (2010).
- V. S. Deshpande, R. M. McMeeking, and A. G. Evans, Proc. Natl. Acad. Sci. USA 103, 14015 (2006).
- W. R. Legant, J. S. Miller, B. L. Blakely, D. M. Cohen, G. M. Genin, and C. S. Chen, Nature Methods 7, 969 (2010a).
- M. F. Coughlin and D. Stamenović, Biophys. J. 84, 1328 (2003).
- M. S. Zand and G. Albrecht-Buehler, Cell Motil. Cytoskel. 13, 195 (1989).
- R. Bar-Ziv, T. Tlusty, E. Moses, S. A. Safran, and A. Bershadsky, Proc. Natl. Acad. Sci. USA 96, 10140 (1999).
- D. Lehnert, B. Wehrle-Haller, C. David, U. Weiland, C. Ballestrem, B. A. Imhof, and M. Bastmeyer, J. Cell Sci. 117, 41 (2004).
- I. B. Bischofs, F. Klein, D. Lehnert, M. Bastmeyer, and U. S. Schwarz, Biophys. J. 95, 3488 (2008).
- C. F. Gerald and P. O. Wheatley, Applied Numerical Analysis (Addison-Wesley Publishing Company, Inc., Reading, 1989).
- R. Fletcher and C. M. Reeves, Comp. J. 7, 149 (1964).
- F. Gittes, B. Mickey, J. Nettleton, and J. Howard, J. Cell Biol. 120, 923 (1993).
- S. Deguchi, T. Ohashi, and M. Sato, J. Biomech. 39, 2603 (2006).
- R. L. Satcher Jr and C. F. Dewey Jr, Biophys. J. 71, 109 (1996).
- K. Luby-Phelps, D. L. Taylor, and F. Lanni, J. Cell Biol. 102, 2015 (1986).
- J. E. Molloy, J. E. Burns, J. Kendrick-Jones, R. T. Tregear, and D. C. S. White, Nature 378, 209 (1995).
- S. Hildebrandt and A. Tromba, The Parsimonious Universe: Shape and Form in the Natural World (Springer, New York, 1996).
- E. H. Lockwood, A book of curves (Cambridge University Press, Cambridge, 2007).
- W. Whewell, Trans. Cambridge Philos. Soc. 8, 659 (1849).
- Q. Wen and P. A. Janmey, Curr. Opin. Solid State Mat. Sci. 15, 177 (2011).
- M. L. Gardel, F. Nakamura, J. H. Hartwig, J. C. Crocker, T. P. Stossel, and D. A. Weitz, Proc. Natl. Acad. Sci. USA 103, 1762 (2006).
- M. Doi and S. F. Edwards, The theory of polymer dynamics (Clarendon, Oxford, 1986).
- J. F. Marko and E. D. Siggia, Macromolecules 28, 8759 (1995).
- K. Kroy and E. Frey, Phys. Rev. Lett. 77, 306 (1996).
- C. Heussinger and E. Frey, Phys. Rev. Lett. 96, 017802 (2006).
- W. R. Legant, A. Pathak, M. T. Yang, V. S. Deshpande, R. M. McMeeking, and C. S. Chen, Proc. Natl. Acad. Sci. USA 106, 10097 (2009b).