Pulling-force generation by ensembles of polymerizing actin filaments

The process by which actin polymerization generates pulling forces in cellular processes such as endocytosis is less well understood than pushing-force generation. To clarify the basic mechanisms of pulling-force generation, we perform stochastic polymerization simulations for a square array of polymerizing semiflexible actin filaments, having different interactions with the membrane. The filaments near the array center have a strong attractive component. Filament bending and actin-network elasticity are treated explicitly. We find that the outer filaments push on the membrane and the inner filaments pull, with a net balance of forces. The total calculated pulling force is maximized when the central filaments have a very deep potential well, and the outer filaments have no well. The steady-state force is unaffected by the gel rigidity, but equilibration takes longer for softer gels. The force distributions are flat over the pulling and pushing regions. Actin polymerization is enhanced by softening the gel or reducing the filament binding to the membrane. Filament-membrane detachment can occur for softer gels, even if the total binding energy of the filaments to the membrane is 100 or more. It propagates via a stress-concentration mechanism similar to that of a brittle crack in a solid, and the breaking stress is determined by a criterion similar to that of the ‘Griffith’ theory of crack propagation.


Introduction
In many cellular processes that require large forces to generate membrane curvature, such as formation of protrusions, endocytosis, and phagocytosis, actin is an essential factor [1]. Bending the membrane requires pushing and pulling forces in balance. Generation of pushing forces by actin polymerization has received substantial quantitative experimental study. For example, polymerization of individual actin filaments [2], and small numbers of actin filaments [3] yields forces on the order of pN. Force measurements on growing branched actin networks in vitro using cantilevers [4,5] found pushing pressures of ≈ 0.2-1.0 × 10 −3 pN nm −2 . Force densities on actinpropelled biomimetic beads [6] reach values of at least 2.5 × 10 −4 pN nm −2 . On the other hand, the processes by which actin polymerization generates pulling forces have received less quantitative study. Experiments on motile fluid vesicles propelled by actin comet tails [7,8] found a force distribution dominated by inward pushing forces on the sides of the vesicle, and directional pulling forces at the rear of the vesicle. The maximum pulling pressure in [7] was about ≈3.5 × 10 −4 pN nm −2 . Measurements of force around podosomes, mechanosensitive adhesion cell structures that exert protrusive forces onto extracellular environment, show that pushing forces from actin polymerization at the core and pulling forces from lateral acto-myosin contractility in the surrounding adhesion ring are required for a single podosome to deform the substrate [9].
Endocytosis in yeast also requires pulling forces from actin. The required magnitude is large, since overcoming the turgor pressure of 0.2 pN nm −2 or more [10] requires a comparable pulling-force density from actin polymerization. Actin patches consisting of an Arp2/3-branched network [11] form during this process. This network constitutes a crosslinked gel whose mechanical properties are not well known. The pulling forces required to initiate invagination are generated only after the arrival of actin [12], suggesting that the network generates them. Furthermore, recent superresolution microscopy studies of the geometry of the process demonstrated accumulation of the membrane-filament binding protein Sla2, within a central dot, surrounded by a ring of the actin nucleator Las17 [13] (a WASP homolog). This suggested a generic mode of pulling-force generation with enhanced actin polymerization in a ring-shaped region creating pulling forces at the center. Deleting the yeast crosslinking protein Sac6, which should reduce the stiffness of the actin gel [14], stops invagination [15]. Reducing the turgor pressure by providing osmotic support across the plasma membrane reduces the requirement for actin filament cross linkers [16], presumably because the force requirement is lowered. These observations, in combination, show that a stiff actin gel is required for robust pulling-force generation. In addition, mutating Sla2 by deleting its actinbinding domain stops the invagination process [15], showing that strong actin-membrane attachments are crucial. Further contrib utions are probably generated by curvature-generating proteins such as clathrin. However, these are not sufficient to drive the process, as shown by the correlated electron-microscopy and light-microscopy studies of [12]. This work showed that no measurable membrane bending occurs without polymerized actin, suggesting that actin polymerization is the dominant factor generating pulling forces.
Although numerous theoretical models have described how actin polymerization generates pushing force [17,18], generation of pulling forces has been studied less extensively. Most studies have assumed that the total force exerted by the actin network on the membrane vanishes, corresponding to overall force balance on the actin network. This assumption is justified by the smallness of the viscous and inertial forces acting on the actin network [19]. Force balance implies that there are two types of filaments, pushing and pulling. Simple calculations based on the surface area of the invagination and the turgor pressure indicated that forces of over 1000 pN were required to drive endocytic invagination in yeast [20]. Subsequent analyses performed by fitting to observed membrane shapes, including force terms from membrane tension, membrane curvature, and curvature-generating proteins, gave estimates of ∼3000 pN for the minimum required actin pulling force [21,22].
A few models have explicitly treated the process by which actin polymerization generates pulling forces during endocytosis in yeast. They have assumed enhanced actin polymerization in a ring-shaped region. Using an actin polymerization rate increasing continuously outwards from the center of the endocytic patch, the continuum-mechanics calculations of [19] showed that even modest actin polymerization forces spread over a large ring can generate a large pulling-force density at the center by a force-amplification process. These calculations were extended in [22] to evaluate the actin growth profile needed to generate the required pulling forces. Subsequent work [23] treating the actin network as a visco-active gel found that it could exert sufficient pulling forces to drive the process. The validity of a continuous deterministic treatment of the discrete stochastic system of filaments and membrane, however, is not clear.
Stochastic simulations of the growth of a rigid 3D actin network [24] during endocytosis in yeast assumed a uniform distribution of the actin nucleator Las17 in a ring-shaped region on the membrane. They also assumed that the filament growth velocity is determined by the average opposing force per filament. However, the distribution of the forces exerted on the membrane was not obtained. Simulations based on realistic dimensions and molecular compositions taken from superresolution experiments found that a 3D branched network of actin filaments can produce forces exceeding 1000 pN, enough to overcome turgor pressure [13]. However, in this model, the actin filament stall forces may have been overestimated [25], which could lead to an overestimate of the pulling force.
Although these models have confirmed the ability of an actin-nucleator ring to generate pulling forces, there has been no systematic study of the mechanisms determining the magnitude of the pulling force, and how this magnitude is affected by key physical properties. By force balance, the pulling force must be limited by the total stall force of the pushing filaments. However, it is neither clear what fraction of this limit can be achieved practically, nor how rapidly the pulling force reaches its maximum value. In addition, there have been no detailed studies of the spatial distribution of the pulling force, explicitly treating stochastic polymerization of individual filaments. Finally, the possibility of pulling filaments detaching from the membrane has not been treated in detail. This process is plausible because of the large magnitude of the pulling force per pulling filament.
The key features that affect pushing-force generation, such as the free-monomer concentration C A and the on/off rate constants k 0 on and k 0 of f , are important for pulling-force generation as well. The single-filament stall force obtained by thermodynamic arguments [26] is where k 0 on and k 0 of f are rates for a free filament not interacting with an obstacle, C A is the actin monomer concentration, δ is the actin step size per added subunit, and θ is the angle of incidence relative to the direction of motion. Thermodynamic analysis shows that the stall force is proportional to the number of filaments and is independent of the geometrical details of the growth process, a result confirmed by simulations [27]. However, it has been suggested that ATP hydrolysis can reduce the stall force [28]. Lateral interactions [29,30] can also affect the stall force. But in branched actin networks such those at endocytic actin patches, the filament spacings of ∼10 nm are large enough that lateral interactions are probably not important. Bending of individual actin filaments, and the elasticity of the actin network on whole, are also likely important for pulling-force generation. Previous studies have shown that fluctuations of the tips or bases of actin filaments are crucial for obtaining rapid network growth under moderate forces [31,32].
Here we calculate the pulling forces that a discrete array of cross-linked actin filaments with varying polymerization properties exerts on a rigid obstacle. We evaluate the total magnitude of the pulling force, its spatial distribution, the dynamics of the force buildup, and the conditions that lead to detachment of the pulling actin filaments from the membrane. The scale and parameters of the model are chosen to correspond to endocytosis in yeast. We assume that there is a group of filaments at the center of the array whose zero-force growth velocity is low because their growing ends are strongly bound to the obstacle. When the remaining 'pushing' filaments grow, the slower-growing filaments inhibit this growth via their indirect linkage to the pushing filaments mediated by the obstacle, and thus exert a pulling force.
We vary several parameters, including the pulling filaments' binding strength and the actin gel stiffness, and evaluate the resulting effects on force generation and actin network deformation. We interpret our numerical results from this complete system of filaments, including 'pushers' and 'pullers', via a mean-force model that treats the pushers and pullers separately. In this model, the growth velocity is determined by the equality of the pusher and puller forces at a common velocity. We calculate the effect of transient attachments of the pushing filaments to the membrane on the magnitude of the pulling force. We evaluate the dependence of the time scale of force generation rate on the gel stiffness. Finally, we investigate how detachment of pullers depends on the gel stiffness.
We find that strengthening the filament-obstacle binding by choosing deeper potential wells for the central filaments decreases the growth rate of the actin network in this region. This increases the total pulling force up to a maximum that becomes the sum of the stall forces of the surrounding pushing filaments when the central filaments do not polymerize at all. However, even very slow polymerization of the central filaments can strongly reduce the pulling force. We also find that the mean-force model accurately predicts the results of full system simulations for the total pulling force. We find that the time required for the maximum pulling force to build up is roughly inversely proportional to the gel stiffness. Transient attachments of the pushing filaments reduce the total pulling force, so maximum force is produced when their potential is purely repulsive. Finally, softening the gel, or weakening the binding of the central filaments, can lead to actin gel detachment from the membrane despite the total binding energy of several hundred k B T or more between the puller filaments and the membrane.

Model
We model the growing actin network at an endocytic site as a 12 × 12 square array of filaments with spacing a, interacting with a flat moving obstacle (see figure 1). The geometry is motivated by the measured architecture [13] of the endocytic actin patch, but to reduce the computational load the 120 × 120 nm square geometry is taken somewhat smaller than the measured circular geometry of radius 100 nm. We treat a fixed number of uncapped filaments rather than treating the dynamics by which filaments are nucleated by Arp2/3 complex and subsequently capped. Thus we model the force-generation properties of the filaments that are uncapped at a given time. The obstacle, corresponding to the combination of the cell wall and membrane, contains a central 6 × 6 patch of a filamentmembrane binding protein such as Sla2 (yellow circles), surrounded by a square band of nucleation promoting factors (NPFs), such as the yeast WASP homolog Las17 (purple circles). As the filaments grow, more rapid growth of the filaments in the outer region causes pulling and pushing forces to act on the gel in the directions indicated by black arrows, deforming the gel. In order to elucidate the physical mechanisms as clearly as possible, we focus on the steady-state force and the buildup to steady state, rather than treating the feedback loops [24] that cause the polymerized-actin count to drop to zero after reaching a peak. We also leave out possible effects of hydrolysis of actin subunits at filaments tips. The possible consequences of this assumption are analyzed in section 2.5.
Each filament tip interacts with the obstacle via a smooth potential U(r), where r is the distance from the filament tip to the obstacle. U(r) can be either purely repulsive or have an attractive well. The 6 × 6 square central filament subset is assumed to have a deep well in U(r), corresponding to filament-membrane binding. The outer filaments have a purely repulsive potential or one with a shallow well. We adopt a picture similar to that of [33] where the filament 'free length' L protrudes beyond a gel region, which we treat as linear and elastic. We define the filament bases (black squares) as being a distance L in from the tip. At this point the filaments are crosslinked, either by Arp2/3 complex or another crosslinker, and the gel begins (see figure 1). For simplicity, we ignore variations in L from filament to filament. The initial filament base positions are staggered randomly within a subunit length interval, so that each filament tip is at a different position relative to the membrane. The filaments polymerize and depolymerize stochastically at rates that depend on r. The filament tips also diffuse rapidly, because of bending fluctuations. We focus on the component of the diffusive motion perpendicular to the obstacle. Diffusive motion of a filament tip parallel to the obstacle does not change the filament-obstacle interaction, because this motion is much smaller in magnitude than the distance (10 nm or more) over which the filament-obstacle interaction changes character. Similarly, diffusive motion of bases describes elastic deformation of the actin gel induced by the forces from the filaments. In treating this effect, crosslinks between bases of adjacent filaments are modeled as springs that constrain the relative motion of filaments in the direction perpend icular to the obstacle. Both the tip and base motions are described using biased Brownian dynamics driven by the filament-obstacle interaction and a linear restoring force. This force is determined by the filament rigidity for tip fluctuations and the actin gel stiffness for the base motions. The obstacle moves stochastically in response to forces from the filaments, via biased Brownian motion. Figure 2 shows the interaction between a single actin filament and the obstacle. To simplify the calculations, we project the growth of the actin filament onto the z direction (the direction of the obstacle motion), with δ cos(θ) being the projected actin step size. We treat only the z-direction growth explicitly. The default filament direction makes an angle θ = 35 • with the normal to the obstacle, consistent with the 70 • Arp2/3 branching angle [34]. The azimuthal angle is unspecified because it does not affect the calculation results. This approach ignores lateral motion of the filament tips along the membrane resulting from polymerization. The possible impact of this approximation is discussed in section 2.5.

Filament-obstacle interaction
The distribution of pulling and pushing forces on the obstacle is determined by the differences in polymerization and bending/deformation between different filaments. To keep track of these differences, we define for each filament a 'height', which determines the filament-obstacle interaction. The height of filament i is determined at a given time step by where n is the number of subunits added to the filament, z tip is the filament tip fluctuation, and z base is the filament base fluctuation. Staggering of filaments is described by the initial base positions h 0 (i) = α i δ cos (θ), where α i is a random number between 0 and 1. This means that even if all the filaments have grown by the same number of subunits, their heights will differ because their bases are at different locations. These correspond to the positions of Arp2/3 branch points where the filaments anchor in the gel. Changes in h result from either polymerization, filament bending, or motion of the filament base. The values of h do not correspond to actual filament lengths, which are not calculated explicitly. Only differences in h from filament to filament are important for determining the forces. The time dependent gap r between the obstacle at position z obst and a given filament tip is We treat the interaction between the obstacle and the filament tip with smooth idealized potential functions, as shown in figure 3, having the form where A, B, C, D, κ 1 , κ 2 , κ 3 , κ 4 , r 1 , and r 2 are constants. Having a non-zero B in equation (4) adds an attractive well to the potential. We refer to this type of potential as a 'simple well'. A 'double-well' potential is obtained  (4) and (5). Blue curve: 'pusher potential' (equation (4)), with A = 54.66 pN · nm, B = 0, and κ 1 = 0.9 nm −1 . Red curve: 'puller potential' having a well depth of 25 k B T , with parameter values (equation (4) by choosing positive C and D in equation (5). In the double-well potential the broad minimum might represent conformational flexibility of a protein binding the actin filament to the membrane, or the presence of two different binding sites [35,36]. The corresponding forces exerted on the obstacle by the filaments are and The pusher filaments in the outer ring of the array have either only repulsive potential terms or a repulsive potential plus one with a shallow 'simple' well. The puller filaments in the 6 × 6 central region have either a deep 'simple well' or a 'double well' potential.

Stochastic treatment of actin polymerization
The actin on-rate, k on , has been defined in most previous models as an average over a time long compared to the time scale of filament-tip and obstacle fluctuations. Here, such an average is not appropriate because of the large force fluctuations that occur over short times from the combination of Brownian obstacle and filament-tip motion, and the rapid variation of the force between the filaments and the obstacle. Therefore we treat instantaneous rates that apply to a particular position of the filament tip relative to the obstacle. Thermodynamic analysis [25] shows that the instantaneous position-dependent rates must satisfy the relationship on and k 0 of f are free-filament on-and off-rates respectively. In order to concretely determine k on and k off , we assume that their magnitudes never exceed the free filament on and off rate values, so that

Stochastic time evolution of obstacle position, filament bending and gel deformation
The obstacle position and filament bending evolve in time according to random thermal forces and deterministic forces from the filament-obstacle interaction potential, as well as linear restoring forces for the filament tips. For conceptual simplicity, we treat the system dynamics using a 'filament-centric' approach, where the actin gel is assumed to be stationary, while the obstacle moves. To treat force generation in the presence of a stationary obstacle, the calculated motion can be reversed by a simple coordinate transformation, assigning the actin gel diffusion coefficient to the obstacle. The diffusive motions of the obstacle and filament bending coordinates are treated by the discrete form of the Langevin equation ( [37], chapter 3): Here D obst is the obstacle diffusion coefficient, D tip is the filament tip diffusion coefficient, F tot is the total force from the filaments acting on the obstacle, F load is the external force applied on the obstacle (used only in the 'Mean-Force' model described below), ∆t is the time step, and k bend is the tip bending stiffness. Further, α and α are random numbers uniformly distributed between − 1 2 and 1 2 , so that < α 2 >=< α 2 >= 1 12 . Displacements in consecutive time steps are uncorrelated. The values of ∆z tip,i are limited to the range defined by a filament being either perpendicular to the obstacle or parallel to it: −δ cos (θ) < ∆z tip,i < δ[1 − cos (θ)]. The use of a uniform distribution for the individual steps is justified, because the central limit theorem guarantees that after many steps the the displacement distribution will approach the Gaussian distribution that characterizes Brownian motion. We find that already after five time steps or 2.5 × 10 −9 s, the rms difference between the calculated displacement distribution and the Gaussian is only 1.3% of the Gaussian peak height. Gel deformation is treated via stochastic motion of the filament bases anchored in the gel. Interactions between different bases include only nearest-neighbor interactions, and are assumed to be proportional to the difference in their z-coordinates. The base displacements follow an equation similar to that of the tips: Here D base is the filament base diffusion coefficient, k elas is the filament-base spring constant embodying the elasticity of the actin gel, and the random number α has properties identical to those of α . The ∆z nn,i term is the difference between the average z base of a filament's neighbors and the value of z base for the filament itself. The ∆z base,i term is also needed, because in its absence the degree of freedom corresponding to uniform motion of the gel surface and the obstacle by the same amount has no restoring force. The ∆z base,i

Choice of parameters
The parameters are given in table 1. As described above, the obstacle diffusion coefficient D obst corresponds to that of the actin gel, D gel . Estimating D obst from the actin gel size of ∼10 2 nm versus 3 nm-size actin monomers with diffusion coefficient D mon = 5 × 10 6 nm 2 s −1 [38], using the inverse proportionality of diffusion coefficients to size, would suggest D obst = D mon /30 = 1.7 × 10 5 nm 2 s −1 . However, using an obstacle diffusion coefficient this large would be extremely computationally demanding, and for this reason we use a smaller value D obst = 10 4 nm 2 s −1 . We believe that this value is large enough to capture the key physical mechanisms since the dimensionless ratio of t dif f = δ 2 /2D obst (the time that takes the obstacle to diffuse the length of one actin step size) to t pol = 1/C A · k 0 on (the time for one free-filament subunit to add), t dif f /t pol = δ 2 C A k 0 on /2D obst = 0.008, is still very small. Here C A is the actin monomer concentration and k 0 on is the free-filament on-rate constant (see table 1). The diffusion coefficients D tip and D base are expected to be much larger than D obst because the moving entities are much smaller than the actin gel. Again, for computational practicality we use values of D tip and D base , that are 10 times larger than D obst , but probably smaller than the physical values. We have evaluated the sensitivity of the results to these two parameters by halving each one separately, and halving them both at the same time. This changed the total pulling force by only about 1%. We have not tested the effect of increasing the parameter values, because of the large computational effort that would be required. But the smallness of the effect resulting from reducing the values suggests that the effect of increasing them would also be moderate unless the behavior changes very abruptly as a function of the diffusion coefficients. We obtain a = 10 nm for the filament spacing from the estimated number of actin filaments 300 [39] at an endocytic patch with radius R = 75 nm.
We obtain the bending spring constant as k bend = 3 k B TL p /L 3 sin (θ) 2 = 4.17 pN nm −1 [40]. Here L p = 17.5 µm is the persistence length of an actin filament [41] and as discussed above L is the free length of the filament beyond the gel surface. We take L to have a typical value of 54 nm, corresponding to the 20 subunits of a typical actin filament in endocytic actin patch. (The estimated number of actin filaments at an endocytic patch is taken to be 300 [39], with a total of 6000 actin monomers [42]). We assume that L remains constant during the polymerization process, as crosslinkers bind to newly grown parts of the actin gel near the membrane.
The actin-gel spring constant k elast is obtained by fitting it to elastic restoring forces calculated for a configuration of filaments where alternating filaments are displaced in opposite directions (see appendix).

Validity of assumptions
Here we discuss the potential impacts of the main simplifying assumptions that we have made:

Assumption of a sharp boundary between pusher and puller filaments.
Because the endocytic protein patches assemble stochastically, the boundary between the pusher and puller filaments will be blurred. In the regions where there are roughly equal numbers of pullers and pushers, the forces will oppose each other and the 'smeared' force density will be reduced relative to the force densities in the strongly pushing or pulling regions. This will reduce the magnitude of the total pulling force that can be obtained.

Ignoring lateral component of polymerization.
This could have at least two effects. First, filament tips of pushers could move into the puller region. This would blur the boundary between pushers and pullers, reducing the total pulling force as described above. Second, if the filament tips are anchored strongly enough in the membrane that lateral motion is inhibited, forces could build up that would slow the polymerization of the pusher filaments. Again, this would reduce the total pulling force.

Treatment of the actin network as an elastic gel.
Viscous flow of the actin network will inhibit its ability to sustain a distribution of pushing and pulling forces, once again reducing the magnitude of the total pulling force. The viscosity of the actin network at endocytic actin patches is not known, but the magnitude of the effect is estimated in the Discussion section.

Assumption of an infinitely hard obstacle.
During endocytosis, the actin gel interacts with the cell membrane, which in turn interacts with the cell wall. We model only the part of the process before the invagination forms. During this time, the membrane is pressed against the cell wall. The force from the membrane is unlikely to deform the cell wall noticeably. This force density is comparable to the turgor pressure, which is about 200 kPa [10]. On the other hand, the Young's modulus of the cell wall is about 110 MPa and its thickness is about 120 nm, suggesting a maximum deformation of about (120 nm)×(200 kPa)/(110 MPa) = 0.22 nm, very small on the scale of the current simulations. In addition, because the membrane is nearly flat, the bending forces that it generates and transmits to the actin gel are very small. Therefore we feel that using an infinitely hard obstacle is a reasonable approximation.

Absence of hydrolysis effects in model.
If subunits at filament tips hydrolyze and release inorganic phosphate before a new subunit is added, the pusher-filament stall force will be reduced [28]. The rates of hydrolysis and phosphate release are not known for filament-tip subunits in endocytic actin patches. However, as the opposing force increases, the rate of subunit addition will slow, increasing the likelihood of release occurring before subunit addition. Therefore the total pulling force, even for very tightly bound pullers, could be significantly below the value predicted from the thermodynamic stall forces of the pusher filaments.

Force generation and gel deformation by a 12 × 12 square array of cross-linked filaments
The simulations for our system of 6 × 6 = 36 pullers and 12 × 12 − 36 = 108 pushers begin with all the filaments having one subunit, and their bases staggered as discussed above. Pushing and pulling forces develop as the pushers and pullers grow at different rates. The total pulling force is defined as the sum of all of the forces on the pulling filaments, and is taken to have a negative sign. We obtain our results from a single run of 10 s, rather than averaging multiple shorter runs. We use this procedure because obtaining reasonable estimates of the steady-state forces requires us, for each simulation run, to go beyond an equilibration time that can be as long as several seconds. For this reason, multiple runs much shorter than ten seconds would not be valuable because they would be dominated too much by the equilibration time. Multiple runs of ten seconds, for each parameter set, are not viable because of the computer time required. However, a single long run will accurately represent the results of shorter runs if the system is ergodic, so that a time average is equivalent to a configuration average. We are not able to perform enough runs to test the ergodicity for our full 144-filament system. However, we have tested it for a smaller 16-filament system with 4 pulling filaments. We compared results of runs of 10 shorter runs of 5 s with the result of one long run of 30 s. We obtained a force of 62 pN from the shorter runs, a force of 66 pN from a run of 30 s, and forces of 62 and 67 pN from two 10 s runs. Thus for this smaller system using a single long run gives average forces accurate to better than 10%. The force fluctuations in the larger 144-filament system are found to be about three times smaller than in the smaller system, so we expect that the average forces in this case are accurate to within a few percent. Figure 4 shows the time course of the total pulling force for systems with different puller-obstacle interaction potentials (see figure 3), treating a range of values of the actin gel stiffness. The total pushing force (not shown) almost entirely balances the pulling force. The total pulling force is thus less than the sum of the pusher stall forces Weak puller-obstacle binding, as in the 5 k B T potential, generates very small forces. The pulling force approaches the limiting value of 760 pN with increasing well depth, but even for potentials as deep as 50 k B T, the pulling force is substantially below the limiting value. Modifying the gel stiffness (halving k elas for 'soft' gels and doubling k elas for 'stiff' gels) might correspond to changing the concentration of actin filament crosslinkers. This does not change the asymptotic total pulling force, but stiffer gels generate large pulling forces earlier in the process. In addition, soft gels can lead to detachment of the pulling filaments from the obstacle, as occurs in the '25 k B T soft' curve at about 1.8 s, where the force suddenly drops to zero. This effect is described in detail in section 3.6.  The total gel deformation is another measure of the force-generating capability of the system. It is defined as the difference between the average base position of all the pullers and that of the pushers. Figure  The time t max that it takes the force to build up to near its maximum value, or equivalently the time required to reach near maximum deformation, can be roughly estimated using dimensional analysis. We take the maximum force F max to be determined by the total gel deformation ∆u z , the actin gel Young's modulus E, and the radius R ( 60 nm in our model) of the gel. Further assuming that F max is proportional to ∆u z E, dimensional consistency requires that F max = (constant) × ∆u z ER. Taking filaments to grow at their zero-force values up until the stall point, we have ∆u z = k 0 on C A δ cos (35 • )t max . Then The inverse proportionality of t max to E is seen in figures 4 and 5, where the force and deformation reach their maximum values for stiff gels faster than for softer gels. The numerical values predicted by equation (15) are not expected to be accurate, but t max is estimated to be about 0.8 s for the case of the medium gel and a 50 k B T puller. By comparison, the time required to reach half-maximum deformation in the simulations is about 1 s.

Effects of transient attachment of pushers to obstacle
It is believed that the WASP family of proteins, and their yeast homologue Las17, create weak transient attachments between filament tips and membrane [45]. We thus calculate how adding a potential well to the pusher potentials affects the magnitude of the pulling force. Figure 6 shows the magnitude of the total pulling force as a function of the pushers' well depth. For each point the mean value is calculated from the  last seven seconds of a ten-second run, to minimize contributions from the equilibration period. The error bars are obtained as the standard deviation of the mean of force values from seven consecutive one-second pieces of these runs. Thus they include both random error and some component of systematic error. For both 25 k B T and 50 k B T pullers, the force drops as the pushers' binding to the obstacle becomes stronger, and the fractional effect is larger for weaker puller potentials. At a well depth of 10 k B T, the drop in total pulling force corresponds to about 1 pN per pusher. Thus the total pulling force is significantly stronger when pushers have purely repulsive interactions with the membrane.

Force distributions
To explore possible spatial variations of the pulling and pushing forces, we present their spatial distribution using heat-map diagrams. Because the simulations are stochastic, the distributions obtained over a finite time interval (figure 7(a)) display noticeable fluctuations.
To show the systematics of the distributions more clearly, we create symmetrized force distribution heatmaps by averaging the filament forces over symmetrically related subsets ( figure 7(b)). For example, the symmetrically averaged force for a filament with coordinates (x, y) includes contributions from filaments at (±x, ±y) and (±y, ±x). We consider the case of 50 k B T pullers, which gives close to the limiting pulling force. We see a fairly flat distribution in the pushing region (blue) balancing the total pulling force in the center (red). Figure 8 shows a cross section of the time-averaged force distribution along a row in the middle of the array, for four different potentials. As expected, the magnitudes of the individual pulling forces are always less than the average maximum force per pulling filament The pulling and pushing force distributions are relatively flat in all cases except the double well potential. We believe this holds because in arrays with single-well pullers, all pulling filaments must grow at the same velocity once steady state is reached; since  force determines growth velocity, all pullers must also experience the same force, and the force distribution is flat. In contrast, double well potentials have a velocity that is nearly force-independent after the pulling force reaches a certain magnitude [25]; therefore a range of forces is compatible with a certain velocity, so equality of growth between pulling filaments does not necessarily lead to a flat distribution of force. This is seen in the 'DW 25 k B T' data points in figure 8. The force distribution extracted from EM images of invagination shapes ( figure 8(b) of [22]) shows a flat distribution over the pulling region similar to our prediction. However, the extracted force profile shows a ≈30% hump in actin force density at the edge of the pushing region (r ≈ 50-80 nm). We do not see this hump in our simulation results. We believe that the hump results from the hemispherical geometry assumed in [22]. This choice makes it easier to shear the gel near the edge of the hemisphere, because it is thinner there. This should in turn reduce the magnitude of the forces. In the present model, we have taken all filaments to be elastically equivalent, which makes the resistance to shearing the same everywhere. Given that the hemispherical shape itself is an idealization, we do not feel that weakening the shear strength near the boundary to mimic the hemispherical shape would render the model more realistic. Figure 9 shows the time evolution of the force distribution for the case of 50 k B T pullers. The 'Early' stage shows the forces averaged over the 1s < t < 2s interval of a 10 s long simulation. As seen in figure 4(b) this stage is well before steady-state, explaining why the distribution has not reached its asymptotic constant value. The filaments at the edge of the pulling region show enhanced forces at this stage. We believe this occurs because at the early stages of the simulation, the forces are not yet strong enough to slow polymerization greatly. Then all the pullers will have added roughly the same number of subunits, and all the pusher filaments will have added a constant number of subunits (larger than the value for the pullers). The corresponding gel deformation will have a constant value in the puller region and a different constant value in the pusher region. Maintaining the difference between these deformations requires a force dipole at the boundary between the pusher and puller regions, which we believe explains the peak in the force seen in the 'Early' results. The 'Late' stage is closer to the asymptotic one, as the network approaches steadystate. The 'Middle' stages has force peaks similar to those seen in the 'Early' results, but also has a bump in the middle, for which we do not have a simple physical explanation.

Mean-force theory
To clarify the mechanisms determining the magnitude of the pulling force, we study a simplified meanforce model based on the force-velocity relations of separate puller and pusher ensembles. In the full system simulations, 108 pushers and 36 pullers exert force on the obstacle simultaneously. The forces exerted by the pullers experience fluctuations in time due to the polymerization dynamics of the pushers, and vice versa. However, on average, the magnitude of the force felt by the pushers equals that felt by the pullers. One would then expect the steady-state force to have the value at which the puller and pusher growth velocities are equal. To explore this hypothesis, we developed a mean-force model based on two force-velocity relations v push (F) and v pull (F). Here F (positive) is the magnitude of the time-averaged total force felt by either the pullers or the pushers. We calculate v push (F) and v pull (F) by performing pusheronly or puller-only simulations using an external force    of magnitude F, with the force pointing in opposite directions in the two cases. The condition determining F is that v push (F) = v pull (F). The value of F satisfying this condition is obtained by linear interpolation from a finite set of force calculations. Figure 10 shows the curves of v push (F) and v pull (F), as well as the comparison between the predicted force from meanforce model with that of the full-system simulations.
As is clear from figure 10(a), slow puller growth will bring the crossing point between the F-V curves in figure 10(a) down and to the right, increasing the pulling force. On the other hand, slowing the growth of pushers will bring the crossing point down and left, reducing the pulling force. As figure 10(b) shows, the mean-force model closely predicts the results of simulations with a full array of filaments, over a range of puller potentials. This result suggests that the time fluctuations of the forces from pushing filaments may not crucially impact the growth velocity of the pullers, and vice versa. This analysis shows that a large pulling force will occur when the puller filaments slow pusher filaments' growth by a large factor. The maximum force occurs when the puller filaments' growth stops completely, stalling the pusher filaments. Figure 11 shows the relationship between total pulling force and average growth velocity of the central filaments. As the plot shows, even a growth rate as low as 10% of the free filament velocity can reduce the generated force substanti ally, by >40%. As described in the Discussion, this is important for ascertaining the effect of viscous flow on force generation.

Effect of mechanical parameters on extent of actin polymerization
We find that softening the gel increases the extent of actin polymerization required to reach a certain force. Figure 12 shows the number of subunits added to all of the filaments in the array during the course of the simulation. We compare different gel stiffnesses and potential depths at a given value (270 pN) of the force, about 80% of the maximum force for 25 k B T pullers.
The total amount of actin polymerization is increased by about 160% for both 25 k B T and 50 k B T pullers going from a stiff gel to a soft one. In addition, ≈40% more actin is polymerized for a medium gel at the same force of 270 pN when the pullers' potential depth is halved from 50 k B T to 25 k B T. This increase is smaller for soft and stiff gels (≈25%). The black 'No obstacle' bar corresponds to the amount of polymerization obtained after 0.5 s (the average time it takes for the force with 25 k B T pullers to reach 80% of maximum) if no pullers are attached to the obstacle.

Obstacle-gel detachment
Our finding of obstacle-gel detachment in some parameter ranges is surprising, given the magnitude of the potential and the number of filaments. Consider the case of a completely rigid, solidly anchored gel. Standard reaction rate theory [46] gives an analytical estimate of the escape time of the obstacle from the potential well: Here w and ∆U are defined in figure 13, where a positive external force of 335 pN (which caused detachment in figure 4 above) pulls on a 6 × 6 array of rigid pullers with 25 k B T potential wells and staggered initial alignment. In this case, w ≈ 5 nm, D obst = 10 4 nm 2 s −1 (as above), and ∆U ≈ 180 k B T. Thus, equation (16) gives an extremely long time of ≈ 4 × 10 75 s for the obstacle to detach from these pullers, so detachment essentially never happens. Consistent with this prediction, when the actin gel is stiff and the binding energy is 25 k B T, our simulations find that obstacle-gel detachment never occurs during 10 s simulations. However, softening the gel or weakening the puller-filament binding can lead to detachment before steady state is reached, as was seen in figures 4 and 5 for the 25 k B T soft-gel case. Thus a soft actin network behaves differently from stiff networks in that it detaches more easily from the obstacle. We evaluated detachment in 40 different runs for the 25 k B T softgel case. We found the distribution to be peaked with an average value of 1.7 s and a standard deviation of 0.3 s. The distribution does not have the exponential form that would result from a stationary Poisson process. We believe that the peaked behavior in this case means that the detachment occurs rather rapidly once the force has reached a sufficient value to drive detachment, which occurs near 1.7 s. To understand the detachment process in more detail, we performed puller-only stimulations with an external force acting on the obstacle. In a particular simulation with 25 k B T pullers in a gel with medium elasticity, with an external pulling force of 335 pN on the obstacle, the obstacle detached from the gel after about 6.5 s. To understand the origin of this effect, we look at the distribution of r, the gap between filament tips and obstacle, during a few time steps right before the rupture happens. Figure 14 shows frames of a heatmap plot of the distribution of r for the pullers. Larger r values (greater stretching) have redder pixels. Frame (a) is at a time well before rupture, to show the baseline appearance of the distribution. Frames (b) through (m) span 12 000 time steps, corresponding to 6 µs. The number of pink pixels increases gradually during this period, indicating the appearance of possible detachment nucleation points around the tip of the actin filaments. The accumulation of these nucleation points eventually spreads from the right bottom corner over the entire pulling region, and the obstacle detaches completely. This is reminiscent of the process of a crack propagating between two dissimilar mat erials due to stress concentration at the edge of the crack.
Some light can be shed on the effect of gel stiffness on the detachment process via the 'Griffith' theory of the critical stress for fracture [47]. The Griffith theory describes how stress concentration at a crack tip aids crack propagation. The stress σ c required to propagate a crack of length l c is given by where γ is the energy density required to break the bonds along the crack. In this criterion, σ c increases with E, consistent with our finding that stiff gels do not detach from the membrane. In applying it to the actin-filament system, we take γ to be the ratio of the puller potential depth to the area a 2 per filament. Since we are studying an incipient crack, we take l c = a = 10 nm. For a puller potential well depth of 25 k B T, using the value E = 0.140 pN nm −2 , we obtain σ c = 0.095 pN nm −2 . By comparison, the stress at a force of 335 pN, for which detachment occurs in the simulations, is σ = 0.17 pN nm −2 . Given the small size of the simulation system and the absence of a welldefined preexisting crack, this level of agreement is reasonable.

Discussion
Our calculations show that actin-based pulling force generation can result from spatially separated ensembles of filaments having different growth velocities. The slower-growing filaments exert pulling forces on the membrane, while faster-growing filaments exert pushing forces. Large pulling forces are generated if the velocity of pullers at a given force is much smaller than that of pushers. When the pullers do not grow at all, the total pulling force is maximized and equals the sum of the stall forces of the pushing filaments. The mechanism explored here for generating slower puller-filament growth is a stronger binding to the membrane, but others may be relevant.
A mean-force model treating pushers and pullers separately with constant forces reveals in more detail how the pulling force is determined by the difference between the puller and pusher force-velocity relations. It accurately predicts the total pulling force from fullsystem simulations, suggesting that time statistics of the force generation are not important for the final results. Our key specific findings are the following:

The pulling force is reduced by growth in the pulling region, viscous deformation of the actin gel, and transient attachment of the pushers to the membrane
When the binding of puller filaments becomes weaker and they grow faster, the pulling force drops rapidly ( figure 4). A similar effect would be expected from viscous flow of the actin gel. Figure 11 shows that a growth rate of 13.6 nm s −1 is sufficient to reduce the pulling force by about 50%. Taking our model system to have a radius of about 60 nm gives a shear rate of about 0.2 s −1 . Since the pulling force per filament is on the order of 20 pN (see figure 8(a)) and the filament spacing is about 10 nm, the stress in the puller region is ∼20 pN/100 nm 2 = 0.2 pN nm −2 . Then substantial force reduction would occur if the actin gel viscosity becomes less than about 0.2 pN nm −2 /0.2 s −1 = 1 pN s nm −2 = 10 6 Pa s. This value is much higher Figure 14. Distribution of r-values for pullers with 25 k B T binding in a medium-stiffness gel. Redder pixels correspond to larger r values. Frame a is the baseline distribution, at a time well before the rupture. Frame (m) is the distribution after the obstacle detaches completely. Frames (b) through (m) span 12 000 time steps, a total of 6 µs. Accumulation of nucleation points around several filament tips (frame (h)), which takes place in a time in the order of microseconds, leads to detachment. than any that have been previously measured. However, the viscous properties of actin networks at the very high actin and crosslinker densities present in endocytic actin patches have not been explored quantitatively. Pulling forces are also maximized when the pushers' interaction potential with the obstacle is monotonically repulsive (figure 6), which gives the largest growth velocity [25]. Transient attachments of the pushers to the obstacle slow their growth and thus inhibit pulling force generation.

Effect of actin gel mechanical properties on pulling-force generation
In the steady-state limit, the magnitude of the pulling force is independent of the actin gel stiffness over the range studied, provided that obstacle-gel detachment does not occur ( figure 4). However, large forces are obtained at earlier times for stiff gels. This suggests that invagination should be slowed by mutations reducing the number of crosslinkers. The endocytic event could also be aborted completely because the actin patch has a finite lifetime, which may be too short to allow the maximum force to build up. This is consistent with requirement of the yeast fimbrin homolog Sac6 for endocytosis in yeast [15,16]. In [15], more than 70% of the endocytic sites in Sac6∆ cells were found to have a flat membrane profile. The majority of the remaining ones invaginated to a distance of ≈100 nm, but then retracted. The 'retraction' phenotype could result from insufficient force building up during the lifetime of the actin patch. The slower buildup of force for the soft gels could also prevent endocytosis if viscous flow of the gel is important. Such flow during the period of force buildup might prevent the force from ever becoming large enough to overcome the turgor pressure.

Force distributions
The steady-state force distribution profiles (figures 7 and 8) reveal a complete force balance between pushing and pulling regions, and fairly constant force densities over these regions. They are roughly consistent with the force profiles obtained from measured membrane deformations [22]. The time evolution of the force distributions shows enhanced forces at the edge of the pulling region at early stages, which are part of a dipole of forces surrounding the interface between the pusher and puller filaments. The force distributions could be measured using a combination of superresolution microscopy and suitably designed molecular force sensors. If a force sensor were inserted into Sla2, and the signal from the force sensor measured using superresolution microscopy, a picture of the force distribution in the central Sla2 region could be obtained and compared to these predictions.

Effect of mechanics and filament-membrane interaction on actin polymerization
We also find that the actin count is increased by either softening the gel or reducing the puller binding energy (see figure 12). We are not aware of data showing the effect of gel softening on the actin count, but an extreme version of reducing the puller binding energy is obtained in Sla2 deletion mutants. In these experiments, extensive actin accumulation has been observed in the form of 'comet tails' [48,49].

Detachment of pullers from the obstacle
This effect can completely disrupt the force-generation machinery. We find that detachment does not occur for completely rigid gels. But for soft gels, it does occur as seen in figure 14 and described by equation (17). The process begins with an initial 'nucleation' event in which one or more filament tips move out of the potential well binding them to the obstacle. This leads to stress concentration, and detachment then spreads over a microsecond time scale, in a mechanism analogous to crack propagation in a solid. These results suggest that softening the actin gel driving endocytosis in yeast by mutations reducing the number of crosslinkers could cause actin gel detachment from the membrane, aborting the endocytic event. This mechanism could provide an alternate explanation of the requirement for Sac6 in endocytosis [15,16]. The detachment mechanism could be distinguished from the direct effect of gel softening by tracking the motion of patches of the actin proxy Abp1, in Sac6∆ cells. If the detachment mechanism operates, the patches will move into the cell rapidly; if the direct effect of gel softening are more important, the patches will remain at the membrane. where we have taken Poisson's ratio σ = 1/2 to correspond to an incompressible actin gel. To obtain the restoring forces for the given displacements, we first express this result in Fourier space: The Young's modulus for the actin gel has been roughly estimated as E = 0.14 pN nm −2 [22]. Thus equation (A.7) obtains the baseline value for the actin gel spring constant as k elas = 0.53 pN nm −1 .