New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Selforganization of actin filament orientation in the dendriticnucleation/arraytreadmilling model

Contributed by Edwin W. Taylor, March 6, 2007 (received for review November 16, 2006)
Abstract
The dendriticnucleation/arraytreadmilling model provides a conceptual framework for the generation of the actin network driving motile cells. We have incorporated it into a 2D, stochastic computer model to study lamellipodia via the selforganization of filament orientation patterns. Essential dendriticnucleation submodels were incorporated, including discretized actin monomer diffusion, MonteCarlo filament kinetics, and flexible filament and plasma membrane mechanics. Model parameters were estimated from the literature and simulation, providing values for the extent of the leading edgebranching/cappingprotective zone (5.4 nm) and the autocatalytic branch rate (0.43/sec). For a given set of parameters, the system evolved to a steadystate filament count and velocity, at which total branching and capping rates were equal only for specific orientations; net capping eliminated others. The standard parameter set evoked a sharp preference for the ±35 degree filaments seen in lamellipodial electron micrographs, requiring ≈12 generations of successive branching to adapt to a 15 degree change in protrusion direction. This pattern was robust with respect to membrane surface and bending energies and to actin concentrations but required protection from capping at the leading edge and branching angles >60 degrees. A +70/0/−70 degree pattern was formed with flexible filaments ≈100 nm or longer and with velocities <≈20% of free polymerization rates.
The polymerization of soluble actin monomers between filament “barbed ends” and the plasma membrane (PM) generates the force of protrusion in cell motility (1, 2). Other proteins required for lamellipodial motility (3) are arp2/3, which nucleates (branches) free barbed ends at ≈70 degrees from existing ones (4); a PMbound activator of arp2/3 (5); ADF/cofilin, which promotes the depolymerization of pointed ends (6) and perhaps debranching reactions (7); and capping protein, a terminator of barbed end growth (8). The generation and persistence of lamellipodia from these elements is described in the “dendriticnucleation/arraytreadmilling” conceptual model (2, 4, 9). This model can be subdivided into three main processes: the kinetics of filament (de)polymerization, branching, and capping; filament–PM interactions, which limit polymerization rates; and the diffusion of actin monomers and other soluble components. Such a system is “complex” in the sense that many copies of each component type interact to exhibit “emergent” system properties not expected from the individual rules of interaction (10). In contrast to “complicated” systems of many dissimilar components with precisely defined interactions, complex systems can selforganize and adapt to environmental change.
An important emergent property is the selforganization of lamellipodial actin filaments into orientations at ±35 degrees with respect to the direction of protrusion (11, 12). Maly and Borisy (11) predicted ±35 or +70/0/−70 degree patterns with a 2D mathematical model based on these dendriticnucleation assumptions. A model by Atilgan et al. (13) allowed 3D branching, but required preferential arp2/3 orientation in the PM for pattern formation. The numerical model described here extends the Maly and Borisy (11) model, removing most of its simplifications and limitations. Reactions are treated stochastically, the elastic properties of the membrane and filaments are included, and the timedependence of the distribution's evolution is obtained. We preserve the assumption that filaments remain oriented in the 2D lamellipodial plane. Our simulation exhibits orientational selforganization and permits the determination of a range of parameter values consistent with pattern stability. We reveal the transient development of the orientation pattern and the approach to steadystate protrusion velocity and filament number.
Several computer models of (rigid) bacterial propulsion by rigid filaments have been proposed (14–16), with those of Carlsson demonstrating a flat force–velocity relationship under autocatalytic branching (15, 16). Mogilner and Oster (17, 18) developed the basic theory of the actinbased elastic Brownian ratchet, which applies smallangle elastic beam theory to thermal filament fluctuations. The present model combines stochastic propulsion and elastic filament models, adding a flexible PM load and thermodynamically realistic filament–membrane interactions.
Assumptions and Methods
An ≈1 μmwide (X) portion of a flexible lamellipodial leading edge (LE) was simulated, with cyclic boundary conditions for all components on the ≈1 μm Y axis edges and a fixed actin monomer concentration [A] _{TE} on the trailing edge. Every filament over the entire lamellipodial thickness was modeled, with 2D positions, orientations, and end states of each filament recorded individually. Over each small time step Δt, the calculation algorithm performed spatially discretized (Fick's law) diffusion and MonteCarlo (stochastic) kinetics calculations for all components, with iterative calculations of PM and filament mechanics as required (Fig. 1 and Tables 1 and 2). Consistent with experimental indications that branching and anticapping mechanisms operate very near the LE, free barbed ends within a Ydistance ε branched new filaments at rate R_{br} and were blocked from the usual capping rate R_{cp} . New filaments deviated from the parent filament barbed end orientation by a normal distribution about the mean branch angle ±(θ _{br} ±σ _{br} ). Polymerization occurred at a rate proportional to the local actin monomer concentration, reduced at the LE by a Boltzmann factor based on the total system energy required for that specific protrusive step. Filaments were assumed to be either rigid or flexible and cantilevered from the nearest branch or pointed end. Values for the filament and PM mechanical properties and most of the reaction rate constants were available. The main unknowns were the effective rates of branching and capping, critical to the development of the filament distribution. With R_{cp} available, ε and R_{br} were determined from simulation (see Fig. 3), completing the standard set of parameters listed in Table 2. Model details and a simulation video depicting individual protrusive steps can be found in supporting information (SI) Methods and Movie 1.
Results and Discussion
SelfOrganization of Filament Orientation Arose Under Any Initial Conditions (ICs).
We first asked under what ICs the model system would organize filament orientation. A set of filaments at random position, orientation, and length has no order, with only the position of the branchinducing and cappingprotective LE to generate directionality. A simulation was run under these ICs (Fig. 2 a), with the branch angle SD (σ _{br} ) equal to zero and otherwise standard parameter conditions (Table 2). Results show the disappearance over 120 sec of total filament length (mass) in most orientations, with one pair of very sharp orientations at ±35 degrees dominating (Fig. 2 b and SI Movie 2). Because the branch angle did not vary, filament populations always occurred in strictly complementary families (related by exactly 70 degrees) and were not able to drift into other orientation families over successive generations of branching. This result shows that from ICs without order, very narrow populations of filaments oriented at 35 degrees show net growth, whereas others show net depolymerization.
When the number of filaments, instead of the total filament length, in each direction is plotted, we see a similar pattern with the addition of backwardfacing filaments at ±105 degrees (Fig. 2 c). Other orientations showed a net decrease in the number of filaments, showing that the total capping rate outpaced the total generation rate for those orientations until their numbers reached zero. Domination by length at ±35 degrees thus did not represent a system in which ±35 degree filaments became particularly long, but instead one in which other filament families were not viable. Filaments that faced backward at ±105 degrees were generated at nearly the same rate as forwardfacing filaments, but rapid capping away from the LE kept their total lengths (and nonproductive monomer depletion) very low (Fig. 2 b), and subsequent debranching and depolymerization reduced their count moderately (Fig. 2 c).
It remained to be shown that the model could adapt from initial orientation patterns devoid of 35 degree filaments into the steadystate (SS) ±35 degree patterns. This is necessary because a feature of motile cells is their ability to reorganize filament orientations in response to changes in the direction of protrusion, without which we would not observe the ±35 degree pattern in cells that turn. Initial conditions of filaments oriented near −70/0 degrees simulated a sudden, 35 degree clockwise turn in protrusive direction from a previous SS pattern (Fig. 2 d). The standard σ _{br} of 7 degrees was used allowing the population to drift into new orientations, but no 35 degree filaments existed initially. Orientations became symmetrical within 1 min, with a broader distribution of orientations centered on ±35 degrees (Fig. 2 e). The SS distribution compared favorably to those measured via digital imageprocessing techniques (Radon transform) from lamellipodial electron micrographs (EMs) (11) (Fig. 2 e, with an image at SS shown in Fig. 2 f). We conclude that ICs of all protrusive simulations with standard parameter values evolve and adapt into ±35 degree patterns consistent with EMs.
Simulation Indicated Plausible Branching Parameters.
The standard parameter set used in Fig. 2 and the rest of this study came from both the literature and simulation. Although experimental data provided relatively direct estimates for N_{fb} , R_{cp} , and [A] values, R_{br} and ε were more difficult to estimate. Modeling allowed us to determine R_{br} , ε, and the filament bending length parameter f_{ts} , consistent both with estimates for N_{fb} , R_{cp} , and [A] and with experimentally accessible values such as orientation pattern and filament length.
A particular N_{fb} is set only indirectly, through filament turnover parameters R_{br} , R_{cp} , and ε. For a given V_{PM} /V_{free} , the profile of N_{fb} with distance from the LE is similar (Fig. 3 a, the exact shape of which will be considered in a subsequent publication). At SS, a fraction of these free barbed ends, f _{ε}, are within the εdemarcated region, and the total rates of filament branching and capping are balanced: (N_{fb} f _{ε}) R_{br} = (N_{fb} [1f _{ε}]) R_{cp} . The SS f _{ε} therefore equals R_{cp} /(R_{br} +R_{cp} ), the setpoint for a negativefeedback loop controlling N_{fb} and V_{PM} /V_{free} : A low V_{PM} relative to V_{free} allows filaments to keep up with the LE more often, raising f _{ε} above the setpoint. The resulting net branching increases V_{PM} /V_{free} . Conversely, a high V_{PM} /V_{free} results in net capping and a decrease in N_{fb} and V_{PM} /V_{free} .
Plots were made of the R_{br} values required to sustain a steady N_{fb} at a given ε, with two values each of N_{fb} , R_{cp} , f_{ts} , and [A] specified and the diffusion coefficient (D) set very high to maintain [A] uniformly (Fig. 3 b). The value of [A] was not a factor in the required R_{br} ; it affected both V_{PM} and V_{free} (= k_{on,brb} [A] δ) equally, and curves are superimposed. (Comparison is therefore by V_{PM} /V_{free} throughout this study.) In contrast, maintaining higher N_{fb} values required a higher R_{br} because the associated increase in V_{PM} /V_{free} diminished f _{ε}. A decrease in R_{cp} was associated with an approximately proportional decrease in required R_{br} , as expected. Filament flexibility was allowed in three parameter sets, and moderately raised the required R_{br} . The most sensitive parameter was ε itself. Decreasing ε dramatically decreased f _{ε}, requiring a compensatory increase in R_{br} .
To narrow the range of acceptable values of ε and R_{br} , the quality of the SS orientation distribution was assessed. Quality was quantified as the mean absolute deviation (α _{d} ) from the mean orientation angle. The α _{d} values for all parameter sets were very similar at the same ε (and R_{br} ) values, with some increase in α _{d} with filament flexibility (Fig. 3 c). Polar histograms of orientation distributions over a range of α _{d} showed focused distributions at high ε and noisy and distorted distributions, incompatible with EM data, at very low ε (Fig. 2 d). These results are compatible with ε values as low as 1 monomer length (δ = 2.7 nm), but not less.
To limit the upper values of ε, we analyzed average filament lengths and orientation pattern development times. Again, because large ε values were associated with low R_{br} , filaments turned over slowly and thus grew very long before being capped. Using EMs, conservative visual estimates of the average length of filaments within 1 μm of the LE ranged from 50 to 500 nm (data not shown). In simulations, the number of filaments was a decaying function of length, with all parameter sets consistent with a maximum 500 nm average length when ε was held to <2.5 δ (Fig. 3 e). Supporting this limit were measurements of the halftime of orientation pattern development (t_{1/2} ) from a 15 degree simulated turn (i.e., ICs of a distribution centered on +20 and −50 degrees), with full recovery defined as a SS α _{d} value. At ε = 3.0 δ, all parameter sets yielded a t_{1/2} of 60 sec or more (Fig. 3 f), corresponding to ≈3 min recovery times for a minor 15 degree turn. In comparison, fibroblasts have lamellipodial protrusion persistence times on the order of 1–2 min before retraction (28). Length and recovery time considerations were thus limited by two criteria to ε < 2.5 δ.
Of related importance to t_{1/2} is the number of branching generations required for network selforganization. The product of t_{1/2} and R_{br} yields the number of generations required for a 50% recovery in α _{d} after a rapid, 15 degree turn (Fig. 2 g). Regardless of parameter conditions, this N_{1/2} value was consistently ≈5 generations for ε = 2.0 δ. Three halftimes (87.5% recovery) required 15 generations. For ε = 3.0 δ, three halftimes still required ≈9 generations. These generations must be largely successive (i.e., parents branching children, then children branching grandchildren) for adaptation to occur.
Based on these studies, we chose a standard set of parameter values consisting of ε = 2.0 δ, R_{cp} = 6.0 /sec, and [A] = 12 μM, with which reasonable protrusion rates of 8 μm/min were achieved. Furthermore, an N_{fb} of 200 fil per μm was compatible with measurements of filament barbed end count (20) and length, and required an R_{br} of 0.43/sec. Filaments were held rigid (f_{ts} = 0) for simplicity except where noted.
The SS Orientation Pattern Is Sensitive to Several Model Parameters.
In the standard model, the same ε value limited capping protection and branching initiation. When branching was allowed at any free barbed end, regardless of position, but capping protection was maintained within ε, orientation distributions remained unchanged (Fig. 4 a). These conditions were effectively similar to standard conditions because the high intrinsic capping rate quickly removed free barbed ends arising beyond ε. (Although we did not simulate branching from the sides of filaments away from their barbed ends, we predict a similar result.) In contrast, when branching was only allowed within ε and capping was allowed anywhere, the orientation pattern changed dramatically. Cases in which N_{fb} and R_{br} were maintained at standard values required a greatly diminished R_{cp} of 0.22 per sec, whereas those with standard N_{fb} and R_{cp} required a greatly elevated R_{br} of 12 per sec. Both cases resulted in distributions with most filament mass facing backward, incompatible with EM results. Protection from capping near the LE is thus essential for the formation of the ±35 degree distribution.
Membrane model parameters were not found to significantly affect orientation distributions over a wide range of values. Given a constant R_{br} , runs varying the resistance to protrusion (γ _{se} ) retained the same V_{PM} and ±35 degree distribution due to the compensating effect of N_{fb} (Figs. 3 a and 4 b). This flat force–V_{PM} relationship is a characteristic of unrestrained autocatalytic branching, observed by Carlsson (16). Membrane flexibility, specified by the bending energy coefficient, κ _{b} , had a significant effect on the LE shape over length scales of <200 nm. Although the pattern is formed with respect to the LE perpendicular, this had no effect on the orientation distribution, with the same pattern attained over two orders of magnitude of κ _{b} variation (Fig. 4 c). We attribute this stability to the LE's horizontal average orientation and rapid fluctuations relative to filament lifetime.
Conversely, filament bending stiffness had a large influence on the orientation pattern. This stiffness is dependent on the third power of effective filament bending length [(f_{ts} l_{ts} )^{3}, with the fraction f_{ts} simulating crosslinking]. Standard parameter simulations with effective filament lengths averaging 0 (f_{ts} = 0, equivalent to rigid) or 41 nm (f_{ts} = 0.25) yielded ±35 degree distributions, but 99 or 139 nm lengths (f_{ts} = 0.50 or 0.75, respectively) resulted in distributions similar to +70/0/−70 degree triplets (Fig. 4 d). These length limits are consistent with Mogilner and Oster's (17) force generation calculations. Triplets were a result of barbedend “splaying” to higher angles under load.
The branch angle (θ _{br} ) itself was not important to the final ±θ _{br} /2 distribution over values of 70–110 degrees, but smaller θ _{br} values resulted in a single uniform mass (Fig. 4 e). Allowing splaying by filament flexibility did not revive the twopeak distribution (data not shown). Note that filaments branching from θ _{br} /2 no longer face backward at θ _{br} = 60 degrees, but rather 30 + 60 = 90 degrees. This allows them an increased opportunity to further branch and disperse the population orientation. Among θ _{br} values that form twopeak distributions, the t_{1/2} from ICs (uniformly distributed over the range of ±90 degrees) is also sharply minimized at θ _{br} = 70 to 80 degrees (data not shown). A 70 degree branch angle is therefore optimized for rapid generation of a stable twopeak orientation pattern.
Changes in V_{free} had no effect on the orientation pattern (data not shown), but the ratio V_{PM} /V_{free} was important. Holding V_{free} constant, the SS V_{PM} was diminished by decreasing R_{br} (this raises the SS f _{ε} requirement; see Fig. 3 a). An R_{br} of 0.06/sec yielded V_{PM} /V_{free} = 0.20 and a ±35 degree distribution (Fig. 4 f). Decreasing R_{br} to 0.03/sec lowered the SS V_{PM} /V_{free} to 0.16 and produced an SS triplet distribution, but also lowered N_{fb} to 80 fil per μm. Raising γ _{se} to restore N_{fb} had no effect on the distribution or V_{PM} /V_{free} (plotted). Note that low R_{br} can lead to long filament lengths and slow selforganization (Fig. 3), and triplet patterns in vivo would be more suggestive of a small ε value (which instead lowers V_{PM} to reach the same required f _{ε}). Note also that, on average, a 70 degree filament advances at V_{free} cos (70) = 0.34 V_{free} , but that stochastically growing filaments bounded ahead by the PM do not maintain a high f _{ε} at this velocity and in fact require V_{PM} /V_{free} < 0.20 to dominate under standard values.
The Evolution of the Final Orientation Pattern and the Development of the SS V_{PM} Were Interdependent.
From ICs of a small number of filaments at random orientation, Fig. 5 a tracks the total number of filaments in each of three complementaryorientation “families” over time. At the low N_{fb} and V_{PM} values early in the simulation, the +70/0/−70 degree triplet multiplied at the highest rate. With increasing V_{PM} , however, the number of filaments decreased in turn for each family until only the ±35 degree pair remained at SS V_{PM} .
Following Maly and Borisy's evolutionary “fitness” concept, fitness = (rate of change of N_{fb} )/(N_{fb} ) for each orientation. The value is positive if a population is growing in number, zero if constant, and negative if decreasing toward extinction. Fig. 5 b shows the fitness landscape for standard branching parameters as a function of φ _{fil} and V_{PM} /V_{free} , with the two symmetric patterns, +70/0/−70 and ±35 degrees, having maximum fitness at low and high velocities, respectively. When either pattern is rotated slightly, the capping rate of the filament with the largest φ _{fil} increases faster than the branching rate of its parent, causing a net decrease in fitness. Because ≈2/3 of new filaments branching from the triplet contribute to the pattern (1/2 each of branches from ±70 degree filaments, and 2/2 branches from 0 degree filaments ≈4 of 6 new filaments), whereas only 1/2 of branches from ±35 degree filaments contribute to their pattern, the triplet has a superior reproductive rate and fitness at low velocities. At higher protrusion rates, this effect is overcome by the low f _{ε} of slow 70 degree filaments.
To identify the SS condition from Fig. 5 b, we note that, under any model context, an increase in SS N_{fb} always increases SS V_{PM} /V_{free} due to the higher net rate of protrusive kinetic events (Fig. 5 c). A typical path through the fitness landscape can therefore be shown as a Darwinian evolution of patterns through the flattened landscape representation in Fig. 5 d. At low N_{fb} , all forwardfacing orientations interact with the same lowvelocity PM and exhibit net branching [f _{ε} R_{br} > (1f _{ε}) R_{cp} for every φ _{fil} ]. This leads to an increase in N_{fb} and V_{PM} and to an environment in which filaments at orientations near 0 and ±70 degrees decrease in count (Fig. 5 d, i). Because any V_{PM} with orientations of positive fitness will ultimately result in an increased total N_{fb} , this cycle will continue (Fig. 5 d, ii) until reaching the equilibrium V_{PM} and its only viable orientations, near ±35 degrees. There, f _{ε} R_{br} = [1f _{ε}] R_{cp} for 35 degree filaments (i.e., fitness = 0), but all other orientations are driven to extinction by net capping. Given different branching parameter values consistent with a low V_{PM} , the triplet reproduces itself at a higher rate than the twopeak distribution, consequently mandating a slightly higher equilibrium V_{PM} and the extinction of the ±35 degree pattern. In either case, populations of orientations reproduce and drive environmental changes (V_{PM} /V_{free} ), that lead to changes in fitness and, ultimately, to a single symmetric pattern and velocity.
Conclusions
We have developed a comprehensive 2D model of lamellipodial protrusion based on the dendriticnucleation/arraytreadmilling mechanism, incorporating diffusion, stochastic kinetics, and elastic filament and PM models. A standard set of model parameter values were determined both directly from the literature and indirectly by using the criteria that model results of filament length, orientation patterns, and development times must be consistent with observed properties. Following these criteria, free barbed ends were protected from capping and allowed to branch at R_{br} = 0.43/sec within a distance ε of 5.4 nm from the leading edge. This conferred the only directionality to the system.
The model accounts for the essential dynamic properties of the network, including a negative feedback loop, controlling the fraction of free barbed ends within ε, that maintained a constant SS protrusion rate V_{PM} regardless of load. Under standard parameter values, the system converged from any ICs to a SS V_{PM} /V_{free} of 0.38 and a free barbed end count N_{fb} of 200 per μm. The filament pattern also selforganized, with only a ±35 degree pattern remaining at this terminal V_{PM} . The system in fact displayed the hallmark adaptation to this pattern from distributions initially devoid of ±35 degree filaments. Any deviation from the symmetrical orientation increased the capping rate at the larger angle more than it increased the branching rate at the smaller angle. Changes in the direction of protrusion were thus corrected for by preferential reproduction of newly symmetrical patterns generated via branching angle “errors.”
Alternate parameter values that resulted in a SS V_{PM} /V_{free} < 0.20 (e.g., very low R_{br} ) resulted in a −70/0/70 degree pattern. Altering V_{free} alone had no effect. Protection from capping within ε was absolutely required for realistic pattern formation, although branching localization was not. The pattern was robust with respect to PM surface and bending energies, but sensitive to filament bending lengths longer than ≈50 nm. These robustness and sensitivity traits describe a selforganizing filament network that resists environmental pressures in maintaining a characteristic orientation pattern and protrusion velocity.
Acknowledgments
We thank the reviewers for valuable comments. This work was supported by a National Institutes of Health Grant GM 62431 (to G.G.B.) and a Northwestern University Pulmonary and Critical Care Division grant (to T.E.S.).
Footnotes
 ^{†}To whom correspondence should be addressed. Email: etaylor3{at}northwestern.edu

Author contributions: T.E.S. and G.G.B. designed research; T.E.S. performed research; T.E.S., E.W.T., and G.G.B. analyzed data; and T.E.S. and E.W.T. wrote the paper.

The authors declare no conflict of interest.

This article contains supporting information online at www.pnas.org/cgi/content/full/0701943104/DC1.
 Abbreviations:
 IC,
 initial condition;
 SS,
 steady state;
 PM,
 plasma membrane;
 LE,
 leading edge;
 EM,
 electron micrograph;
 RF,
 rigid and flexible (beambending) filament model;
 fil,
 filament.
 © 2007 by The National Academy of Sciences of the USA
References

↵
 Miyata H ,
 Hotani H
 ↵
 ↵

↵
 Mullins RD ,
 Heuser JA ,
 Pollard TD
 ↵

↵
 Carlier MF ,
 Laurent V ,
 Santolini J ,
 Melki R ,
 Didry D ,
 Xia GX ,
 Hong Y ,
 Chua NH ,
 Pantaloni D

↵
 Blanchoin L ,
 Pollard TD
 ↵
 ↵
 ↵

↵
 Maly IV ,
 Borisy GG

↵
 Verkhovsky AB ,
 Chaga OY ,
 Schaub S ,
 Svitkina TM ,
 Meister JJ ,
 Borisy GG
 ↵

↵
 Alberts JB ,
 Odell GM
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Boal D

↵
 Isambert H ,
 Venier P ,
 Maggs AC ,
 Fattoum A ,
 Kassab R ,
 Pantaloni D ,
 Carlier MF

↵
 Pollard TD
 ↵
 ↵

↵
 Schafer DA ,
 Jennings PB ,
 Cooper JA
 ↵
 ↵