Elastic Anisotropy Governs the Decay of Cell-induced Displacements
Page of Elastic Anisotropy Governs the Decay of Cell-induced Displacements
Shahar Goren , Yoni Koren , Xinpeng Xu , Ayelet Lesman School of Mechanical Engineering, Faculty of Engineering, Tel-Aviv University, Israel Physics Program, Guangdong Technion—Israel Institute of Technology, Shantou, Guangdong 515063, People’s Republic of China *Correspondence should be addressed to Ayelet Lesman ([email protected]) and Xinpeng Xu ([email protected])
Significance
Tissues are made up of cells and an extracellular matrix (ECM), a cross-linked network of stiff biopolymers. Cells actively alter the ECM structure and mechanics by applying contractile forces, which allow them to sense other distant cells and regulate many tissue functions. We study theoretically the decay of cell-induced displacements in fibrous networks, while quantifying the changes in the elastic properties of the cell's local environment. We demonstrate that cell contraction induce an anisotropic elastic state, i.e. , unequal principal elastic moduli, in the ECM which dictates the slow decay of displacements. These observations suggest a new mechanical mechanism through which cells can mechanically communicate over long distances, and may provide biomaterials design parameters to guide morphogenesis in tissue engineering.
Abstract
The unique nonlinear mechanics of the fibrous extracellular matrix (ECM) facilitates long-range cell-cell mechanical communications that would be impossible on linear elastic substrates. Past research has described the contribution of two separated effects on the range of force transmission, including ECM elastic non-linearity and fiber alignment. However, the relation between these different effects is unclear, and how they combine to dictate force transmission range is still elusive. Here, we combine discrete fiber simulations with continuum modeling to study the decay of displacements induced by a contractile cell in fibrous networks. We demonstrate that fiber non-linearity and fiber reorientation both contribute to the
Page of strain-induced anisotropy of the elastic moduli of the cell’s local environment. This elastic anisotropy is a “lumped” parameter that governs the slow decay of the displacements, and it depends on the magnitude of applied strain, either an external tension or an internal contraction as a model of the cell. Furthermore, we show that accounting for artificially-prescribed elastic anisotropy dictates the displacement decay induced by a contracting cell. Our findings unify previous single effects into a mechanical theory that explains force transmission in fibrous networks. This work provides important insights into biological processes that involve the coordinated action of distant cells mediated by the ECM, such that occur in morphogenesis, wound healing, angiogenesis, and cancer metastasis. It may also provide design parameters for biomaterials to control force transmission between cells, as a way to guide morphogenesis in tissue engineering. Introduction
Cells in tissues are surrounded by a cross-linked network of semi-flexible biopolymers, such as collagen and elastin, known as the extracellular-matrix (ECM) (1, 2). Cells actively adhere and are mechanically connected to the ECM, which enables them to sense and respond to their mechanical microenvironment by applying active contractile forces (3). These cellular forces can alter the structure and mechanical properties of the ECM in the proximity of the cell (4). As such, cell-ECM mechanical interactions play important roles at the cellular level, such as in migration (5, 6), proliferation (6–8), differentiation (5–7, 9, 10), cancer invasion (11–13) and mineral deposition (14).
When cells contract their microenvironment, they cause substantial displacements and structural changes that can reach a distance of tens of cell diameters away (15–18). Such phenomena support long-range cell-cell mechanical communication, a process that can mechanically couple distant cells and coordinate processes such as capillary sprouting (19) and synchronous beating (20). The long-range transmission of cellular forces is usually attributed to the unique nonlinear mechanics and fibrous nature of the ECM (15–17, 21–29). It is currently clear that multiple physical effects (28, 30), including fiber stiffening (15, 31), buckling (17, 22) and collective network responses ( e.g. , fiber alignment (18, 25, 31))
Page of play central roles in enhancing force transmission in the ECM. Here we show how these effects that were previously considered separately are interlinked and contribute to a more general mechanical mechanism. Specifically, fibrous networks are known to show substantial stiffening under different kinds of tensile loading (30, 32–35), and contracting cells have been shown to stiffen their surrounding environment when embedded in fibrous gels (15, 34, 36). Stiffening of networks can result either from the stretch-stiffening of individual fibers, as described by the semi-flexible polymer model (2, 37, 38), or from collective network alignment due to strain-induced fiber reorientation (23–25, 39). In various biological contexts, it was frequently shown that fibers align in the vicinity of contracting cells and in the matrix between neighboring cells (12, 13, 18, 23, 26, 27, 31, 36, 40, 41). The tendency of the ECM to stiffen and align under tension was shown to facilitate long-range transmission of internal cellular forces (15, 28, 29, 39, 42). In particular, Hall et al. (16) demonstrated that cell contractility induces local collagen alignment which was shown to be highly correlated to enhanced force transmission. Wang et al. (26) devised a continuum model that incorporates the effect of alignment, and managed to reproduce the nonlinear properties of collagen networks and the increased force transmission. Importantly, however, the isotropic stiffening alone is not enough to account for far-reaching forces observed in experiments, and other mechanisms must be involved (27). For example, buckling of fibers under compression results in softening of fibrous network when compressed over some small critical strain (2, 30, 43–45). In previous works (17, 21, 22, 26, 28) including ours (17, 28), fiber buckling was shown to enhance force transmission. Ronceray et al. (22) has shown, for example, that cell contracting in a simulated fibrous networks creates a buckled region near the cell where the decay of force is very slow. Indeed, fiber buckling has been observed near contracting cells in both experiments (46, 47) and network simulations (22, 48), which provide direct evidence for its involvement. Although the effects of buckling and stiffening have been clearly observed, it is still unknown if they are related phenomena. Taking into account that both buckling (21) and stiffening (15, 31 ) lead to changes in the elastic properties of the fibrous environment, may hint into the underlying mechanism in controlling force transmission. Despite the understanding of the effects of fiber nonlinearity (buckling, stiffening) and fiber reorientation to enhancing force transmission, they were previously considered separately. The manner these effects are related and combined to determine the slow decay of cell-induced displacements has yet to be quantified and is still unclear (28, 29). This work combines discrete fiber simulations with continuum modeling to study how the displacements induced by a contracting cell are transmitted in fibrous networks and quantifies how the nonlinear mechanics of individual fibers and collective fiber alignment combine Page of to increase the elastic anisotropy of the cell's local environment. We demonstrate that when a cell pulls on the network, the elastic modulus along the radial direction becomes higher than the modulus along the angular direction, leading to a substantial increase in the elastic anisotropy of the network. This newly established elastic state dictates the decay of displacements in cross-linked fibrous networks. We thus propose a new theory of strain-induced elastic anisotropy through which cells can maintain long-range mechanical communications with one another in fibrous networks consisting of various types of fibers with different properties. We have further validated this mechanism of anisotropy-facilitated displacement transmission by studying the decay of displacements in intrinsically pre-designed anisotropic networks. Our findings provide a new strategy in controlling force transmission in biomaterials by altering the material elastic anisotropy. Results and Discussion Discrete fiber simulations of a contracting cell in an isotropic network
Using a simplified finite element model of a contracting cell in a two-dimensional isotropic fibrous network (Fig. 1A), we explore how the nonlinear mechanical properties of the fibers and the strain-induced fiber alignment affect the decay of displacements. The cell is represented by a circular cavity that is isotropically contracting in the fibrous network. The fibrous networks are designed to be isotropic (in fiber orientation), and homogenous (in fiber density) at the scale of a cell (network construction and nonlinear properties are described in the Supporting Information, Sec. I and II, respectively). Four fiber models are used to compare different mechanical behaviors of the ECM biopolymers: linear, compression-buckling, tension-stiffening and a model including both buckling and stiffening (Fig. 1B). In all cases, nodes are modeled as freely rotating hinges ( i.e. , neglecting the rotational resistance of crosslinkers), and only stretching energy is considered. Here we focus on networks with a high connectivity around eight, where the deformations are found to be near-affine (see section III in Supporting Information for simulation evidence of the nearly-affine network deformation). This is in contrast to previous works that concentrated on the classification of different mechanical regimes of fibrous networks, according to connectivity, bending rigidity and internal pre-stresses (39, 49–53). Note that high-connectivity networks can represent fibrous biopolymer gels that are more vigorously cross-linked (54–56) or combined with synthetic gels (57). We take the advantage of the near-affine deformation in such high- connectivity networks to develop
Page of and compare with continuum models, which allows us to gain deeper insights into the associated mechanisms. To further simplify the system, we perform the simulations in two-dimensions, which has been shown to qualitatively capture all of the main mechanical behaviors of fibrous ECM, while being substantially simpler computationally (2, 22, 58). Two decaying regimes of displacements induced by a contracting cell
To simulate the decay of network displacements induced by a contracting cell, we embed a cell represented by a circular cavity of radius 𝑅 cell in the fibrous network. The cell contracts isotropically by radial displacement 𝑈 cell at its edge (Fig. 1A). From these simulations, we obtain the network displacement 𝑈 as a function of the distance 𝑅 from the center of the cell (Fig. 2A). We consider cell contractions, 𝒞 ≡ −𝑈 cell /𝑅 cell , ranging from 2% to 50%. Distances are normalized by 𝑅 cell (𝑅̃ = 𝑅/𝑅 cell ) and displacements are normalized by 𝑈 cell (𝑈̃ = 𝑈/𝑈 cell ). Since our network is isotropic at the undeformed state and because all the simulated fiber models behave linearly at small strains, we expect to find 𝑈̃~𝑅̃ −1 for small strains, as in linear isotropic elastic continuum (59). We indeed observe this scaling at the far-field 𝑅̃ ≫ 1, where strains are small (Fig. 2C). In contrast, close to the cell, up to a distance of around cell , the displacement decays more slowly. In this near-field regime, strains are high such that fiber alignment and fiber nonlinearities (stiffening/buckling) play substantial roles in slowing the decay of displacements. We note that the slow decay of displacements in the near-field regime can be well fitted by an effective power law 𝑈̃~𝑅̃ −𝑛 , with 𝑛 ≤ 1 for all fiber model types and cell contractions (Fig. 2C). Moreover, we find that as the cell pulls more on the network, the fitted exponent 𝑛 in the near-field decreases, giving rise to longer range transmission of displacements (Fig. 2D). Remarkably, even for networks of linear fibers, there is a monotonic decrease of 𝑛 as the cell contracts more. This is in contrast to linear isotropic continuum in which we expect a constant 𝑛 = 1 , and the decay of displacement is independent of the cell contraction and the elastic properties of the linear medium (59). The range of displacement transmission is further enhanced ( i.e. , 𝑛 decreases more with increasing cell contraction) by introducing fiber buckling and stiffening (inset of Fig. 2D). Interestingly, the exponents 𝑛 for all types of fibers falls to a single master curve (Fig. 2D) when plotted as a function the normalized cell contraction , 𝒞/𝒞 crit , in which 𝒞 crit is the critical cell contraction, a characteristic parameter of the network, at which the near-field region starts to be evident, with n = 0.9 , as previously defined by Xu et al. (28). This indicates that 𝒞/𝒞 crit , is an essential
Page of “emergent” dimensionless parameter involved in the long-range transmission in fibrous networks; it measures the degree of the nonlinearity in the cell-contracted network. The increase in the range of displacement transmission by the nonlinear near-field region can be also quantified by introducing an effective cell contraction, 𝑈̃ eff , through fitting the far-field displacement by 𝑈̃ = 𝑈̃ eff ∗ 𝑅̃ −1 (60). 𝑈̃ eff reflects the degree of cell contraction that would induce the same magnitude of displacements far from the cell in a linear isotropic elastic continuum (Fig. 2E). 𝑈̃ eff increases with contraction and is enhanced when nonlinear mechanical properties are introduced. We find that for strong cell contraction, 𝑈̃ eff scales linearly with the normalized cell contraction, 𝒞/𝒞 crit , for all types of fibers simulated in this work. This unified linear relation for all fiber models is consistent with our theory prediction (see Eq. (S36) in Sec. VI in Supporting Information) and further justifies the importance of the normalized cell contraction,
𝒞/𝒞 crit as an essential “emergent” dimensionless parameter involved in the long-range transmission. Overall, these results indicate that mechanical and geometrical nonlinearities give rise to two power-law regimes with distinct scaling, and nonlinear properties of individual fibers play important roles mainly in the near-field regime. These findings extend on previous studies where multiple scaling regimes were identified (22, 28), by systematically quantifying the effects of the strength of cell contraction and the mechanical properties of fibers on the identified regimes.
Strain-induced elastic anisotropy: effects of fiber alignment, stiffening and buckling
Our next goal is to understand why the near-field regime shows a slower decay of displacements, and what mechanism dictates the slope of the decay. We hypothesize that when cells pull on the network, the elastic modulus along the radial direction becomes much higher than the modulus along the angular direction ( i.e. , the network becomes elastic anisotropic), and that this new elastic state dictates the slow decay of displacements. Our motivation originates from the established theory of anisotropic continuum: when a circular cavity contracts isotropically in a linear anisotropic elastic continuum, the displacements decay as
𝑈̃~𝑅̃ −√𝐸 /𝐸 (21, 61) (Supporting Information, Sec. VI), where 𝐸 and 𝐸 are the stiffness moduli in the radial and angular directions, respectively. Recent experiments by Han et. al. (34) indeed demonstrated that stiffness becomes anisotropic near contracting cells, with a higher extent of stiffening in the radial relative to the angular direction. Moreover, in our simulations, radially aligned fibers are stretched, whereas angularly aligned fibers are typically compressed (Fig. 2B). In the case of nonlinear Page of elastic fibers, this can clearly lead to elastic anisotropy in the network, as radial stretched fibers will become stiffer than the angular compressed fibers. We first carry out simulations to study how an external strain induces elastic anisotropy in a bulk fibrous network. We generate a rectangular piece of the network and subject it to uniaxial stretch (Fig. 3A, B). We then apply an additional infinitesimal uniaxial stretch to obtain the linear (zero-strain) longitudinal modulus 𝐸 . We find that 𝐸 ≈ 𝜇𝜙 f with 𝜇 being the modulus of the fibers and 𝜙 f being the fiber volume fraction. This agrees with the prediction of classical cellular solid theory (44, 62), justifying our methods of measuring elastic moduli. We then measure the elastic modulus in the transverse direction, and verify that it is equal to the longitudinal value, indicating elastic isotropy of the network at infinitesimal strain (Fig. 3C, for 𝜖 ≈ 0 ). For each mechanical fiber model, we stretch this bulk network externally over a range of tensile strains 𝜖 , and allow the network to reach equilibrium, leaving free boundary conditions in the transverse direction. The network then assumes an aligned, anisotropic state (Fig. 3B). Upon each external strain 𝜖 , we measure the elastic moduli 𝐸 and 𝐸 along the longitudinal and transverse directions, respectively, by applying an additional infinitesimal tensile stress along each direction (Supporting Information, Section III and Fig. S5). The elastic modulus ratio , 𝐸 /𝐸 , measuring the elastic anisotropy of the network, systematically decreases with increasing external strain (Fig. 3C), showing that the fibrous networks become increasingly anisotropic as they are stretched. Note that this strain-induced anisotropy occurs in networks composed of all model fibers, even for networks composed of linear fibers. This indicates that strain-induced elastic anisotropy is driven not only by mechanical nonlinearities of individual fibers, but also by geometrical anisotropy such as collective fiber alignment. In addition, we find in Fig. 3C that the data points of 𝐸 /𝐸 for all types of fibers can be approximately fitted by a master curve if plotted versus the normalized strain, 𝜖/𝜖 crit . Here 𝜖 crit is the critical strain, a characteristic parameter of the network, in which 𝐸 /𝐸 becomes smaller than 0.9, and the strain-induced elastic anisotropy is significant. The normalized strain, 𝜖/𝜖 crit , similar to the normalized cell contraction, 𝒞/𝒞 crit in cell contraction networks, is another “emergent” dimensionless parameter for the fibrous network upon external uniaxial strain; it measures the degree of the nonlinearity in the externally-stretched network and determines the strain-induced elastic anisotropy.
Collective fiber alignment measured by nematic order parameter: a characteristic parameter of network elastic anisotropy
Page of We next aim to quantify the elastic anisotropy that develops in the fibrous network due to cell-induced strains, in the circular geometry as depicted in Fig. 1A. For this purpose, we use fiber alignment as a geometric measure of elastic anisotropy. We assume that if two different networks (both isotropic at undeformed state) are composed of the same type of model fibers and have the same degree of fiber alignment or nematic order parameter (NOP), they should be in the same mechanical state; i.e. , they have the same elastic anisotropy. We then approximate the elastic anisotropy of each local area of the network in the cell-contracting system by that of the uniaxially-stretched bulk networks with the same NOP, as measured in Fig. 3C. Thus, we next measure the NOP in both the bulk uniaxial and cell-contraction systems, and use it as a mapping parameter to estimate the elastic anisotropy of the cell system from the measurements of the bulk system. For rectangular networks under uniaxial strain, we measure the orientation of each fiber with respect to the direction of the applied external strain (Fig. 4A). Before stretching, the network is isotropic and the distribution of fiber angles is uniform without any preferred orientation. After stretching, this distribution shifts towards the direction of the applied stretch (Fig. 4B), indicating a clear collective fiber alignment. The degree of fiber alignment can be measured by the NOP in two dimensions (23):
𝑆 = 〈cos(2𝜃)〉 . S ranges from −1 to , where a value of corresponds to an isotropic network, a value of corresponds to a network with fully aligned fibers along the strain direction, and a value of −1 corresponds to a network of transversely aligned fibers. The network comprised of linear fibers shows a linear increase of 𝑆 with strain, while nonlinear networks exhibit stronger alignment with strain (Fig. 4C, inset). To understand this behavior, we use affine approximation to derive the analytical relation between alignment and external strain, and obtain: 𝑆 = (1 + 𝜈)𝜖 , with 𝜈 = −𝜖 /𝜖 being the Poisson ratio (Supporting Information, Section IV). We find that the theoretical prediction of S is in good agreement with the measured 𝑆 from the bulk simulations (Fig. 4C), despite possible local non-affine deformations in our fibrous networks (Supporting Information Section III). Thus, a larger Poisson ratio is directly coupled to enhanced alignment. Indeed, nonlinear fiber models are characterized by a larger Poisson ratio that increase with the external strain, whereas linear fiber networks have a relatively constant Poisson ratio ( 𝜈~ 0.4 ) (Supporting Information, Fig. S6). S therefore depends linearly on the applied strain only for linear networks. To quantify fiber alignment ( i.e., S) in the cell contraction system, we measure the angles relative to the radial direction (Fig. 4D). This is in accordance with the evaluation of alignment described in recent studies (24, 48). 𝑆 is averaged over fibers in annular domains around the cell. Note that although Page of undeformed networks are macroscopically isotropic, they can locally have some weak preferred fiber orientation. We thus measured the degree of alignment in respect to the undeformed state: ∆𝑆 =𝑆 deformed − 𝑆 undeformed for each annular domain around the cell. As shown in Fig. 4E, ∆𝑆 decays with a power-law of ∆𝑆~𝑅 −2 at distances far from the cell, independent of the degree of cell contraction. This is expected, as 𝑆 is linearly proportional to strain for small strains (63), and small radial strain in the far-field regimes decays as −2 (because strains are the first derivatives of displacements). In the near-field regime, we observe that 𝑆 increases linearly with cell contraction and remarkably, it is almost independent of the fiber model type (Fig. 4F), in comparison to the case of uniaxially-stretched networks in rectangular geometry. To understand this better, we derive an analytical prediction for the dependence of 𝑆 on cell contraction 𝒞 in the near-field, using affine displacement approximation with a 𝑈̃ = 𝑅̃ −𝑛 scaling (Supporting Information, Section IV). We obtain the relation: 𝑆 = 𝑄(𝑛, 𝑅̃ ∗ )𝒞 , where 𝑛 is the near-field power-law exponent and 𝑅̃ ∗ is the size of the near-field region normalized by the cell radius. The exact form of 𝑄(𝑛, 𝑅 ∗ ) is given in (Supporting Information, Section IV). We find that the dependence of 𝑄 on 𝑛 and 𝑅 ∗ is rather weak and 𝑄(𝑛, 𝑅 ∗ ) is almost constant for all of our simulations. This analytical relation fits nicely to the simulation results (dotted curve in Fig. 4F(, which suggests that fiber alignment near the cell can be well predicted from affine theory, and has robust trends independent of fiber mechanics. To conclude, we now have a full description of the network alignment in both the uniaxially-stretched bulk and cell simulation networks. Cell-induced elastic anisotropy dictates the slow displacement decay in the fibrous network
As the cell pulls on the matrix, the network alignment in the near-field increases, and the elastic anisotropy is expected to rise. For each magnitude of cell contraction, we estimate the elastic anisotropy, /𝐸 , in the near-field based on the measurements in uniaxially-stretched bulk networks composed of the same model fibers and having the same NOP (Fig. 5A, schematics). We can now plot the displacement transmission power-law exponent 𝑛 against the elastic anisotropy in the near-cell area. We find a linear relation that is similar for all networks types, independent of particular fiber mechanical properties (Fig. 5B). Note that this unified linear relation takes into account the unified master curve for 𝑛 versus normalized cell contraction 𝒞/𝒞 crit (Fig. 2D). Furthermore, the obtained unified linear relation indicates that once the elastic anisotropy is set by
𝒞/𝒞 crit , the displacement decay will be dictated accordingly, independent of the fiber mechanical model. Thus, the cell-induced elastic anisotropy is a
Page of “lumped” parameter, that takes into account the effects of both fiber nonlinearities and fiber alignment, and which governs the slow decay ( 𝑛 < 1 ) of cell-induced displacements in fibrous networks. Controlling the range of displacement transmission by modifying the intrinsic elastic anisotropy
The change of force transmission by elastic anisotropy can be seen more explicitly by looking at the decay of displacements in an intrinsically anisotropic fiber networks with pre-defined tunable elastic anisotropy. An intrinsically anisotropic fiber network can be constructed by introducing an angle-dependent fiber modulus, for example, 𝜇 𝑓 = 𝜇 cos 𝜃 + 𝜇 sin 𝜃 , where 𝜇 𝑓 is the axial fiber modulus, 𝜃 is the fiber orientation (measured relative to the radial direction), 𝜇 and 𝜇 are two variable stiffnesses (Fig. 6A, inset). Note that in the undeformed state, the network is elastically anisotropic but geometrically isotropic, with a uniform fiber angle distribution. In this case, calculations from affine deformation predicts that the ratio between the two principal moduli takes the form of 𝐸 /𝐸 = (𝜌𝜇 + 𝜇 )/(𝜇 +𝜌𝜇 ), where 𝜌 = 5.0 for uniform distribution of fiber orientation [Supporting Information, section V]. We carry out bulk uniaxial simulations to measure 𝐸 /𝐸 for different predefined ratios 𝜇 /𝜇 , and obtain 𝜌 = 4.1 by least-squares fitting (Fig. 6A). Note that in contrast to the strain-induced anisotropy that is determined by the normalized external strain, 𝜖/𝜖 crit , or normalized internal cell contraction, 𝒞/𝒞 crit , the elastic anisotropy in intrinsically anisotropic networks is determined by design based on predefined fiber anisotropy, 𝜇 /𝜇 . We then study the decay of displacements induced by a contractile cell in such intrinsically anisotropic networks. In very small cell contractions, we can consider only the effect of pre-defined elastic anisotropy, while neglecting the effect of fiber alignment. We find that the power-law exponent shows a very good linear proportionality to the elastic anisotropy, 𝑛 ≈ √𝐸 /𝐸 for all values of 𝑛 (Fig. 6B). When the radial direction is stiffer than the angular direction, ( 𝐸 < 𝐸 ) we obtain 𝑛 < 1 , i.e., the displacements decay slowly, in which case, the range of cell-cell communications mediated by the matrix is considerably enhanced. In contrast, when the radial direction is softer than the angular direction ( 𝐸 > 𝐸 ), we obtain 𝑛 > 1 , i.e., the displacement decays faster than in linear isotropic elastic medium. In this case, the range of cell-cell communication is restricted. This indicates that the transmission of cellular forces in fibrous networks and hence, the efficiency of matrix-mediated cell-cell communications, can be reprogrammed by modifying the network anisotropic properties. Page of In summary, we show that the unified mechanism of elastic anisotropy governs the long-range transmission of cell-induced displacement in networks consisting of various types of fibers with different mechanical properties. These findings show that the two effects that were previously considered separately: fiber non-linearity (buckling/stiffening) and fiber reorientation (or alignment), are related and both contribute to the network’s elastic anisotropy, which is a key strain-dependent “lamped” parameter that dictates the decay of the displacements.
The validity of this unified mechanism of elastic anisotropy is shown in two different cases: cells that contracts strongly enough to generate a local anisotropic environment and cells that contracts weakly in an existing pre-designed anisotropic network. In the former case, the elastic anisotropy is induced by cell contraction and is determined by the normalized cell contraction,
𝒞/𝒞 crit . In the latter case, the elastic anisotropy is intrinsic, pre-designed and is determined by the fiber elastic anisotropy, 𝜇 /𝜇 . In both cases, the elastic anisotropy dictates the decay of cell-induced displacements. In this work, we only focused on fibrous networks of high connectivity in two dimensions and neglected the rotational resistance of cross-linkers and the potential feedback of cells to the change in their mechanical microenvironment (3). All these aspects of cell-cell and cell-matrix interactions should be investigated more carefully in the future by extending the continuum theory and finite element simulation presented in this work. Still, this work relates and unifies previous single effects into a general mechanism of elastic anisotropy that explains how (strongly or weakly) contracting cells can maintain long-range mechanical communications with one another. Our findings provide new insights into various biological processes that involve cell-ECM interactions, ECM-mediated cell-cell interactions and for providing a new scaffold design parameter for tissue engineering application.
Numerical Methods
We use Matlab to create the network geometry and architecture, and the finite element software Abaqus/CAE 2017 (
Dassault Systèmes Simulia ) to model the network mechanics and definitions. Truss elements are used to represent the ECM fibers. The response of five networks in different randomnesses is averaged to reduce the effect of the network geometry on the obtained results. The software’s implicit static solver is used to simulate the cell contractions for the four fiber models, up to 50%, and the bulk uniaxial and bi-axial stretching for the two fiber models that exhibit strain stiffening. In the case of high
Page of tensile strains, we use the implicit dynamic solver for the linear and buckling fiber models, with sufficient damping and loading time that is divided into a large number of increments in order to obtain a quasi-static solution. We ensure that the kinetic energy of the whole model declines to a very small fraction of the internal energy at the end of each stretching (less than 0.001%). ACKNOWLEDGMENTS.
We thank Samuel A. Safran from the Weizmann Institute of Science and Yair Shokef from Tel-Aviv University for their useful comments and suggestions. This work was supported by Guangdong Technion–Israel Institute of Technology, by the Israel Science Foundation (1474/16) and Israel Science Foundation- Israeli Centers for Research Excellence (1902/12). 1. Frantz C, Stewart KM, Weaver VM (2010) The extracellular matrix at a glance.
J Cell Sci
Rev Mod Phys
Soft Matter
Proc Natl Acad Sci
Biophys J
Burns
Differentiation
J Cell Biol
Cell
Cell Motil
Proc Natl Acad Sci
Cancer Res
Proc Natl Acad Sci
Calcif Tissue Int
Page of Contractile Cells to Communicate Local Position and Orientation.
PLoS One
Proc Natl Acad Sci
J R Soc Interface
Nature
J Cell Sci
Nat Phys
J Mech Phys Solids
Proc Natl Acad Sci U S A
PLoS One
Biophys J
J Biomech Eng
Biophys J
Biophys J
Phys Rev E
Phys Biol
Sci Rep
6. 31. Ma X, et al. (2013) Fibers in the Extracellular Matrix Enable Long-Range Stress Transmission between Cells.
Biophys J
Nature
New J Phys
Proc Natl Acad Sci
J Phys Chem B
Biophys J
Page of stiffening of cross-linked networks. Biophys J
Soft Matter
Sci Rep
5. 41. Sopher RS, et al. (2018) Nonlinear Elasticity of the ECM Fibers Facilitates Efficient Intercellular Communication.
Biophys J
Phys Rev E - Stat Nonlinear, Soft Matter Phys
Biomech Model Mechanobiol
Phys Rev E
Soft Matter
Soft Matter
Phys Biol
Bull Math Biol
Soft Matter
Nat Phys
Phys Rev Lett
Biophys J
Phys Rev E
Acta Biomater
Contin Mech Thermodyn . 56. Khalily MA, Goktas M, Guler MO (2015) Tuning viscoelastic properties of supramolecular peptide gels via dynamic covalent crosslinking.
Org Biomol Chem
Nat Phys . 58. Head DA, Levine AJ, MacKintosh FC (2003) Deformation of Cross-Linked Semiflexible Polymer Networks.
Phys Rev Lett
Theory of elasticity / by L. D. Landau and E. M. Lifshitz ; translated from the Russian by J. B. Sykes and W. H. Reid ed Lifshitz EM (Evgenii M (Pergamon,
Page of London). linear ela. 60. Shokef Y, Safran SA (2012) Scaling Laws for the Response of Nonlinear Elastic Media with Implications for Cell Mechanics.
Phys Rev Lett
Theory Of Elasticity Of An Anisotropic Body . 62. J. A. A. Gibson L, F. Ashby M (1988)
Cellular Solids: Structure And Properties
63. Fischer-Friedrich E (2018) Active Prestress Leads to an Apparent Stiffening of Cells through Geometrical Effects.
Biophys J
Page of Figure 1: Cell contraction simulation. (A) Illustration of a contractile cell, simulated by a circular cavity, embedded in an isotropic, homogeneous fibrous network. Isotropic inward displacements are imposed at the boundary of the cell, and the outer boundary of the network is held to be fixed (not shown). (B) Force-strain curves for individual fibers used in the simulations. Four different mechanical models are used: linear, stretch-stiffening, compression-buckling, and nonlinear fibers with both buckling and stiffening.
Page of Figure 2: Decay of displacements induced by isotropic cell contraction. (A) Color map of displacements induced by 50% cell contraction. (B) Color map of strains of individual fibers. Radially aligned fibers are stretched (red) while angularly aligned fibers are compressed (blue). (C) Normalized displacement 𝑈̃ as a function of the normalized distance 𝑅̃ from the cell center, for the four fiber models, with 40% cell contraction. Two power-law regimes can be identified. (D) Near-field effective power-law exponent n, versus normalized cell contraction, 𝒞/𝒞 𝑐𝑟𝑖𝑡 . Inset: n decreases linearly with contraction for linear fibers, due to increased fiber alignment. Fiber nonlinearities facilitate the decrease in n. Error bars are averaged over 5 network realizations. (E) Effective cell displacement 𝑈̃ eff versus normalized cell contraction, 𝒞/𝒞 𝑐𝑟𝑖𝑡 . 𝑈̃ eff increased linearly with contraction for all fiber models (inset). Page of Figure 3: Elastic anisotropy induced by external uniaxial tension. (A) Schematic illustration of the rectangular network upon external uniaxial stress. (B) Strain color-map: The fibers along the axis of strain are stretched (red) while transverse fibers are compressed (blue). (C) The ratio of the two principal moduli, 𝐸 /𝐸 , measuring the elastic anisotropy, versus normalized cell contraction, 𝜖 /𝜖 𝑐𝑟𝑖𝑡 . At a given tensile pre-strain 𝜖 , additional infinitesimal tensile strains are applied in both directions to measure the respective elastic moduli, 𝐸 (longitudinal) and 𝐸 (transverse) (see illustration in a). The modulus ratio, 𝐸 /𝐸 , deviates from 1 when the network is strained, and becomes close to zero for high strains. All the three effects – fiber alignment, buckling and stiffening, contribute to this strain-induced anisotropy. Inset shows 𝐸 /𝐸 versus cell contraction (without normalization). Page of Figure 4: Strain-induced collective alignment of fibers. (A) Schematic illustration of uniaxially-stretched network. The orientational angle, 𝜃, of each fiber is defined with respect to the strain axis. (B) Histogram of the distribution of individual fiber angles. (top) Fiber angle distribution in an undeformed isotropic network. (bottom) Fiber angle distribution for a linear fiber network upon 0.52 uniaxial strain; collective fiber alignment along the strain direction is indicated with more fiber fraction towards zero angle. (C) Nematic order parameter S of the stretched network as a function of the applied uniaxial strain 𝜖 . All fiber models follow the affine prediction 𝑆 = 1/2(1 + 𝜈)𝜖 . (D) Schematic illustration of a cell-contracted network. The orientational angle, 𝜃, of each fiber is defined with respect to the radial direction. (E) Nematic order parameter S of the contracted network (of linear fibers) as a function of 𝑅̃ . S decays with ~𝑅 −2 in the far field with a slower decay in the vicinity of the cell. 𝑆 increases with cell contraction. (F) Averaged S in the near-field region as a function of cell contraction. S increases linearly with cell contraction, and fits well to the affine theory (dashed black line). Error bars were calculated based on the average of 5 network realizations. Page of Figure 5: Elastic anisotropy dictates the decay of displacements in fibrous networks. A Schematic illustration of the transformation from bulk to cell contraction simulations: If two different networks (both isotropic at undeformed state) are composed of similar fibers have the same nematic order parameter, they are assumed to be in the same mechanical state, i.e., have the same elastic anisotropy, 𝐸 /𝐸 . (B) The near-field power- law exponent n for the decay of the displacement is plotted as a function of the square root of the network anisotropy √𝐸 /𝐸 . We found 𝑛 ∝ √𝐸 /𝐸 which is consistent with the theoretical prediction for an anisotropic continuum material. Page of Figure 6: The decay of cell-induced displacements in intrinsically anisotropic networks. (A) An intrinsically anisotropic network is constructed by introducing an orientation-dependent fiber modulus 𝜇 𝑓 =𝜇 cos (𝜃) + 𝜇 sin (𝜃) (inset). The elastic anisotropy of the network is predicted by affine theory, as 𝐸 /𝐸 = (𝜌𝜇 + 𝜇 )/(𝜌𝜇 + 𝜇 ) (dashed line) with 𝜌 ≈ 4.1 (black dots, simulation data) which are closed to the theory prediction 𝜌 = 5.0 . (B) The near-field power-law exponent, 𝑛, is plotted as a function of network anisotropy at infinitesimal cell contractions; a very good linear fitting to 𝑛 = √𝐸 /𝐸 (dashed line) is obtained in both regions: 𝐸 > 𝐸 (with n>1 indicating fast displacement decay) and 𝐸 > 𝐸 (with n<1 indicating slow displacement decay). Page of Supporting Information
I. Mechanical properties of four model fibers .............................................................................. 22
II. Numerical construction of an isotropic and homogeneous fibrous network ............................ 23
III. Nonlinear elastic responses of fibrous networks to external uniaxial stretch ......................... 24
1. Nearly-affine deformations of the network ........................................................................... 24
2. Stretch-induced anisotropy and Poisson ratio of the network .............................................. 25
IV. Strain-induced collective fiber alignment ............................................................................... 27
1. Fiber alignment induced by external uniaxial stretch ........................................................... 27
2. Fiber alignment induced by internal cell contraction ........................................................... 28
V. Anisotropic elasticity of an intrinsically anisotropic network ................................................. 29
VI. Continuum theory: decay of displacements induced by cell contraction ............................... 30
1. Linear anisotropic elastic continuum .................................................................................... 30
2. Strain-induced anisotropic nonlinear elastic continuum ....................................................... 31 In this supporting information, we provide all the details about the method of discrete fiber simulations of fibrous networks, quantification of elastic anisotropy and fiber alignment, and the simplified continuum model we are using. I. Mechanical properties of four model fibers
Each individual fiber is modeled as a one-dimensional linear truss element undergoing uniaxial tension or compression. The four material models used to simulate the mechanical elastic behavior of the individual fibers ( Fig. 1A in the main text) are modeled as follows: linear fibers have an elastic modulus of , both tensile and compressive. For buckling fibers, the elastic modulus is ten times smaller in compressive strains larger than 2%. Stiffening fibers in tension are characterized by elastic modulus, which increases exponentially in tensile strains above 2% and is constant in lower strains. Buckling-stiffening fibers combine both softening and stiffening of the two last mentioned material models,
Page of respectively. The elastic modulus assigned to each of the individual fibers exhibit buckling and stiffening is given by 𝐸 = { 𝜌𝐸 ∗ , 𝜀 < 𝜀 𝑏 𝐸 ∗ , 𝜀 𝑏 < 𝜀 < 𝜀 𝑠 𝐸 ∗ 𝑒𝑥𝑝[(𝜀 − 𝜀 𝑠 )/𝜀 ], 𝜀 > 𝜀 𝑠 (S1) while the values used in the model are 𝐸 ∗ = 11.5𝑘𝑃𝑎 as the given elastic modulus, 𝜀 𝑏 = 2% is the minimal compressive strain for which buckling occurs, 𝜀 𝑠 = 2% is the tensile strain above which strain stiffening occurs, 𝜀 = 5% is the strain-stiffening coefficient and 𝜌 = 0.1 is the buckling ratio, which is the ratio between the compressive and tensile elastic moduli. This mechanical behavior is largely based on computational models described in previous studies (1). II. Numerical construction of an isotropic and homogeneous fibrous network
First, we randomly insert some nodes in two dimensions throughout a circular domain of radius 𝑅 b by following a uniform distribution. Pairs of nodes are then connected by single fibers generated following a minimum cost algorithm as follows. The probability that a potential fiber connecting two neighboring nodes would be generated is determined by a cost function, P . In every iteration, each node looks at the thirty nearest neighbors in its vicinity as potential fibers to be generated. The cost function associated with each potential fiber is given by: 𝑃 = 𝑁 𝑖𝑗 + 𝑎𝐶 𝑖𝑗 + 𝑏(𝐴 𝑖 + 𝐴 𝑗 ) (S2) with the two constants 𝑎 and 𝑏 determining the relative importance of each term. A potential fiber can be generated only if its cost function, P, is negative. The three terms in Eq. (S2) are explained as follows. (i) 𝑁 𝑖𝑗 , is the degree of nearness of each node i to the other node j . For example, if node 𝑖 is the second nearest neighbor of node 𝑗 , and node 𝑗 is the third nearest neighbor of node 𝑖 , then 𝑁 𝑖𝑗 = 2 + 3 = 5 . (ii) 𝐶 𝑖𝑗 = 𝑐 𝑖 + 𝑐 𝑗 − 2𝑐 opt takes into account of the current connectivity numbers 𝑐 𝑖 and 𝑐 𝑗 of nodes 𝑖 and 𝑗 , respectively, with 𝑐 opt being a chosen optimal connectivity. We can set the Fig. S1: Schematic illustration of local geometrical consideration in generating new fibers. For a new candidate fiber between nodes, i, and j, we define 𝐴 𝑖 = 𝛽 𝑖 −𝛼 𝑖 , where 𝛼 𝑖 and 𝛽 𝑖 are the maximum free angles before and after adding the new fiber. If 𝐴 𝑖 is negative, the fiber will have a lower cost. The same is done for node j. Page of average connectivity of the network by tuning 𝑐 opt . The higher 𝑐 opt is, the smaller 𝐶 𝑖𝑗 becomes and so the probability for creating the fiber gets higher. (iii) The last term relates to the maximum free angle of each node (as shown in Fig. S1), representing the local geometry around the nodes. Here 𝐴 𝑖 is the difference in maximum free angle after and before adding the new potential fiber to node i . A negative 𝐴 𝑖 increases the probability that a fiber would be generated. In summary, in every iteration, each node is connected by a new fiber to the node with the minimal (negative) cost (in comparison to the other twenty-nine nearest neighbors). The iteration process finishes when no new fiber can be added. The major advantages of this method of network construction are as follows. (i) The generated networks are therefore always isotropic. Nodes tend to be connected by fibers that are oriented in a wide range of directions, in comparison to other methods where many fibers that are aligned mostly in some particular directions. This increases the isotropy and homogeneity of the networks. (ii) It allows for control of geometrical features like pore size and connectivity in an elegant way. (iii) It does not depend on macroscopic length scales such as cell and system sizes, but only depends on the density of nodes. In addition, for the simulations of networks that are contracted by an embedded cell, we remove all fibers or fiber sections inside the cell area, not affecting the fibers outside. The cell edge remains linked with the surrounding network. This guarantees that the procedure of embedding the cell in the network does not change the geometry around it. For bulk simulations with a network stretched by uniaxial stresses, a rectangular section of network is chosen and the rest of the network is removed. III. Nonlinear elastic responses of fibrous networks to external uniaxial stretch
1. Nearly-affine deformations of the network
Our networks are classified as high-connectivity, stretch dominated networks, and are therefore expected to deform in an affine manner. For affine deformations, the coordinated of the nodes are transformed by 𝑥 → 𝑥(1 + 𝜖 ); 𝑦 → 𝑦(1 + 𝜖 ) , where 𝜖 , 𝜖 are the macroscopic principal strains. From this it can be proved that the strain and orientation of each fiber following deformation, depend only on its initial orientation and the macroscopic network principal strains, and are given by: 𝜖 𝑓 = √(1 + 𝜖 ) cos (𝜃 ) + (1 + 𝜖 ) sin 𝜃 − 1 (S3) cos (θ) = (1 + ϵ ) cos (𝜃 )(1 + ϵ ) cos (𝜃 ) + (1 + ϵ ) sin (𝜃 ) (S4) Page of where 𝜖 𝑓 and 𝜃 are the strain and final orientation of the fiber, 𝜃 is the initial fiber orientation. The following figures compare our simulation results to these predictions. The figures are shown only for linear fibers and stiffening-buckling fibers. Figures S2-S4 indicate that our network deforms in a way that is close to affine over all strains and fiber models. Large strains and fiber nonlinearity increase the local deviations from the affine predictions, but the overall affine trend is roughly preserved. Notice that in general, the average strains in the network are smaller than the affine predictions, leading to lower deformation energies. Notice also the large Poisson effect for nonlinear fibers.
2. Stretch-induced anisotropy and Poisson ratio of the network
The initially (before deformation) isotropic and homogeneous fibrous networks show highly nonlinear elastic responses to the applied external uniaxial stretch. On the one hand, the Poisson ratio, ν , depends significantly on the applied uniaxial strain, 𝜖 (See Fig. S5), except for the networks of linear fibers where ν ≈ 0.4 . ν increases from its linear value ~0.4 to be even larger than 1.6, beyond the upper limit, 1.0, of Poisson ratio in 2D linear isotropic, elastic continuum. Note that the volume change is described by strain 𝜖 + 𝜖 = 𝜖 (1 − 𝜈 ) ; that is the network volume increases if 𝜈 <1 and decreases if Poisson ratio is greater than 1 for the networks of stiffening fiber, i.e., 𝜈 >1. This indicates the presence of nonlinearity as well as fiber dilution/ densification in the deformed fibrous networks. Fig. S2: Comparison of affine predictions (Eqs. (S3) and (S4)) to simulations, for Linear fibers under small strains. Fig. S3: Comparison of affine predictions (Eqs. (S3) and (S4)) to simulations, for Linear fibers under large strains. Fig. S4: Comparison of affine predictions (Eqs. (S3) and (S4)) to simulations, for Stiffening-Buckling fibers under large strains.
Page of On the other hand, the initially isotropic network will become anisotropic in elasticity when it is stretched strongly enough. Fig. S6 shows that the two principal moduli, 𝐸 and 𝐸 , become unequal over some threshold strain, representing elasticity anisotropy induced by the applied strain, 𝜖 . Note that the longitudinal modulus, 𝐸 , increases with strain in all the cases. The curve of linear fibers overlaps that of buckling fibers, while the curve of stiffening fibers overlaps that of stiffening-buckling fibers; the slope of later curve is much larger than the former one due to stretch-stiffening nonlinearity. In comparison, the transverse modulus, 𝐸 , decreases with strain in all the cases except for stiffening fibers. The curve of buckling fibers overlaps that of stiffening-buckling fibers, the slope magnitude of which is much larger than that of linear fibers due to the compression-buckling nonlinearity in the transverse directions. Interestingly, the transverse modulus, 𝐸 , of stiffening fibers increases (not decreases as the other fiber types) with strain, 𝜖 . However, when we plot 𝐸 /𝐸 as a function of the applied strain, 𝜖 (as shown in Fig. 3 in the main Fig. S5. Two principal moduli, 𝐸 and 𝐸 , along the longitudinal and transverse directions, respectively, as a function of the applied strain, 𝜖 . Elastic anisotropy, represented by the inequality of 𝐸 and 𝐸 , is shown to be induced by the applied external uniaxial strain. Fig. S6: Poisson ratio, 𝜈 , as a function of applied uniaxial strain, 𝜖 for fibrous networks composed of four different model fibers. Page of text), we find 𝐸 /𝐸 is always smaller than 1.0 for all the fiber types, although the magnitudes of 𝐸 /𝐸 for given strain depend significantly on the particular fiber type. IV. Strain-induced collective fiber alignment
1. Fiber alignment induced by external uniaxial stretch
The external uniaxial stress applied on the fibrous networks will induce anisotropy not only in elasticity but also in geometry -- nematic order accompanying collective fiber alignment. To quantify the collective uniaxial fiber alignment in 2D networks upon external uniaxial stretch, we employ the 2D nematic order parameter:
𝑆 ≡ 〈cos(2𝜃)〉 = ∫ 𝑝(𝜃)cos(2𝜃) 𝑑𝜃 (S5) where 𝜃 ∈ [0, 𝜋) is the fiber orientation with respect to the direction of the applied stretch and 𝑝(𝜃) is the probability distribution function of the fiber orientation. For initially (before stretched) isotropic networks, we assume the fiber angle distribution to be uniform, i.e., 𝑝(𝜃) = 1/𝜋 , and the corresponding nematic order parameter is 𝑆 = 0 . For deformed fibrous networks with two principal strains, 𝜖 and 𝜖 , if we assume affine deformations of each fibers, we then obtain the new angle , 𝜃, of the fiber with initial angle, 𝜃 , is given by 𝜃 = atan (tan(𝜃 ) 1 + 𝜖 ), (S6) which in the case of small strains ( 𝜖 , 𝜖 ≪ 1) can be expanded to the first order in strains as 𝜃 ≈ 𝜃 + 𝜖 − 𝜖 ) . (S7) Moreover, we can get the probability distribution of fiber orientation in deformed networks as 𝑝(𝜃) = 𝑝(𝜃 ) 𝑑𝜃 𝑑𝜃 ≈ [1 + (𝜖 − 𝜖 ) cos(2𝜃)], (S8) from which we obtain 𝑆 = 1𝜋 ∫(𝜖 − 𝜖 ) cos (2𝜃) 𝑑𝜃 = 𝜖 − 𝜖 Particularly, in the case of uniaxial deformation with 𝜖 = −𝜈𝜖 , we obtain 𝑆 = (1 + 𝜈)2 𝜖 . (S10) Note that Poisson ratio is usually not a constant but a nonlinear function of strain 𝜖 as shown in Fig. S5. Therefore, the nematic order parameter induced by strain is generally a nonlinear function of strain as Page of shown in Fig. 4C (in the main text) and is linearly proportional to strain only if Poisson ratio is a constant (also see Fig. 4C). We have compared the theoretical prediction Eq. (S10) with simulation results; very good agreement is achieved as shown in Fig. 4C. There, we plot as a function of 𝜖 and find that all the data obtained from four types of fiber networks fall onto one single master curve as indicated by Eq. (S10).
2. Fiber alignment induced by internal cell contraction
Let’s now consider the collective fiber alignment in fibrous networks that is induced by internal cell contraction. In the near-field region, the decay of the cell-induced displacement follows a power law 𝑢 = 𝑈 cell ( 𝑅𝑅 cell ) −𝑛 , (S11) with 𝑅 cell and 𝑈 cell < 0 being the cell radius and its contracting displacement, respectively, from which we obtain the two principal strains as 𝜖 = 𝜕𝑢𝜕𝑅 = −𝑛𝑈 cell 𝑅 −𝑛−1 𝑅 cell−𝑛 , 𝜖 = 𝑢𝑅 = 𝑈 cell 𝑅 −𝑛−1 𝑅 cell−𝑛 . (S12) In this case, we obtain from Eq. (S9) that the nematic order parameter, S , of small section of the network at 𝑅 , is given by 𝑆(𝑅) = − 1 + 𝑛2 𝑈 cell 𝑅 −𝑛−1 𝑅 cell−𝑛 . (S13) Averaging S over the near-field region from the cell boundary 𝑅 = 𝑅 cell to some characteristic length
𝑅 =𝑅 ∗ , we obtain the averaged near-field S as 〈𝑆〉 = ∫ 𝑆(𝑅) 𝑅 ∗ 𝑅 cell 𝑅𝑑𝑅∫ 𝑅 𝑅 ∗ 𝑅 cell 𝑑𝑅 = 𝑄(𝑛, 𝑅 ∗ ) 𝑈 cell 𝑅 cell , (S14) with 𝑄(𝑛, 𝑅 ∗ ) ≡ − 1 + 𝑛1 − 𝑛 𝑅̃ ∗1−𝑛 − 1𝑅̃ ∗2 − 1 , (S15) and 𝑅̃ ∗ ≡ 𝑅 ∗ /𝑅 cell . We have also compared the theoretical prediction (S14) with simulation results. All the data obtained from four different types of fiber networks fall on one curve which agrees with (S14) very well as shown in Fig. 4F in the main text. Page of V. Anisotropic elasticity of an intrinsically anisotropic network
We construct a homogeneous fibrous network composed of linear fibers with intrinsic anisotropy by following the same method as explained in Sec. I with uniform fiber angle distribution, i.e., 𝑝(𝜃) = 1/𝜋 for 𝜃 ∈ [0, 𝜋) , and we introduce an orientation-dependent stiffness, 𝜇 𝑓 , to each fiber as 𝜇 𝑓 = 𝜇 cos (𝜃) + 𝜇 sin (𝜃), (S16) where 𝜃 is the angle of the fiber with respect to radial direction of the cell, 𝜇 and 𝜇 being the two extrema of 𝜇 𝑓 along the radial (with 𝜃 = 0 ) and angular (with 𝜃 = 𝜋/2 ) directions, respectively. That is, we have constructed a fibrous network with intrinsic anisotropy in elasticity but not in geometry (without collective fiber alignment). If the fibers in the network are deformed affinely, the strain of a fiber with orientation 𝜃 is: 𝜖 𝑓 = √(1 + 𝜖 ) cos (𝜃) + (1 + 𝜖 ) sin (𝜃) − 1, (S17) where 𝜖 , 𝜖 are the two principal strains along the radial and angular directions, respectively. Then the corresponding deformation energy of the fiber is 𝐹 𝑓 (𝜃, 𝜖 , 𝜖 ) = 𝑙 𝑓 𝜖 𝑓2 = 𝑙 cos (𝜃) + 𝜇 sin (𝜃)) (√(1 + 𝜖 ) cos (𝜃) + (1 + 𝜖 ) sin (𝜃) − 1) , (S18) with 𝑙 being the undeformed length of the fiber. For small strains, 𝐹 𝑓 (𝜃, 𝜖 , 𝜖 ) can be expanded in 𝜖 , 𝜖 to the second order as 𝐹 𝑓 ≈ 𝑙 cos 𝜃 + 𝜇 sin 𝜃) [1 + (1 + 𝜖 ) cos 𝜃 + (1 + 𝜖 ) sin 𝜃− 2 (1 + 𝜖 cos 𝜃 + 𝜖 sin 𝜃 + 𝜖 𝜃 − cos 𝜃) + 𝜖 𝜃 − sin 𝜃)− 𝜖 𝜖 cos 𝜃sin 𝜃)]. (S19) For very small deformation, fiber angle distribution stays to be more or less uniform; we can integrate 𝐹 𝑓 (𝜃) over fiber angle to get the total deformation energy as a function of 𝜖 , 𝜖 : 𝐹(𝜖 , 𝜖 ) ≈ 𝑁 𝑓 𝜋 ∫ 𝐹 𝑓 (𝜃)𝑑𝜃= 𝑁 𝑓 𝑙
4𝜋 [(5𝜋8 𝜖 + 𝜋8 𝜖 + 𝜋2 𝜖 𝜖 ) 𝜇 + (𝜋8 𝜖 + 5𝜋8 𝜖 + 𝜋2 𝜖 𝜖 ) 𝜇 ], (S20) Page of with 𝑁 𝑓 is the total number of fibers. From 𝐹(𝜖 , 𝜖 ), we can calculate the two principal (linear) elastic moduli: 𝐸 = 𝜕 𝐹𝜕𝜖 = 𝑁 𝑓 𝑙
16 (5𝜇 + 𝜇 ), (S21a) 𝐸 = 𝜕 𝐹𝜕𝜖 = 𝑁 𝑓 𝑙
16 (𝜇 + 5𝜇 ), (S21b) from which we obtain 𝐸 𝐸 = (5𝜇 + 𝜇 )(𝜇 + 5𝜇 ). (S22) We fitted the simulation results using the following expression 𝐸 𝐸 = (𝜌𝜇 + 𝜇 )(𝜇 + 𝜌𝜇 ), (S23) with 𝜌 ≈ 4.1 , which is very close to the prediction 𝜌 = 5.0 in Eq. (S22) from affine theory. VI. Continuum theory: decay of displacements induced by cell contraction
1. Linear anisotropic elastic continuum
In this section, we calculate the decay of displacements induced by a contracting cell in linear anisotropic elastic medium in two dimensions. To be specific, we consider a circularly contracting cell in 2D linear spherically isotropic elastic medium, in which there are three independent material parameters, two principal Young’s moduli, 𝐸 and 𝐸 (along its radial and tangential axes respectively), and one Poisson’s ratio 𝜈 or 𝜈 with 𝜈 /𝐸 = 𝜈 /𝐸 . The corresponding elastic energy density is given by 𝐹 = 12 𝐸 𝜈 𝜖 + 𝐸 𝜈 𝜈 𝜖 𝜖 + 12 𝐸 𝜈 𝜖 , (S24) from which we obtain the stress-strain relations as 𝜎 𝑟 = 𝐸 𝜈 𝜖 𝑟 + 𝐸 𝜈 𝜈 𝜖 𝜃 𝜎 𝜃 = 𝐸 𝜈 𝜈 𝜖 𝑟 + 𝐸 𝜈 𝜖 𝜃 , (S25) or equivalently, 𝜖 𝑟 = 1𝐸 𝜎 𝑟 − 𝜈 𝐸 𝜎 𝜃 𝜖 𝜃 = − 𝜈 𝐸 𝜎 𝑟 + 1𝐸 𝜎 𝜃 . (S26) Page of For the displacement fields generated by a circularly contracting cell in 2D infinite elastic medium, we have circular symmetry and the two principal strains are given by 𝜖 𝑟 = 𝑑𝑢𝑑𝑅 ; 𝜖 𝜃 = 𝑢𝑅 , (S27) where 𝑢 is the displacement along the radial direction. The 2D mechanical equilibrium condition is given by 𝑑𝜎 𝑟 𝑑𝑅 + 𝜎 𝑟 − 𝜎 𝜃 𝑅 = 0. (S28)
Substituting Eq. (S27) into Eqs. (S25) and (S28) subsequently, we obtain the equilibrium equation for displacement field 𝑢(𝑅) as 𝑑 𝑢𝑑𝑅 + 1𝑅 𝑑𝑢𝑑𝑅 − 𝐸 𝐸 𝑢𝑅 = 0. (S29) The general solution of Eq. (S29) is 𝑢 = 𝑎𝑅 𝑛 + 𝑏𝑅 −𝑛 , (S30) with 𝑛 = √𝐸 /𝐸 . (S31) Particularly, for cells in an infinite elastic medium with boundary conditions, lim 𝑅→ ∞ 𝑢 = finite and 𝑢(𝑅 = 𝑅 cell ) = 𝑈 cell , we obtain Eq. (S11) (see Sec. IV).
2. Strain-induced anisotropic nonlinear elastic continuum
More generally, if the fibrous network that is contracted by a circular cell can be simply mapped to be a spherically isotropic elastic medium as 𝜎 𝑟 = 𝐸̃ 𝜖 𝑟 + 𝐸̃ 𝜈 𝜖 𝜃 𝜎 𝜃 = 𝐸̃ 𝜈 𝜖 𝑟 + 𝐸̃ 𝜖 𝜃 , (S32) where 𝐸̃ , 𝐸̃ , and 𝜈 are the mapping parameters that usually also depend on strains. Substituting it into Eq. (S28), we can obtain a similar equation as (S29) for the highly nonlinear near-field region, in which we can assume that the radial moduli 𝐸̃ and Poisson’s ratio 𝜈 changes very slowly and increases with cell contractions. In this case, the decay of near-field displacement follows an effective power law 𝑢̃ ~1/𝑅̃ 𝑛 with exponent, 𝑛 , determined by relative magnitude between two principal moduli. This is consistent with the simulation data showing a good power-law decay of cell-induced displacements in the near-field region (see Fig. 2B in the main text). We match the far-field solution and the near-field solution given by Page of 𝑢̃ far = 𝑈̃ eff /𝑅̃, (S33) 𝑢̃ near = 𝑎𝑅̃ 𝑛 + 𝑏𝑅̃ −𝑛 , (S34) respectively, using the continuity of the displacement and stress at 𝑅 = 𝑅 𝑇 . Here 𝑈̃ eff is the normalized effective far-field displacement, 𝑈̃ eff /𝑅̃ 𝑇2 = 𝒜 −1 with 𝒜 = 𝒞/𝒞 cr being the normalized cell contraction and measuring the nonlinearity of the cell-contracted network; we obtain 𝑎 = 1 − 𝑏 and 𝑏 =𝒜 −1 𝑛+12𝑛 𝑅̃ 𝑇𝑛+1 with the matching radius 𝑅̃ 𝑇 satisfying the equation 𝑇1+𝑛 − 1 − n2n 𝑅̃
𝑇1−𝑛 − 𝒜 = 0. (S35)