A Two-Step Biopolymer Nucleation Model Shows a Nonequilibrium Critical Point
Alexander I. P. Taylor, Lianne D. Gahan, Rosemary A. Staniforth, Buddhapriya Chakrabarti
AA Two-Step Biopolymer Nucleation Model Shows a Nonequilibrium Critical Point
Alexander I. P. Taylor, ∗ Lianne D. Gahan,
1, 2
Rosemary A. Staniforth, † and Buddhapriya Chakrabarti ‡ Department of Molecular Biology and Biotechnology, University of Sheffield, Sheffield S10 2TN, UK. Department of Physics and Astronomy, University of Sheffield, Sheffield S3 7RH, UK. (Dated: December 18, 2019)Biopolymer self-assembly pathways are central to biological activity, but are complicated by theability of the monomeric subunits of biopolymers to adopt different conformational states. As a re-sult, biopolymer nucleation often involves a two-step mechanism where the monomers first condenseto form a metastable intermediate, and this then converts to a stable polymer by conformational re-arrangement of its constituent monomers. While existing mathematical models neglect the dynamicsby which intermediates convert to stable polymers, experiments and simulations show that thesedynamics frequently occur on comparable timescales to condensation of intermediates and growth ofmature polymers, and thus cannot be ignored. Moreover, nucleation intermediates are responsiblefor cell toxicity in numerous diseases, such as Alzheimer’s, Parkinson’s, and prion diseases. Dueto the relationship between conformation and biological function, the slow conversion dynamics ofthese species will strongly affect their toxicity. In this study, we present a modified Oosawa modelwhich explicitly accounts for simultaneous assembly and conversion. To describe the conversiondynamics, we propose an experimentally motivated initiation-propagation (IP) mechanism in whichthe stable phase arises locally within the intermediate, and then spreads through additional con-version events induced by nearest-neighbor interactions, in a manner analogous to one-dimensionalGlauber dynamics. Our mathematical analysis shows that the competing timescales of assemblyand conversion result in a nonequilibrium critical point, separating a regime where intermediatesare kinetically unstable from one where conformationally mixed intermediates can accumulate. Inturn, this strongly affects the rate at which the stable biopolymer phase accumulates. In contrastto existing models, this model is able to reproduce commonly observed experimental phenomenasuch as the formation of mixed intermediates, and abrupt changes in the scaling exponent γ , whichrelates the total monomer concentration to the accumulation rate of the stable phase. Our workprovides the first step towards a general model of two-step biopolymer nucleation, which can be usedto quantitatively predict the concentration and composition of biologically crucial intermediates. I. INTRODUCTION
Biopolymer formation is essential to life. Uncon-trolled self-association of biological molecules to formnon-covalent polymers is implicated in diseases such asAlzheimer’s [1], Parkinson’s [2], and sickle-cell anemia [3].In many cases, the monomeric subunits of biopolymersadopt distinct conformational states. As a result, nucle-ation often proceeds via a two-step mechanism, in whichthe monomers first condense into a metastable interme-diate, and the stable polymer develops by conformationalconversion of the monomers within this phase [4–10]. In-termediates have diverse morphologies (Fig. 1(a)), de-pending on the monomer interactions which drive con-densation. Polymeric intermediates are commonly ob-served in experimental and computational studies ofamyloid formation [4, 11–13], and convert to stable poly-mers via reorganization of their secondary or tertiarystructure. In addition, diverse biopolymers includingactin and prions are known to nucleate in spheroidalintermediates [5–7, 10, 14–16]. Biopolymer nucleationintermediates are implicated in cell death and resulting ∗ ataylor8@sheffield.ac.uk † r.a.staniforth@sheffield.ac.uk ‡ b.chakrabarti@sheffield.ac.uk disease [17, 18], and an understanding of their formationis crucial for therapeutic development.In order to accurately predict the rate and mecha-nism of two-step biopolymer nucleation, theoretical mod-els must realistically represent the process by which in-termediates convert to stable polymers. Existing ana-lytical models treat conversion as an instantaneous pro-cess [19–23], neglecting the dynamics by which the sta-ble polymeric phase spreads through the intermediate.However, there is no general reason why the orderingtransition should be sufficiently rapid, or cooperative, forthis to be valid. In contrast, experimental and computa-tional studies often point to progressive ordering on slowtimescales [7, 15, 24–27], indicating that the stable phasearises locally and spreads through the intermediate. Thismeans that condensation and ordering can occur concur-rently [15, 24–26]. In addition, simulations support anautocatalytic mechanism by which conformationally or-dered monomers promote conversion of their neighbors[24, 28, 29]. This raises the possibility that condensationand ordering may be kinetically coupled, with larger in-termediates converting at a higher rate [20, 22]; however,the kinetic prerequisites for such a mechanism have notbeen rigorously determined. Thus, in order to describethe time-evolution of intermediate populations and elu-cidate the role of autocatalysis in two-step nucleation,there is an urgent requirement for a model which ex-plicitly considers the conversion of intermediates by an a r X i v : . [ q - b i o . B M ] D ec autocatalytic initiation-propagation (IP) mechanism.In this paper, we address these issues. We developa minimal model of two-step nucleation (Fig. 1(b)-(f))in which polymeric intermediates undergo an autocat-alytic conversion process analogous to one-dimensional(1D) Glauber dynamics [30]. Through a combinationof analytical calculations and stochastic simulations [31],we show that this simple addition to Oosawa’s classi-cal model of nucleated polymerization [32, 33] producesrich, nonlinear dynamics. In autocatalytic systems, anonequilibrium critical point arises, separating a regimewhere intermediates do not accumulate from one wherethe two-step nature of nucleation becomes apparent. Inthe vicinity of this critical point, nucleation dynamicsresemble a second-order phase transition similar to thetwo-dimensional (2D) Ising model [34, 35]. It is clearfrom our derivation that this phase transition does notrequire a particular dimensionality, and we expect ourresults to be applicable to morphologically diverse in-termediates. Thus, we propose that models of explicitconformational conversion have the capacity to reconcileapparently contradictory interpretations of biopolymernucleation, and resolve important unanswered questionsregarding the origin and nature of metastable intermedi-ates. II. OUTLINE OF THE MODEL
In this section, we outline the mathematical basisof our two-step biopolymer nucleation model. Thismodel combines the established formalism for nucle-ated biopolymer self-assembly [32, 33] with a proposedinitiation-propagation (IP) mechanism of conformationalconversion, in which the stable phase arises locallywithin the intermediate and spreads autocatalyticallyby nearest-neighbor interactions (Fig. 1). In this re-spect, our IP mechanism resembles 1D Glauber dynamics[30]. Although the proposed dynamics affect all polymer-ized monomers, we are particularly interested in theireffects on macroscopic, observable nucleated polymer-ization. Therefore, we concentrate on the aspects ofthese dynamics which are most relevant to the biopoly-mer nucleation and growth, and leave a more comprehen-sive treatment of the internal composition of developingbiopolymers as a subject for future work.
A. Microscopic processes
We consider nucleated polymerization of a freemonomer population maintained at a fixed concentration m , far from equilibrium. This situation arises physiolog-ically, when there is a monomer source, and is also validfor experimental systems on timescales shorter than thetypical timescale τ inf for significant monomer depletionto occur, which we discuss in greater detail in Sec. II. B.Each biopolymer is treated as a 1D lattice of monomers FIG. 1. Schematic of the proposed two-step biopoly-mer nucleation model. (a) Biopolymer nucleation pathwaysfrequently involve diverse metastable intermediates, whichprogress to stable polymers by conversion of monomer sub-units. (b-f) Schematic representation of our model, account-ing for explicit growth and conversion of linear intermediatesas described in Sec. II. A. Color scheme for monomers: green,unassembled; red, metastable (type-1); blue, stable (type-2). with length j . A polymerized monomer is in either themetastable (type-1) or stable (type-2) state (Fig. 1(b)-(f)). Polymerization is initiated at a rate v n = k n m j ,and each polymer initially contains j ≥ k c , or k c + k p if adjacent to a type-2 monomer (Fig. 1(e) & (f)).Thus, k p represents the excess rate due to autocatalysis.For the purpose of this study, it is not necessary to con-sider cases where a type-1 monomer is bounded by twotype-2 monomers, as such cases do not affect the macro-scopic kinetics of nucleated polymerization.Each polymer can elongate by reversible monomer ad-dition at a single ‘active’ end (Fig. 1(c) & (d)). Strongbias towards elongation at a single end is often observedexperimentally [36–38], and is caused by the fact thatmost biopolymers exhibit structural differences betweentheir ends, which result in contrasting physicochemicalproperties. Although inclusion of this bias simplifies ouranalysis, it is not expected to qualitatively alter the re-sults. Therefore, we fully expect our theory to generalizeto more unusual cases [39] where both ends elongate atsimilar rates. It is also important to note that the for-ward and reverse elongation rates must exhibit an equalbias, in order to satisfy detailed balance when m is atthe solubility limit. From this point onwards, we sim-ply refer to the forward process as elongation, and thereverse process as dissociation.During elongation, the conformational state of themonomer being incorporated into the polymer is tem-plated by that of the monomer already at the active end[5, 40] (Fig. 1(c) & (d)). Therefore, if the monomer atthe active end is in the type-1 state, then elongationproceeds via addition of type-1 monomers with a for-ward rate k +1 m and a reverse rate k − . Similarly, if themonomer at the active end is in the type-2 state, thentype-2 monomers are added with forward and reverserates k +2 m and k − , respectively. Thus, polymers ex-hibit one of two distinct growth modes, type-1 and type-2, corresponding to the conformation of the monomer atthe active end and the resulting mechanism of elonga-tion. Growth mode switching occurs when the monomerat the active end converts, causing the emergence of anelongating type-2 phase. Experimental results typicallysupport k +2 (cid:29) k +1 [4, 41], meaning that the switchingrate determines the lifetime of intermediates with a hightype-1 structure content, as well as the observed macro-scopic kinetics of nucleated polymerization. B. Characteristic timescales
Before deriving master equations to predict the macro-scopic kinetics of nucleated polymerization, it is informa-tive to consider the key conformational transitions that anascent polymer undergoes, and the timescales on whichthey occur (Fig. 2). There are three key events that leadto emergence of a type-2 growth mode: (i) initiation ofthe polymer at a time t , with all constituent monomersinitially in the type-1 state; (ii) non-autocatalytic con-version of the first polymerized monomer at time t c >t , causing the formerly ‘pure’ intermediate to becomea mixed intermediate; and (iii) growth mode switch-ing at time t s ≥ t c , which coincides with conversionof the monomer at the active end (Fig. 2(a)). If thefirst monomer to convert is situated at the active end,then t s = t c ; otherwise, additional conversion eventsare required for growth mode switching to occur, sothat t s > t c . If k p (cid:28) k c , these additional conversionevents will be primarily non-autocatalytic and uncorre-lated with one another; if k p (cid:29) k c , conversion of the firstmonomer will stimulate autocatalytic conversion of itsneighbours, causing the type-2 state to propagate fromits point of origin to the active end. We also identifythree important time intervals: (i) the lifetime of a pureintermediate, τ c = t c − t ; (ii) the propagation time τ p = t s − t c , which is the time for a mixed intermedi-ate to adopt a type-2 growth mode; and (iii) the totalswitching time τ s = t s − t = τ c + τ p , which is the totaltime taken for an intermediate to acquire a type-2 growthmode (Fig. 2(a)). In addition, we introduce a single do-main wall propagation time τ (cid:48) p , which is distinct from τ p and will be discussed in greater detail in Sec. II. C.We also consider two timescales related to the overallself-assembly of populations of mature polymers: (i) thecharacteristic self-assembly time τ char , and (ii) the inflec-tion time τ inf (Fig. 2(b)). Both are common experimen- FIG. 2. Characteristic timescales of two-step biopolymernucleation. (a) Schematic showing an example of how thelength ( j ; black line) and AEM domain size ( x ; red line) ofa single polymer can change over time. The three key events(Sec. II. B.) of initiation ( t = t ), conversion of the firstmonomer ( t = t c ), and growth mode switching ( t = t s ) are an-notated with diagrams representing the change in conforma-tional state of the polymer, with the following color scheme forconstituent monomers: green, unassembled; red, metastable(type-1); blue, stable (type-2). (b) Schematic showing themacroscopic self-assembly kinetics of a typical polymer popu-lation (Sec. II. B.), and relevant characteristic times. Curvesrepresent the proportions of monomers which are unassembled( m ( t ); green), or incorporated into polymers with a type-1( M ( t ); red) or type-2 ( M ( t ); blue) growth mode. tal observables which are derived from phenomenologicalanalysis of biopolymer self-assembly curves, and can beused to derive mechanistic conclusions [42, 43]. Specifi-cally, the time-dependent mass of mature polymers fol-lows a convex profile under circumstances where m isapproximately constant. If monomer depletion can oc-cur, the self-assembly curve is sigmoidal and approachesa constant limit due to exhaustion of the free monomer[42] (Fig. 2(b)). For the purpose of this study, we identifymature polymers as species with a type-2 growth mode,and denote the effective concentration of monomer in-corporated into such species as M ( t ). The self-assemblytime τ char is the time required for M ( t ) to reach an arbi-trary, often small proportion of the initial free monomerconcentration m , such that M ( τ char ) = Cm and C isa constant in the interval (0 ,
1) [42] (Fig. 2(b)). Whilea number of different values of C are used, for small C these different characteristic times typically producesimilar scaling behaviors. The inflection time τ inf isthe time at which M ( t ) has its inflection point, suchthat ∂ t M ( τ inf ) = 0 [43] (Fig. 2(b)). Because inflectionis a consequence of monomer depletion, τ inf provides ageneral timescale on which monomer depletion occurs.Therefore, the assumption that m is approximately con-stant implies that τ s (cid:28) τ inf , meaning that M ( t ) is locallyconvex. C. Growth mode switching as an absorbingboundary problem
The ability of type-2 monomers to autocatalyticallyinduce conversion of their neighbors allows the type-2 state to propagate to the active end from a distantsite of origin. Thus, the rate of growth mode switch-ing depends not only on the conformational state of themonomer at the active end, but also the length distribu-tion of the domain of type-1 monomers adjacent to theactive end, which we term the ‘active end metastable’(AEM) domain. When the size x of this domain becomeszero, growth mode switching occurs (Fig. 2(a)). We areprimarily interested in far-from-equilibrium cases where m (cid:29) k − /k +2 , so that growth mode switching has a neg-ligible reverse process and can be treated as an absorb-ing boundary problem. Initially, the polymer is formedwith all monomers in the type-1 state, so that x = j for t ≤ t < t c . Following the first conversion event at t = t c ,the polymer acquires a mixed composition, so that x < j for t ≥ t c . For t > t c , the AEM domain can grow due toelongation at the active end (Fig. 1(c)), or contract dueto either dissociation of monomers from the active end(Fig. 1(d)) or conversion of monomers within the AEMdomain (Fig. 1(e) & (f)). Growth mode switching occurswhen the AEM domain size hits an absorbing boundaryat x = 0, at the time t s .The fact that both autocatalytic and non-autocatalyticconversion events may occur during propagation leads toa distinction between the total propagation time τ p , andthe single domain wall propagation time τ (cid:48) p (Fig. 2(a)).While autocatalytic conversion is restricted to the nearestneighbors of type-2 monomers, non-autocatalytic conver-sion is not. Thus, non-autocatalytic conversion of type-1monomers within the AEM domain may cause a new do-main wall to form, superseding the existing AEM domainwall. While τ p denotes the total time between the onsetof propagation and growth mode switching, τ (cid:48) p denotesthe time taken for a single domain wall to successfullypropagate to the active end, without being superseded.Therefore, τ (cid:48) p only applies to propagation of the last AEMdomain wall to form before growth mode switching. Ifthe original domain wall successfully propagates to theactive end, then τ (cid:48) p = τ p ; otherwise, τ (cid:48) p < τ p . D. General master equations
We consider the multivariate probability mass function(PMF) p j,x ( t ; j , t ), which gives the probability that apolymer will have total length j , AEM domain size x , anda type-1 growth mode at time t . The master equation for p j,x ( t ) has the following form: dp j,x ( t ) dt = δ j,j δ x,j δ ( t − t )+ k c j (cid:88) y = x +1 p j,y ( t ) − k c xp j,x ( t )+ k p [ θ ( j - x -2) p j,x +1 ( t ) − θ ( j - x -1) p j,x ( t )]+ k +1 m [ p j − ,x − ( t ) − p j,x ( t )]+ k − [ p j +1 ,x +1 ( t ) − p j,x ( t )] , (1)where δ j,j and δ ( t − t ) are the Kronecker and Diracdeltas, respectively, and θ ( z ) is a step function such that θ ( z ) = 1 if z ≥ θ ( z ) = 0 if z <
0. As with previ-ous models [43–45], j is taken to be a minimum polymerlength below which rapid disassembly almost always oc-curs; therefore, an absorbing boundary exists such that p j − ,x ( t ) = 0. In addition, to account for irreversiblegrowth mode switching, we impose the absorbing bound-ary condition p j, ( t ) = 0. While the first term in Eq. (1)accounts for initiation at t = t , the terms proportionalto k c , k p , k +1 m , and k − account for non-autocatalyticconversion, autocatalytic conversion, elongation, and dis-sociation, respectively.The terms accounting for non-autocatalytic conversionintroduce a jump process to Eq. (1), which complicatesefforts to satisfy the boundary conditions. Nonetheless,analytical solutions can be obtained in specific cases. Inthe case where x = j , corresponding to the length distri-bution of intermediates with a purely type-1 composition,the summation disappears and the problem becomes sig-nificantly more tractable. We use this special case toobtain (cid:104) τ c (cid:105) in Sec. III. B. For mixed intermediates, thetime-dependence of AEM domain size can be solved bysumming over the polymer length. Let φ x ( t ; x c , t c ) rep-resent the probability that a mixed intermediate has anAEM domain of size x and a type-1 growth mode at time t , given the boundary condition that x = x c when t = t c .We can determine φ x ( t ) by summing p j,x ( t ) over j : φ x ( t ) = ∞ (cid:88) j = j ∨ x +1 p j,x ( t ) , (2)where j ∨ x + 1 denotes the maximum of j and x + 1.The corresponding master equation is as follows: dφ x ( t ) dt = δ x,x c δ ( t − t c ) + k c ∞ (cid:88) y = x +1 φ y ( t ) − k c xφ x ( t )+ k +1 m [ φ x − ( t ) − φ x ( t )]+ ( k − + k p )[ φ x +1 ( t ) − φ x ( t )] . (3)As with related theories [42–44], we consider the casewhere k − is small, allowing losses due to destabilizationof polymers with length j < j to be neglected. The ab-sorbing boundary at x = 0 corresponds to growth modeswitching, and is retained. Eq. (3) is the discrete analogof a jump-diffusion equation: ∂φ ( x, t ) ∂t = δ ( x − x c ) δ ( t − t c ) + k c (cid:90) ∞ x φ ( y, t ) dy − (cid:16) k c x + µ ∂∂x − D ∂ ∂x (cid:17) φ ( x, t ) , (4)where µ = k +1 m − k − − k p and D = ( k +1 m + k − + k p ) / D ≥ − D ≤ µ ≤ D . Because Eq. (3)can be written in an equivalent form to Eq. (4), µ and D are also informative in a discrete context. We definean effective P´eclet number x c | µ | /D for propagation ofthe type-2 state to the active end. Propagation is pri-marily diffusive when x c | µ | /D <
1, and advective when x c | µ | /D >
1. In the advective case, the AEM domaintends to contract if µ <
0, and expand if µ >
0. Thus, µ is the expected rate of AEM domain expansion due tothe competing effects of type-1 growth and type-2 prop-agation, and this competition provides the basis of thecritical behavior we observe in this paper. III. CRITICAL PHENOMENA IN GROWTHMODE SWITCHING
In this section, we investigate how the competing ef-fects of propagation and elongation on AEM domain sizealter the lifetime of intermediates. Using a combina-tion of stochastic simulations and analytical theory, weuncover the existence of a nonequilibrium critical pointresulting from transcritical bifurcation of an order pa-rameter P e around the point µ = 0. This critical pointshares many common features with equilibrium second-order phase transitions, including a jump in ∂ µ P e , criticalslowing, and diverging susceptibility to autocatalytic ef-fects. These effects will be shown in the following sectionto significantly affect the macroscopic kinetics of two-stepbiopolymer self-assembly. A. Existence of a nonequilibrium critical point
To evaluate how the competing effects of propagationand elongation on AEM domain size alter the lifetimeof intermediates, we calculate the mean switching time (cid:104) τ s (cid:105) . For this purpose, we use a rejection-free kineticMonte Carlo (rfKMC; see Supplemental Material at [31])algorithm to determine hitting times at the boundary x = 0 in Eq. (1). The rfKMC algorithm converges to theexact solution for the dynamics as the number of poly-mers N → ∞ [46–48], and its application to biochemi-cal problems including biopolymer self-assembly is well- established [47, 49, 50]. In this case, the rfKMC algo-rithm is preferable to deterministic methods as it samplesbroad polymer size distributions more efficiently. For thepurpose of the simulations, we set k − = 0 so that thereis a single absorbing boundary, although our subsequentmathematical analysis (Sec. III. C. & D.; [31]) is moregenerally applicable to small, non-zero k − . In addition,while we have carried out simulations for diverse j , theresults are qualitatively similar in all cases [31]; therefore,we focus on the j = 2 case as a representative example.To explore the system’s limiting behavior in a 2D space,it is convenient to normalize both the timescales and rateconstants with respect to k c . Therefore, we use relativetime k c t and sample points in ( k p /k c , k +1 m /k c ) space byvarying a pair of bounded, mechanistically significant pa-rameters, µ/D and ζ = k p / ( k c + k p ). While µ/D definesthe tendency for expansion or contraction of the AEMdomain (Sec. II. D.), ζ is the proportion of the conver-sion rate adjacent to a domain wall which is attributableto autocatalytic effects. Thus, ζ provides a measure ofthe importance of autocatalysis in conformational con-version of the assembled monomers.In Fig. 3(a), we show how the normalized mean growthmode switching time k c (cid:104) τ s (cid:105) varies as a function of themodel parameters ( µ/D, ζ ). We use these data to plota nonequilibrium phase diagram, which we divide intothree different dynamical regimes ([31]; Fig. 3(c)): (A)random, (B) cooperative, and (C) competitive. In therandom switching regime, elongation and autocatalysisdo not occur, so switching depends entirely on non-autocatalytic conversion of the monomer at the activeend. In the cooperative regime ( ζ → µ/D (cid:47) k c j . This means that massaccumulation and autocatalytic conversion cooperate tobring about a change in growth mode. In the compet-itive regime, which arises for µ/D (cid:39)
0, the effects ofelongation on AEM domain size outweigh those of au-tocatalytic conversion. Therefore, the AEM domain un-dergoes net expansion, which inhibits propagation of thetype-2 state to the active end. This causes an increasein k c (cid:104) τ s (cid:105) , allowing accumulation of large, metastable in-termediates with a mixed composition. In this regime,mass accumulation and autocatalytic conversion can besaid to compete with one another. Although regimes A-C are continuous with one another, we observe a singularpoint at ( µ/D, ζ ) = (0 , (cid:104) τ s (cid:105) with respect to µ/D diverges in a manner resemblinga second-order phase transition (Fig. 3(a) & (b), 5(a)). FIG. 3. Two-step biopolymer nucleation exhibits distinct dynamical regimes. (a) Nonequilibrium phase diagram superimposedon a heat map of mean growth mode switching time k c (cid:104) τ s (cid:105) (color scale). The solid and dashed lines are contours which separatethree dynamical regimes (A, random switching; B, cooperative; C, competitive). At the critical point (CP), ∂ (cid:104) τ s (cid:105) /∂ ( µ/D )undergoes a jump. (b) A close-up of the heat map of k c (cid:104) τ s (cid:105) at high ζ values (0 . ≤ ζ < . µ/D or ζ close to thecritical point causes a rapid change in the switching time. (c) Schematics of regimes A-C, with color representing monomerconformation: green, unassembled; red, metastable (type-1); blue, stable (type-2). B. Lifetime of pure intermediates
The switching time (cid:104) τ s (cid:105) must reflect a similar criti-cal point in the lifetime of pure intermediates (cid:104) τ c (cid:105) , orthe propagation time of mixed intermediates (cid:104) τ p (cid:105) , since (cid:104) τ s (cid:105) = (cid:104) τ c (cid:105) + (cid:104) τ p (cid:105) . The fact that the critical point oc-curs on the µ = 0 line, along which the competing effectsof autocatalytic conversion and elongation on AEM do-main size are balanced, indicates that an effect on prop-agation is responsible for the observed criticality. Toascertain whether (cid:104) τ c (cid:105) also plays a role in the criticalpoint, we solve Eq. (1) in the case x = j to obtaina moment-generating function (MGF) for τ c . As oursimulations were carried out in the k − = 0 case, weapply the same condition when determining (cid:104) τ c (cid:105) . TheMGF of the first-conversion time is given by (cid:104) e − sτ c (cid:105) =1 − s (cid:80) ˆ p j,j ( s ), where ˆ p j,j ( s ) = (cid:82) ∞ t p j,j ( t ) e − s ( t − t ) dt isthe Laplace transform of p j,j ( t ). By taking the Laplacetransform of Eq. (1) in the case x = j and rearranging,we obtain the following result [31]: (cid:104) e − sτ c (cid:105) = 1 − sk c e k +1 m kc (cid:16) k +1 m k c (cid:17) − (cid:0) s + k +1 m kc + j (cid:1) × γ (cid:16) s + k +1 m k c + j , k +1 m k c (cid:17) , (5) where γ is a lower incomplete gamma function. Thus, (cid:104) τ c (cid:105) = 1 k c e k +1 m kc (cid:16) k +1 m k c (cid:17) − (cid:0) k +1 m kc + j (cid:1) × γ (cid:16) k +1 m k c + j , k +1 m k c (cid:17) . (6)Eq. (6) decreases monotonically with k +1 m /k c , and hasthe limits (cid:104) τ c (cid:105) → ( k c j ) − as k +1 m /k c →
0, and (cid:104) τ c (cid:105) → k +1 m /k c → ∞ . This behavior reflects the fact thatelongation increases the number of monomers within thepolymer which can non-autocatalytically convert, reduc-ing the time until the first monomer converts. Thus,if growth does not occur ( k +1 m /k c → k c j . Conversely, if growthis rapid, then the length goes to infinity and the timeuntil the first conversion event becomes negligible. Im-portantly, Eq. (6) does not predict any form of critical-ity, and the fact that it does not depend on k c meansthat it cannot be expressed as a function of µ/D and ζ alone. This prediction is confirmed by analysis of k c (cid:104) τ c (cid:105) values from our rfKMC trajectories, which we present inFig. 4(a). Given that k +1 m /k c = [ ζ (2+ ξ )] / [(1 − ζ )(2 − ξ )],the phase behavior of k c (cid:104) τ c (cid:105) in Fig. 4(a) is accounted forby Eq. (6) alone, with no critical point. Thus, changes in (cid:104) τ p (cid:105) rather than (cid:104) τ c (cid:105) must be responsible for the observedcriticality. FIG. 4. The critical point in (cid:104) τ s (cid:105) = (cid:104) τ c (cid:105) + (cid:104) τ p (cid:105) is caused by a similar critical point in (cid:104) τ p (cid:105) , but not (cid:104) τ c (cid:105) . (a) Heat map of k c (cid:104) τ c (cid:105) ,the normalized mean lifetime of pure intermediates. The phase behavior is described by Eq. (6), and does not exhibit a criticalpoint at ( µ/D, ζ ) = (0 , k c (cid:104) τ p (cid:105) , the normalized mean time for mixed intermediates to undergo growth modeswitching. Mixed intermediates are kinetically unstable ( k c (cid:104) τ p (cid:105) = 0) along the line between ( µ/D, ζ ) = ( − ,
1) and (0 , ∂ (cid:104) τ p (cid:105) /∂ ( µ/D ) = 0 along this line as well. A second-order critical point (CP) occurs at ( µ/D, ζ ) = (0 , ∂ (cid:104) τ p (cid:105) /∂ ( µ/D ) jumps to a non-zero value. C. Nature of the critical point
The critical point in (cid:104) τ s (cid:105) must be caused by an abruptchange in the manner in which the type-2 state propa-gates to the active end. This deduction is supported byanalysis of k c (cid:104) τ p (cid:105) values from our rfKMC trajectories,which we present in Fig. 4(b). While k c (cid:104) τ c (cid:105) is unaffectedby the critical point, k c (cid:104) τ p (cid:105) exhibits second-order criticalbehavior similar to k c (cid:104) τ s (cid:105) . Specifically, mixed interme-diates are kinetically unstable ( k c (cid:104) τ p (cid:105) = 0) along a linefrom ( µ/D, ζ ) = ( − ,
1) to (0 , k c (cid:104) τ p (cid:105) is constant alongthis line, it follows that ∂ (cid:104) τ p (cid:105) /∂ ( µ/D ) = 0 along this lineas well. At the critical point at ( µ/D, ζ ) = (0 , ∂ (cid:104) τ p (cid:105) /∂ ( µ/D ) = 0, allowing k c (cid:104) τ p (cid:105) totake nonzero values to the right of this point. To un-derstand this result, let us consider the propagation ofa single AEM domain wall from its site of origin to theactive end. For this purpose, we solve Eq. (3) in thecase where k c x (cid:28) k p , so that the terms proportional to k c become negligible. This describes the propagation ofa single domain wall when the probability of additionalnon-autocatalytic conversion events is negligible ( ζ → x remains finite. We obtain φ x ( t ; x c , t c ) using latticeFourier transforms and the image method (see Supple-mental Material at [31]): φ x ( t ) = e − D ( t − t c ) (cid:16) D + µ D − µ (cid:17) x − xc × (cid:8) I | x − x c | [ ν ( t − t c )] − I | x + x c | [ ν ( t − t c )] (cid:9) , (7)for x ≥ t ≥ t c , where I x ( t ) is a modified Besselfunction of the first kind and ν = (cid:112) D − µ . The switching rate, which is equal to the hitting rate atthe boundary x = 0, is given by h ( t ; x c , t c ) = ( k − + k p ) φ ( t ; x c , t c ). We define an escape probability P e ( x c ) =lim t →∞ (cid:80) x φ x ( t ; x c , t c ), which is the probability that agiven AEM domain wall does not reach the active endas t → ∞ . In the ζ → x goes to infinity. In this case, additional non-autocatalytic conversion events become non-negligible,allowing propagation to be interrupted by formation of anew domain wall closer to the active end. Thus, P e ( x c )is the probability that the AEM domain does not dis-appear without additional non-autocatalytic conversionevents occurring. From final value theorem, we find that P e ( x c ) = 1 − lim s → + ˆ h ( s ; x c ). The Laplace transformˆ h ( s ; x c ) = (cid:82) ∞ t c h ( t ; x c , t c ) exp[ − s ( t − t c )] dt has the follow-ing form [31]: ˆ h ( s ; x c ) = (cid:16) s + 2 D − r ± D + µ (cid:17) x c , (8)where r ± = ± (cid:112) s + 4 Ds + µ . As shown in Fig. 5(b),Eq. (8) delineates a parabola-like curve, which opens tothe right and whose upper ( r − , unphysical) and lower( r + , physical) branches meet at the point ( s ∗ , ˆ h ∗ ), where s ∗ = ν − D ≤
0. When µ (cid:54) = 0, both branches interceptthe ˆ h axis; when µ = 0, ( s ∗ , ˆ h ∗ ) touches the ˆ h axis and thebranches exchange intercepts, causing an abrupt changein ∂ µ P e ( x c ): P e ( x c ) = (cid:40) , µ ≤ , − (cid:0) D − µ D + µ (cid:1) x c , µ > . (9)In the ζ → k c (cid:104) τ s (cid:105) ∼ lim j →∞ (cid:8) (cid:80) j − x c =0 [1 − P e ( x c )] (cid:9) − = P e (1) [31]. Thus, Eq. (9) explains thecritical point in k c (cid:104) τ s (cid:105) , and the agreement between ourcalculated scaling law and the simulated result is shownin Fig. 5(a). The emergence of the critical point in the ζ → H [35]. Therefore, we propose that theappearance of nonzero P e ( x c ) for µ > P e ( x c ) is an orderparameter and ( µ, ζ ) play an analogous role to ( T, H ) inthe Ising model. This transition is associated with theemergence of a regime where mixed intermediates havea finite lifetime, and their continued formation allows anonequilibrium steady-state population to persist.
D. Diverging propagation times and susceptibilityat the critical point
As further evidence for this critical point, we observethat the moments of the propagation time diverge when( µ/D, ζ ) = (0 ,
1) (Fig. 5(c)). As discussed in Sec. II. C., τ (cid:48) p is the time taken for a single domain wall to success-fully propagate to the active end, without being inter-rupted by formation of a new domain wall closer to theend (Fig. 2(a)). In the ζ → n th -order momentis given by (cid:104) ( τ (cid:48) p ) n (cid:105) = ( − n lim s → + ( ∂ ns ˆ h ) / ˆ h . The upperand lower branches of ˆ h meet at s ∗ = 0 when µ = 0,so the derivatives ∂ ns ˆ h diverge and the moments behavecorrespondingly. More generally,lim ζ → (cid:104) τ (cid:48) p (cid:105) = x c | µ | , (10a)lim ζ → Var( τ (cid:48) p ) = 2 x c D | µ | , (10b)which strongly resembles the diverging magnetic suscep-tibility χ H = ∂ H M near the Curie point. Because τ (cid:48) p is the time interval in which non-autocatalytic conver-sion events can interrupt propagation, a connection to asusceptibility of the form χ ζ ∝ ∂ ζ P e is highly intuitive.To test this connection, we determine analytic boundsfor the susceptibility in the ζ → χ ζ = − ∂ ζ P e = k p ∂ k c P e , where the signof our definition ensures a positive value, but does notaffect the presence or absence of a divergence. This isessentially the effect of ζ on the probability that prop-agation will be interrupted by a non-autocatalytic con-version event in the AEM domain. By differentiating Eq.(3) with respect to k c , we obtain a master equation whichcan be used to determine ∂ k c φ ( t, k c ; x c , t c ) [31]. In turn,this can be used to assess the impact of k c on the hittingrate and P e . We solve this master equation to obtainupper and lower bounds for ∂ k c φ ( t, k c ; x c , t c ): χ ζ = C ( µ, D, x c ) k p (1 − P e )2 (cid:104) τ (cid:48) p (cid:105) , (11a)1 ≤ C ( µ, D, x c ) ≤ x c (cid:18) D | µ | (cid:19) + (cid:18) D | µ | (cid:19) . (11b) Thus, there is a divergence in χ ζ at the critical pointwhich has a scaling law between | µ | − and | µ | − and isclosely related to the divergence of (cid:104) τ (cid:48) p (cid:105) and Var( τ (cid:48) p ). E. Relationship to a transcritical bifurcation
The exchange of intercepts by the branches of ˆ h ( s ; x c )closely resembles a transcritical bifurcation (Fig. 5(b)),raising the question of whether P e ( x c ) can be written as afixed point of the corresponding normal form. We beginby noting that the escape probability can be expressedas a limit P e ( x c ) = lim t →∞ Φ( t ; x c , t c ) of the probabilityΦ( t ; x c , t c ) = (cid:80) x φ x ( t ; x c , t c ) that the AEM domain wallhas not yet propagated to the active end, under condi-tions where k c is sufficiently small that additional non-autocatalytic conversion events are negligible. Because P e ( x c ) = 1 − [1 − P e (1)] x c , a bifurcation in P e (1) mustcause a similar behavior in P e ( x c ) [31]; therefore, we fo-cus on the case where x c = 1. Using the conservationlaw ∂ t Φ( t ; x c , t c ) = δ ( t − t ) − h ( t ; x c , t c ), we re-expressEq. (8) in the form: s ˆΦ( s ; 1) − s = µs ˆΦ( s ; 1) − k +1 m s ˆΦ( s ; 1) . (12)Either by applying final value theorem and the powerrule of limits [51] to Eq. (12), or by taking the inverseLaplace transform and evaluating the convolution thatcorresponds to the ˆΦ( s ; 1) term, it is possible to obtainthe following result in the t → ∞ limit [31]: d Φ( t ; 1 , t c ) dt ≈ µ Φ( t ; 1 , t c ) − k +1 m Φ( t ; 1 , t c ) . (13)This is the normal form of a transcritical bifurcation. Thephysical and non-physical branches of Eq. (8) correspondto the regions of Eq. (13) where ∂ t Φ( t ; 1 , t c ) decreasesor increases with Φ( t ; 1 , t c ), respectively, and the ex-change of intercepts between these branches correspondsto an exchange of stability between the fixed points ofΦ( t ; 1 , t c ). Solutions to Eq. (13) provide insights into thediverging mean and variance of τ (cid:48) p at the critical point,as the hitting rate can be derived from Φ( t ; 1 , t c ) accord-ing to the formula h ( t ; 1 , x c ) = − ∂ t Φ( t ; 1 , t c ) [31]. When µ (cid:54) = 0, the solution is a logistic function, so that thehitting rate is exponentially bounded. However, when µ = 0, the solution is ∼ ( k +1 m t ) − , meaning that thehitting time distribution becomes heavy-tailed and itsmoments diverge. This is an example of critical slow-ing, and is caused by the disappearance of the linearterm in ∂ t Φ( t ; 1 , t c ). In equilibrium systems, a relation-ship between critical slowing and diverging susceptibilitycan be established by considering the Landau potential V (Φ), which satisfies the law ∂ t Φ = − ∂ Φ V (Φ). Boththe relaxation time following a small variation, and thesusceptibility of Φ to a linearly coupled external field,are inversely proportional to ∂ V (Φ). At the criticalpoint, disappearance of the linear term in ∂ t Φ means that ∂ V (Φ) = 0, so that both the relaxation time and suscep-tibility diverge. In our system, the connection between FIG. 5. A continuous phase transition arises in the limit ζ = k p / ( k c + k p ) →
1. (a) Dependence of k c (cid:104) τ s (cid:105) on µ/D for varying k p . Color scale: solid black, k p = 0; gray, k p = k c ; purple, k p = 10 k c ; blue, k p = 10 k c ; green, k p = 10 k c ; amber, k p = 10 k c ;red k p = 10 k c ; dashed black, asymptotic scaling law from Eq. (9), showing emergence of a critical point at µ/D = 0. (b)Origin of the criticality, showing an exchange of intercepts by the physical and non-physical branches of ˆ h ( s ; x c ) in the ζ → µ = − D ; black, µ = 0; red, µ = D . (c) Divergence of the mean ( k p (cid:104) τ (cid:48) p (cid:105) /x c ; dashed) and variance( k p Var( τ (cid:48) p ) /x c ; solid) of the single domain wall propagation time around the critical point in the ζ → k p and x c . critical slowing and the susceptibility is more direct, asthe mean propagation time defines a lower bound for thesusceptibility. Nonetheless, the principle is similar, inthat the exchange of stability by the fixed points causesdisappearance of the linear term in Eq. (13), allowingthe susceptibility to another parameter to diverge.It is possible to obtain an effective potential of theform V (Φ) = − µ Φ / k +1 m Φ / t ; 1 , t c ). Due to the irreversible natureof growth mode switching, Φ( t ; 1 , t c ) can only decrease;thus, the regions of Eq. (13) for which ∂ t Φ( t ; 1 , t c ) > h ( s ; 1)situated at s < IV. EFFECTS OF CRITICALITY ON THEMACROSCOPIC SELF-ASSEMBLY KINETICS
Autocatalytic biopolymer systems with high ζ will ex-hibit a nonlinear response to µ near the critical point.In this section, we evaluate the effects of this critical-ity on the macroscopic self-assembly kinetics. We beginby estimating likely ζ values based on typical biologicalfree energy scales. We find that stabilization of confor-mational transition states by even modest free energydifferences ( < k B T ) can result in highly autocatalytic conformational conversion. Thus, different biopolymersare expected to exhibit the full spectrum of autocatalyticbehaviors predicted by our model. In systems with high ζ , the critical point has a profound effect on the macro-scopic kinetics, inducing rapid changes in the size andabundance of steady-state populations of intermediates,and the concentration-dependent scaling behavior of theformation of mature polymers. These effects are notpredicted by existing local equilibrium models [52–54],and indicate that nonequilibrium criticality may explainwidespread experimental phenomena [55, 56] which one-step nucleation models have not been able to address. A. Estimation of ζ from Kramers theory Critical behavior arises in the autocatalytic ( ζ → k c is non-zero and k p is fi-nite, meaning that a value of ζ = 1 is not strictly pos-sible, ζ can be asymptotically close to this limit. Real-istic ζ values can be estimated from reaction rate the-ory. Let us consider the conversion rate of a single type-1 monomer in the AEM domain, situated n monomersfrom the active end. We represent the conformation ofthe monomer with a single continuous degree of freedom, q n ; however, we note that our analysis can be general-ized to cases where monomers have additional slow de-grees of freedom without changing the conclusions. Thefree energy G ( q n ; q n − , q n +1 ) of a polymerized monomerdepends on its own conformation q n , and a term ac-counting for interactions with the adjacent monomers q n − and q n +1 . The monomer at position n has freeenergy minima corresponding to the type-1 and type-2 states, which are situated at q n = p and q n = p .These minima are separated by a barrier at q n = p ‡ ,which corresponds to a conformational transition state.The fact that conversion has a negligible reverse pro-0cess entails that G ( p ; q n − , q n +1 ) (cid:29) G ( p ; q n − , q n +1 ),and the effectively two-state nature of our model en-tails that G (cid:48)(cid:48) ( p ; q n − , q n +1 ) and G (cid:48)(cid:48) ( p ; q n − , q n +1 ) areboth large. Therefore, the conformation q n (cid:48) of all non-converting monomers n (cid:48) (cid:54) = n is assumed to be close to p or p . Monomers situated closer to the active end( n (cid:48) < n ) are part of the AEM domain, and must be inthe type-1 state ( q n (cid:48) ≈ p ); monomers situated furtherfrom the active end ( n (cid:48) > n ) may be in either the type-1 or type-2 state ( q n (cid:48) ≈ p or q n (cid:48) ≈ p ). Although itis possible to write a Langer equation [57, 58] for theconversion rate which accounts for continuous conforma-tional variation in the monomers at positions n (cid:48) (cid:54) = n ,the results are hard to interpret. The main reason forthis is that such a model predicts dynamics based on the j -dimensional free energy function (cid:80) n G ( q n ; q n − , q j +1 ),which has 2 j − relevant saddle points. Instead, we sim-plify the analysis by neglecting variation in q n (cid:48) about thepoints p or p . We write the free energy of the monomerat position n as a pair of functions corresponding to thetwo possible conformations of the monomer at position( n + 1): G ( q n ) = G ( q n ; p , p ) , (14a) G ( q n ) = G ( q n ; p , p ) . (14b)The monomer conformation exhibits ovderdampedLangevin dynamics, so that the corresponding conforma-tional probability density functions ρ ( q n , t ) and ρ ( q n , t )satisfy the following Fokker-Planck equation: ∂ρ i ∂t = ∂∂q n D i (cid:20) ∂ρ i ∂q n + βρ i dG i dq n (cid:21) (15)where i = 1 , n + 1), β = 1 /k B T , and D i ( q n ) is the confor-mational diffusion coefficient. The conversion rate can beestimated from G i ( q n ) and D i ( q n ) using Kramers theory[58, 59]. The non-autocatalytic conversion rate is givenby the passage rate at the barrier ( q n , q n +1 ) = ( p ‡ , p ): k c = βD ( p ‡ )2 π (cid:2) G (cid:48)(cid:48) ( p ) | G (cid:48)(cid:48) ( p ‡ ) | (cid:3) / e − β ∆ G ‡ , (16)where ∆ G ‡ = G ( p ‡ ) − G ( p ). Similarly, the autocat-alytic conversion rate is given by the passage rate at thesaddle point ( q n , q n +1 ) = ( p ‡ , p ): k c + k p = βD ( p ‡ )2 π (cid:2) G (cid:48)(cid:48) ( p ) | G (cid:48)(cid:48) ( p ‡ ) | (cid:3) / e − β ∆ G ‡ , (17)where ∆ G ‡ = G ( p ‡ ) − G ( p ). Therefore, the definition ζ = k p / ( k c + k p ) gives the result: ζ = 1 − D ( p ‡ ) (cid:2) G (cid:48)(cid:48) ( p ) | G (cid:48)(cid:48) ( p ‡ ) | (cid:3) / D ( p ‡ ) (cid:2) G (cid:48)(cid:48) ( p ) | G (cid:48)(cid:48) ( p ‡ ) | (cid:3) / e − β ∆∆ G ‡ , (18)where ∆∆ G ‡ = ∆ G ‡ − ∆ G ‡ is the reduction in thefree energy barrier of conversion due to autocatalytic ef-fects. It is not immediately obvious whether the pre-exponential ratio would have a value greater or less than 1. While ordering of the ( n + 1) monomer might be ex-pected to reduce conformational diffusion, the simulta-neous loss of complementarity between the monomers’structures might exert the opposite effect. The effectson G (cid:48)(cid:48) ( p ) and G (cid:48)(cid:48) ( p ‡ ) are also hard to gauge, and maybe biopolymer-specific. In the absence of further infor-mation, we set the value of the pre-exponential ratioto 1, so that ζ = 1 − exp ( − β ∆∆ G ‡ ), and we exam-ine the relationship between ζ and ∆∆ G ‡ . Eq. (18)predicts that even high ζ values will be associated withrelatively modest free energy differences; for example, ζ = 0 . ⇔ k p ≈ k c implies ∆∆ G ‡ = 6 . k B T ,which is well within the range of typical biological freeenergy scales at T = 310 K . For comparison, formationof a single hydrogen bond in protein has a free energychange of ∼ k B T [60, 61], and estimated free energybarriers for protein folding are mostly in the 2 − k B T range [62], accounting for the effects of large numbers ofcompeting, non-additive interactions. In addition, a re-cent experimental study of biopolymer nucleation by theA β peptide indicated that autocatalytic interactions be-tween polymers reduce the nucleation free energy barrierby ∆∆ G ‡ = 19 k B T [63]. Although this form of au-tocatalysis is different from the form that we considerhere, as it depends on interactions between polymersrather than within a polymer, the two processes mayshare a common physical basis. It is interesting to notethat an equivalent ∆∆ G ‡ value in our model would give ζ ≈ − − ⇔ k p ≈ k c , resulting in highly autocat-alytic behavior. Thus, we believe that values of ζ closeto the ζ = 1 limit are highly plausible, meaning thatbiopolymer systems will frequently exhibit high levels ofautocatalysis and resulting critical behavior. B. Nonlinear response to perturbations in theautocatalytic limit
We now investigate the response of the macroscopicself-assembly kinetics to variations in µ/D , when ζ ishigh. To do so, we calculate the time-evolution of theprincipal moments of the polymer length distribution.Let f i,j ( t ) represent the concentration of polymers withgrowth mode i = 1 , j at time t . In parallelwith the notation used by [43, 44], we define the followingprincipal moments: P i ( t ) = ∞ (cid:88) j = j f i,j ( t ) , (19) M i ( t ) = ∞ (cid:88) j = j jf i,j ( t ) . (20)The zeroth-order moment P i ( t ) is the total concentrationof polymers with growth mode i , and acts as an effectivenormalization constant for f i,j ( t ). The first-order mo-ment M i ( t ) is the effective concentration of monomers1incorporated into polymers with growth mode i , or the‘polymer mass-concentration’. The mean polymer length L i ( t ) can be obtained by the normalization: L i ( t ) = P i ( t ) M i ( t ) . (21)Let us begin by examining the effect of the critical pointon M ( t ) and L ( t ). The length distribution f ,j ( t ) of thepopulation of type-1 species contains contributions fromindividual polymers initiated at all t < t . If nucleatedpolymerization begins when t = 0, then f ,j ( t ) = v n j (cid:88) x =1 (cid:90) t p j,x ( t ; t ) dt (22)where v n = k n m j is the nucleation rate. Thus, P ( t ) = v n ∞ (cid:88) j = j j (cid:88) x =1 (cid:90) t p j,x ( t ; t ) dt , (23) M ( t ) = v n ∞ (cid:88) j = j j j (cid:88) x =1 (cid:90) t p j,x ( t ; t ) dt . (24)In Sec. III. A., we obtained numerical solutions for p j,x ( t ; t ) using an rfKMC algorithm. We will now usethese solutions to predict M ( t ) and L ( t ), by applyingEq. (23-24). As we are specifically interested in theresponse of the system to varying µ/D at high ζ , wewill focus on cases where µ/D varies through the range[ − , ζ is fixed at the value ζ = 0 . ⇔ k p ≈ k c . This value is associated with a comparativelysmall ∆∆ G ‡ , and so is physically very plausible (Sec.IV. A.).In Fig. 6, we present the effect of µ/D on the time-evolution of αM ( t ) and L ( t ), where α = k c /v n is anormalization constant facilitating the use of relativetime k c t . The integral (cid:82) ∞ p j,x ( t ; t ) converges, allow-ing both αM ( t ) and L ( t ) to approach steady-statevalues αM ss = α lim t →∞ M ( t ) (Fig. 6(a) & (b)) and L ss = lim t →∞ L ( t ) (Fig. 6(b)) in the t → ∞ limit.These nonequilibrium steady states reflect the concen-trations and size of intermediates likely to be observedin physiological contexts, if the free monomer concentra-tion and chemical environment are stable on timescales (cid:29) τ steady , where τ steady is the approximate timescaletaken to reach the steady state. In experimental contexts,the pre-steady-state behavior is likely to dominate in theself-assembly lag phase if τ inf (cid:28) τ steady , and the steady-state behavior is likely to dominate if τ inf (cid:29) τ steady .We will evaluate the effects of the critical point on boththe pre-steady-state kinetics and the steady-state kinet-ics. The pre-steady-state kinetics appear linear whenplotted on double-logarithmic axes (Fig. 6(a)), indicat-ing that they are described by a power law of the form αM ( t ) ∝ t d . Power laws are often observed in biopoly-mer self-assembly kinetics where t is small, and d can be interpreted as the number of sequential processes con-tributing to biopolymer mass accumulation; thus, nucle-ation without polymerization gives d = 1, and nucleatedpolymerization gives d = 2. To assess whether the pre-steady-state kinetics are affected by the critical point,we calculated d = ∂ ln M ( t ) /∂ ln t | k c t =0 . for varying µ/D (Fig. 6(b)). Although d increases with µ/D , reflect-ing the increasing importance of polymerization in massaccumulation, there is no clear response to the criticalpoint. This is unsurprising, given that the critical pointresults from changes in growth mode switching, which isnegligible in the pre-steady-state.In contrast, the steady-state kinetics are strongly af-fected by the critical point. We begin by examining theeffect of the critical point on τ steady . For the purpose ofthis analysis, we define τ steady such that αM ( τ steady ) =0 . αM ss . As shown in Fig. 6(a), τ steady initially decreaseswith increasing µ/D , as both elongation and propagationcooperate to accelerate growth mode switching. How-ever, this trend reverses close to the critical point, causinga delay in the approach to the steady-state and allow-ing greater accumulation of polymer mass. This delaycorresponds to the emergence of a population of large,conformationally mixed intermediates. As a result, both αM ss and L ss exhibit a highly nonlinear response to µ/D around the critical point (Fig. 6(a) & (b)). In physiolog-ical and experimental contexts, this means that smallvariations in the free monomer concentration or solventcomposition have the potential to induce large changes inthe populations of biopolymer nucleation intermediates.In addition, the total rate of growth mode switching isstrongly dependent on the size and nature of the interme-diate populations, so the increase in τ steady will stronglyaffect the accumulation of mature polymers with a type-2growth mode. C. Low concentration-dependence is a dynamicalsignature of a nonequilibrium critical point
As discussed in Sec. II. B., a common experimentaldescriptor of the macroscopic self-assembly kinetics isthe characteristic time τ char at which M ( τ char ) = Cm ,where C is an arbitrary constant (Fig. 2(b)). Theconcentration-dependence of τ char can be described usinga scaling exponent γ , which has the following definition: γ = − ∂ ln τ char ∂ ln m . (25)In cases where C is sufficiently small (0 < C (cid:47) . γ [42]. For example, the Oosawa model predicts γ = n c / n c is the nucleation order, and the nucleation-conversion-polymerization (NCP) model, which assumesinstantaneous growth mode switching, predicts γ = n c / τ steady on the macroscopic self-assembly kinetics, we pre-dict γ values based on our simulations. We begin by2 FIG. 6. Varying µ close to the critical point ( k p = 10 k c ) results in a highly nonlinear response of the macroscopic kinetics.(a) Time-evolution of the mass of species with a type-1 growth mode, αM ( t ), for varying µ/D . Color scale: red, µ/D = − − .
9; yellow, − .
5; green, −
1; cyan, − .
5; blue, 0 (starred, critical); indigo, 0 .
5; purple, 1; gray, 1 .
5. The dashed lineplots the ( k c t, αM ( t )) at which αM ( t ) = 0 . αM ss , for varying µ/D . (b, left axis) Effect of µ/D on the steady-state mass( αM ss ; solid black) and length ( L ss ; dashed black) of species with a type-1 growth mode. (b, right axis) Effect of µ/D on thescaling exponents d (red) and γ (blue). (c) Dependence of γ on µ/D , for varying ζ . Color scale: red, ζ = 0; amber, ζ = 0 . ζ = 0 .
99; blue, ζ = 0 . ζ = 0 . ζ = 0 . observing that: df ,j ( t ) dt = v n ( k − + k p ) (cid:90) t p j, ( t ; t ) dt + k +2 m [ f ,j − ( t ) − f ,j ( t )] − k − [ f ,j +1 ( t ) − f ,j ( t )] . (26)Our rfKMC simulations were carried out in the k − → P ( t ) = v n t − P ( t ) . (27)To simplify our analysis, let us consider the case where k +2 m (cid:29) M ss /τ s , so that the predominant contribu-tion to the type-2 polymer mass comes from elongation,rather than conversion of intermediates. Therefore, M ( t ) ≈ k +2 m (cid:90) t P ( τ ) dτ. (28)In the pre-steady-state ( t (cid:28) τ steady ), P ( t ) is approx-imately linear, so both P ( t ) and M ( t ) will exhibithigher-order scaling with time. In the steady-state ( t (cid:29) τ steady ), P ( t ) is approximately constant, meaning that P ( t ) and M ( t ) will exhibit linear and quadratic scal-ing with time, respectively. It is more convenient toexpress concentrations in relative units than it is tochoose arbitrary values of the rate constants. There-fore, we calculate α (cid:48) M ( t ) using Eq. (23, 27-28), where α (cid:48) = ( k +1 /k p ) j +1 ( k c /k n k +2 ) is a normalization constant,and set C = 0 . /α (cid:48) . From the definition of µ and D , wefind that ∂ ln (cid:0) µD (cid:1) ∂ ln m = 4 − (cid:0) µD (cid:1) (cid:0) µD (cid:1) . (29)Thus, under circumstances where the rate constants arekept fixed, variation in µ/D can be attributed to changes in the free monomer concentration m . This means thatwe are able to determine γ by calculating τ char for vary-ing, finely spaced µ/D values, estimating the logarithmicderivative with respect to µ/D from the correspondingfinite difference equation, and then multiplying by theright-hand side of Eq. (29). In Fig. 6(c), we present theresults obtained for various ζ . At low ζ , γ monotoni-cally transitions between the limits lim µ →− γ = 1 andlim µ → γ = 2 /
3, which agree with the predictions of [33]and [21] discussed above. However, a minimum in γ ap-pears at high ζ , which starts at the point µ/D = 2 andmoves towards µ/D = 0 as ζ →
1. This depression occursbecause µ/D scales with m , so that increasing m causesa reduction in growth mode switching, inhibiting fur-ther acceleration of type-2 mass accumulation. In highlyautocatalytic systems ( k p > k c ), the effect is severeenough to cause a temporary reversal in concentration-dependence. As the minimum γ decreases, the drop in γ becomes increasingly sharp about the point µ/D = 0.Both the low concentration-dependence ( γ <
0) andsharp drop in γ predicted by our model are commonly ob-served in experimental studies, but cannot be explainedby existing models [52–54]. Furthermore, while exist-ing models invoke local equilibria to explain low γ [52–54], our model predicts a low γ purely as a result of thedynamics. We thus propose that a low concentration-dependence in biopolymer self-assembly may representthe dynamical signature of nonequilibrium critical point. V. DISCUSSION
Biopolymer self-assembly is a fundamental mechanismof biological organization, and plays a key role in dis-ease. Experimental and computational studies have re-vealed that monomers often populate metastable confor-mational states during the process of self-assembling toform a biopolymer [40, 64–66]. As a result, nucleation of3a biopolymer often proceeds via one or more metastableintermediates, whose constituent monomers are confor-mationally distinct from those within the mature poly-mer. In addition, monomers within the same interme-diate often exhibit considerable conformational hetero-geneity. This complexity has challenged our ability toextract detailed mechanistic information at the molecu-lar level. Current mathematical models either neglect theconformational heterogeneity of the assembled monomersentirely, or assume that monomers within the same inter-mediate rapidly assume the same conformational state,so that intermediates with a mixed composition are notobserved. However, mixed intermediates are commonlyobserved in experiments and simulations, indicating thatmathematical models which neglect these dynamics areonly applicable under highly simplified experimental con-ditions [64]. Moreover, this conformational diversity isresponsible for many of the most important biological be-haviors of biopolymers [67–69]. One notable example isthe role of intermediates in the assembly of amyloid fibresfrom protein monomers. The conformational heterogene-ity of amyloid self-assembly intermediates allows them tointeract with a diverse range of other biomolecules [70].This disrupts cellular function and is causative in themajority of degenerative conditions of aging, includingAlzheimer’s and Parkinson’s diseases [16, 18, 71–74].In this paper, we address this issue by introducing aminimally complex process by which the self-assembledmonomers are able to exhibit conformational heterogene-ity, and change conformation over time. This simple ad-dition causes our model to exhibit rich, critical dynam-ics. As a result, the macroscopic self-assembly kinet-ics display features such as a highly nonlinear depen-dence of steady-state intermediate populations on theself-assembly conditions, and loss or even inversion of theconcentration-dependence of the self-assembly rate of thestable phase. These phenomena are commonly observedin experimental studies [55, 56], but are not explained byexisting models. Thus, for the first time, the initiation-propagation (IP) model proposed here allows quantita-tive predictions to be made regarding the concentrationof conformationally heterogeneous intermediates.
A. Competing dynamics in morphologically diverseintermediates
While polymeric intermediates are commonly observedin biopolymer self-assembly, a diverse range of otherintermediate morphologies are also observed, includingsmall prenucleation clusters with an apparently ran-dom organization of constituent monomers [75, 76], largespheroidal intermediates [52, 77, 78], and species possess-ing fractal geometry which is evident at a macroscopicscale [79]. In addition to applying to polymeric interme-diates, we expect the conclusions of our study to gen-eralize to these other diverse geometries. To make thisgeneralization, we note that propagation involves incor- poration of monomers into a growing stable phase byautocatalytic conversion of condensed monomers in theintermediate phase. Thus, the dimensionality of propa-gation is defined by that of the stable phase. In an n -dimensional intermediate, propagation will manifest aslinear growth of the stable phase through the interme-diate, with the growth mode switching when an auto-catalytic end reaches the solvent. This form of prop-agation is well-supported by experimental observationsof spheroidal intermediates [6, 14, 15], which are oftenlarge enough to be observed by light microscopy tech-niques. In these studies, the stable polymeric phase ap-pears to originate within the intermediate by conforma-tional rearrangement of a subset of the monomers, andthen emerge from the surface of the intermediate. Whenthe active end reaches the solvent, the dramatic changein the physicochemical environment, and the availabilityand conformational state of nearby monomers, would beexpected to induce a marked change in the growth behav-ior of the polymeric phase analogous to the switch frompropagation ( k c + k p ) to type-2 growth ( k +2 m ) behaviorin our polymer model.The fact that polymerization is likely to manifest aslinear growth-like behavior in multidimensional interme-diates means that we expect IP models involving suchintermediates to reduce to pseudo-1D absorbing bound-ary problems. A transition from a regime where P e = 0to one where P e > x , ie.the case where the intermediate is so large as to behavelike a continuous phase, rather than a discrete collectionof monomers. In this limit, we obtain the classical resultsfor first-passage of Brownian motion with drift [81], anda jump in ∂ µ P e ( x c ) occurs about the point µ = 0 [31].Besides the discrete polymer example considered in thispaper, IP mechanisms with a number of other intermedi-ate geometries are likely to be mathematically tractable.For example, in the case of spherical intermediates witha continuous-like composition, the problem considered inSec. II. can be recast as a continuous advection-diffusion(Eq. (4)) along a 1D trajectory towards a spherical ab-sorbing boundary. In the simplified case of radial propa-gation, we have performed a quick calculation to estimatethe rate at which large spherical intermediates give riseto a growing stable polymers [31]. While this rate is pro-portional to the volume of the intermediates in the coop-erative regime, it is proportional to their surface area inthe competitive regime. This behavior closely resemblesthe way in which initiation sites for successful propaga-tion events are restricted to the ends of polymeric in-termediates in the competitive switching regime. Thus,a more general principle may exist, in which the com-petition between growth and propagation dynamics cre-ates two distinct regimes. In one, a polymerizing stablephase may originate from anywhere within the volumeof the intermediate; in the other, the sites where such4a phase may originate are restricted to the boundary ofthe intermediate. The fact that a critical point separat-ing these regimes appears to be a widespread feature oftwo-step biopolymer nucleation pathways underlines themore general principle that the prerequisites for criticalbehavior are relaxed in far-from-equilibrium systems.In addition to being compatible with a wide range of in-termediate morphologies, our IP mechanism also predicts significant morphological diversity. The incorporation ofa single alternative monomer conformation means thatpolymers of a given length may exhibit a large numberof different sequences of type-1 and type-2 monomers.Because the conformation of a monomer affects its me-chanical properties, we would expect these different se-quences to result in diverse morphologies. In this study,we have primarily focused on the mechanism of growthmode switching, and have ignored aspects of the confor-mational conversion dynamics which are not relevant tothis problem. However, in future studies it may be possi-ble to rigorously predict changing polymer morphologiesby obtaining analytical solutions to describe the progres-sive conversion of monomers throughout the polymer.As an example of how polymer composition might af-fect morphology, we observe that monomers in the type-1and type-2 states are likely to exhibit differing degrees ofconformational ordering. It is most often assumed thatthe intermediate is less ordered than the stable species, sothat conformational conversion can be regarded as an or-dering transition. In this case, the conformational degen-eracy of the type-1 state would be expected to result ingreater flexibility for type-1 domains. As the monomerswithin a polymer converted to the type-2 state, we wouldexpect to see a progressive increase in the persistencelength l p . However, in cases where the persistence lengthof the type-1 and type-2 monomers are markedly dif-ferent, presence of a small number of type-1 monomerswould be expected to strongly affect the overall morphol-ogy of the polymer. Thus, polymers with a small but sig-nificant type-1 content might exhibit radically differentmorphologies from mature polymers, as is often observedexperimentally [66, 82, 83].It would also be possible to introduce further con-formational diversity to our model, by introducing ad-ditional metastable conformations for the polymerizedmonomers. While analysis of such a model would not betrivial, additional conformational freedom and the pres-ence of other slow dynamics might be expected to yieldeven richer behaviors than those presented here. How-ever, it is useful to note that structural data from solid-state nuclear magnetic resonance (NMR) studies sug-gest that monomers within amyloid fibrils occupy a lim-ited range of stable or metastable conformations [84, 85],while high-resolution microscopy data have revealed a di-verse range of morphologies [82, 86]. Therefore, even inamyloid systems where the free monomer is often highlyconformationally degenerate, the range of accessible con-formations seems to remain limited within the polymer.Thus, a simple model in which polymerized monomers occupy one of two conformations may be sufficient to ac-count for much of the observed morphological diversity. B. Explicit conversion dynamics are essential to ageneral two-step biopolymer nucleation model
The inclusion of explicit conformational conversion dy-namics also allows our IP model to exhibit a generalitynot seen in other two-step nucleation models. Propa-gation of the stable phase is not assumed to be instan-taneous, and can occur on timescales similar to eitherpolymerization of the intermediate or non-autocatalyticconversion (Fig. 2). As a result, the dynamical behaviorof our model depends on two key parameters, µ/D and ζ .The first of these, µ/D , describes the competing effects ofpropagation and polymerization on the time-dependenceof growth mode switching. When µ/D <
0, the stablephase propagates through the intermediate faster thanthe intermediate phase can polymerize, causing acceler-ated switching. When µ/D >
0, polymerization is rapidand reduces the likelihood of successful propagation. Thesecond parameter, ζ , compares the timescale for propa-gation due to autocatalytic conversion with the timescalefor non-autocatalytic conversion. When ζ = 0, nearest-neighbor interactions between monomers do not promoteconformational conversion, so that the rate of growthmode switching is independent of the length of the poly-mer. When ζ = 1, nearest-neighbor interactions stronglypromote conversion, so that the timescale for propagationis much faster than that of non-autocatalytic conversion,and the rate of growth mode switching is accelerated.Existing models which neglect the conversion dynam-ics are not able to consider the variation accounted forby µ/D and ζ , and thus make various assumptions re-garding the rate at which the new growth mode emerges.While certain models [20, 22] have correctly supposeda connection between intermediate size and the rate ofemergence of a growing stable phase (eg. k c (cid:104) τ s (cid:105) ≈ /j ),the lack of competing dynamics means they are restrictedto a limited range of points close to the top left ( ζ → µ/D <
0) of our phase diagram (Fig. 3(a)). In contrast,other models have assumed that there is no relationshipbetween intermediate size and the switching rate [19, 21](ie. k c (cid:104) τ s (cid:105) ∝ ζ = 0); and when growth is suf-5ficiently rapid that growth mode switching is inhibited( µ/D → ∼ k B T ) would result in a significant ζ value,and in the vast majority of our phase space some level ofreduction in k c (cid:104) τ s (cid:105) is observed (Fig. 3(a)).We have seen that existing models of two-step biopoly-mer nucleation describe a limited range of possible sce-narios, which are restricted to points close to three lines( ζ → µ/D < ζ = 0; and µ/D →
2) on thelimits of our phase diagram (Fig. 3(a)). As a result,these models exhibit markedly different behaviors whichappear difficult to reconcile. By explicitly consideringthe dynamics by which conversion of monomers withina polymeric intermediate leads to emergence of a newgrowth mode, we are able to explore a much broaderregion of parameter space, and show that these appar-ently contrasting scenarios are in fact different limits ofthe same underlying reaction scheme. In addition, wefind that these previously unexplored regions contain richphase behavior which is not anticipated based on the be-havior observed in the limits. Therefore, our model hasthe capacity to explain a wide variety of ‘poorly behaved’experimental or physiological systems which are not de-scribed by existing models, such as those with mixed in-termediates or low concentration-dependence (Sec. IV.B. & C.). Moreover, by classifying existing two-stepnucleation models in the context of our phase diagram(Fig. 3(a)), we have shown that these models occupy dis-parate regions of phase space which can only be unitedby considering the competitive dynamics accounted forby µ/D and ζ . Therefore, while our model still does notdescribe some of the more complex scenarios that havebeen suggested [23], we have shown that explicit conver-sion dynamics are essential features of a general modelof two-step biopolymer nucleation. We thus believe thatour study represents a crucial step towards such a model. C. Biological importance of the nonequilibriumcritical point
In addition to generalizing existing models, the inclu-sion of an explicit conformational conversion process al-lows our model to predict a nonequilibrium critical point.While critical points have previously been described inbiopolymer systems, the possibility of nonequilibriumcriticality has largely been overlooked. Instead, the mostcommonly suggested critical point is the critical micel-lar/oligomer concentration (CMC) [52, 77, 87], which isthe concentration at which the free monomer and inter-mediate are in local equilibrium, such that m = k − /k +1 .CMC models typically assume that formation of interme-diates causes rapid depletion of the free monomer, so that k − /k +1 represents an upper limit on the free monomerconcentration; as a result, saturation behavior occurs.In this paper, we consider a different scenario. We fo-cus primarily on the case where the free monomer de-pletes slowly, as occurs in physiological circumstanceswhere there is a monomer source, or in many experi-mental situations [88–91]. Thus, while our basic model(Fig. 1) has the required processes to produce a CMC-type effect, we are primarily interested in the case wherethe free monomer is supersaturated with respect to theintermediate ( m > k − /k − ), and thus out of local equi-librium. Under these circumstances, a second criticalpoint emerges at m = ( k − + k p ) /k +1 (Fig. 3 & 5), ata higher concentration than the CMC. On our nonequi-librium phase diagram (Fig. 3a), the CMC would mani-fest as a vertical line (ie. constant µ/D ) situated to theleft of the nonequilibrium critical point. The differencein concentration between these two points will be givenby k p /k +1 . Although we set k − = 0 in our simulationsto focus primarily on the nonequilibrium critical point,our mathematical analysis is valid for more general k − subject to the condition that destabilization of interme-diates is negligible. Thus, the CMC and the nonequi-librium critical point described in this paper are com-patible so long as the concentration difference betweenthe two is sufficiently large. This situation is highly rel-evant to nonequilibrium scenarios where the dissociationof monomers from the intermediates is relatively slow, asoccurs in Alzheimer’s and Parkinson’s diseases [89–91].In cases where k p (cid:28) k − , the CMC and nonequilibriumcritical point will become close to one another on thephase diagram, and the presence of two non-negligibleabsorbing boundaries in Eq. (1) will further complicatethe dynamics of self-assembly and conversion. A rigor-ous analysis of this effect is both mathematically involvedand outside the scope of this paper, and so has been leftas a subject for future work.It is interesting to note that in cases where physio-logical and experimental systems exhibit slow variationin free monomer levels over time, the phase behavior ofbiopolymer self-assembly will vary accordingly. Experi-mental systems are likely to exhibit a monotonic deple-tion of the free monomer. In this case, µ/D will decreaseover time, causing the system to pass the nonequilib-rium critical point en route to the CMC. The forma-tion of stable polymers due to growth mode switchingwill then cause further monomer depletion, resulting ineither destabilization or switching of the remaining in-termediates. Thus, accelerated growth mode switchingwill occur when the free monomer concentration passesthe nonequilibrium critical point at m = ( k − + k p ) /k +1 ,with destabilization of the remaining intermediates oc-curring shortly after. This sort of behavior has been ob-served for the amyloid- β (A β ) peptide. Under conditionsof high supersaturation, Chimon et al. [7] described theformation of large, spheroidal intermediates by A β , sim-ilar to those discussed in Sec. V. A. Over time, theseintermediates progressively acquired secondary structure6content similar to the stable polymeric phase, but grow-ing, stable polymers only appeared following depletion ofthe free monomer. This behavior is strongly indicativeof a transition from the competitive switching regime tothe cooperative regime, as described by our model.The biological implications of the nonequilibrium crit-ical point extend beyond the biopolymer self-assemblyprocess itself. The rapid change in steady-state behaviorwhen µ/D ≈ D. Relationship to autocatalytic secondarynucleation
The model developed here can be applied to any two-step biopolymer nucleation process. A fundamental as-pect of our model is its ability to predict mixed inter-mediates consisting of self-assembled monomers in twodistinct conformations, resulting in diverse intermediatemorphologies. Conversion of a self-assembled monomerfrom the metastable type-1 state to the stable type-2state occurs either non-autocatalytically, with rate k c , orautocatalytically with rate k c + k p . Thus, the stable phaseis able to initiate locally within the intermediate, andthen propagate autocatalytically to the active end. This type of intramolecular structural propagation has beenextensively documented in the biochemical literature inmodels of multi-subunit proteins [7, 15, 24, 26, 27, 94],and is responsible for the nonequilibrium criticality thatour model exhibits.A different type of autocatalysis which is often dis-cussed in models of biopolymer nucleation is secondarynucleation, the process by which the surface of an ex-isting stable biopolymer catalyzes the formation of ad-ditional biopolymers by heterogeneous nucleation [95–98]. Thus, while the proposed IP mechanism can beconsidered a form of intra-polymer autocatalysis, sec-ondary nucleation is an example of inter-polymer auto-catalysis. Secondary nucleation has been used to explainfeatures such as exponential ( M ( t ) ∝ exp( κt )) ratherthan quadratic ( M ( t ) ∝ t ) accumulation of the stablephase, and the formation of significant concentrations oftoxic intermediates during the ‘growth phase’ (ie. when t ≈ τ inf ) [97]. Given the prominence of secondary nucle-ation models in the protein polymerization literature, itis pertinent to question how they are mechanistically re-lated to our IP model. At the microscopic level, the pro-cess of secondary nucleation results from the alignmentof free monomers to the surface of the mature biopoly-mer [99]. In the context of our model, this alignmentcould promote initiation of type-1 polymers which thendetach from the mature polymer and independently un-dergo growth mode switching. This would manifest as anincrease in the initiation rate v n . There is also the possi-bility that interactions between conformationally orderedmonomers in the stable phase and disordered monomersin the nascent intermediate promote ordering of the lat-ter. This would create a more tightly coupled process inwhich conversion and propagation are accelerated, result-ing in an increase in k c and k p . Thus, it is possible thatthe IP mechanism and secondary nucleation share a com-mon mechanistic basis, and are simply different examplesof a common tendency for autocatalytic conformationalchange induced by interactions between monomers. Allthese effects are easy to incorporate into the frameworkof our IP model, and we anticipate that future studieswhich combine the IP mechanism with inter-polymer au-tocatalysis may be able to gain deeper insights into themicroscopic mechanism of secondary nucleation. E. Competing dynamics can explain lowconcentration-dependences
One of the key experimental validators for biopolymerassembly models is their ability to accurately describethe concentration-dependence of the process and relatethis to theory [5, 23, 32, 33, 40, 55, 56, 97, 100, 101].The Oosawa model [32, 33] was first to ascribe a physicalmeaning to the value of γ (Eq. (25)), the scaling expo-nent for the concentration-dependence of accumulationof the stable polymer phase (Sec. IV. C.). The Oosawamodel predicts that γ = n c /
2, where n c is the effective7order of nucleation and often represents the minimumsize for a polymer to be stable or metastable (ie. j ) [42].While n c is often interpreted as reflecting the number ofmonomers within the critical nucleus ( n ∗ = n c −
1) [102],which is the least stable species during the self-assemblypathway, this explanation has been challenged by thewidespread observation of γ < γ would imply n ∗ <
1, which is physically implausible[56]. In addition, many biopolymer systems with higher γ still exhibit much lower concentration-dependences thanwould be expected, based on biophysical predictions oftheir n ∗ [20, 23, 103]. Some authors have proposedthat saturation effects, such as the presence of a sat-urable catalytic surface which induces heterogeneous nu-cleation, may explain the low concentration-dependences[53, 54, 104]. Indeed, γ < γ appears to saturate with increasingmonomer concentration. However, under some exper-imental conditions the concentration-dependence dropsfurther to γ < . γ without invoking a saturation ef-fect. In the NCP model, low γ values are explained in thecontext of a ‘cascade nucleation system’, where sequentialconversion of N − n c / ( N + 1) ≤ γ ≤ n c /
2. Inthis model, although the number of intermediate typescan be infinite, their conversion occurs as a single-stepprocess without considering the underlying dynamics,and the intermediates do not exhibit a growth process.As a result, the conversion mechanism of an NCP in-termediate corresponds to the scenario predicted at thepoint ( µ/D, ζ ) = ( − ,
0) on our nonequilibrium phasediagram (Fig. 3(a)).The IP mechanism proposed here generalizes the con-version mechanism considered in the NCP model to ex-plicitly consider the dynamics by which intermediatesgrow and individual monomers convert. As a result, ourmodel makes predictions regarding the pre-steady andsteady-state scaling behaviors of intermediates and ma-ture polymers which can easily be generalized to the caseof multiple intermediates in the same manner as the NCPmodel. However, the inclusion of additional dynamics al-lows our model to explore a broader region of parameterspace, leading to our observation of sharp changes in γ near the critical point (Fig. 5(b) & (c)). This results in lo-cal depression of the concentration-dependence, which ismore pronounced and occurs more abruptly at higher ζ .High ζ values correspond to systems in which autocatal-ysis plays a dominant role in conformational conversion,and can be achieved on comparatively modest free energyscales (Sec. IV. A.). This reduction in concentration-dependence allows our far-from-equilibrium model to vi-olate constraints on the concentration-dependence of ex- isting models over a broad range of µ/D values, and thiseffect is likely to be further enhanced when saturationeffects are also active. Under many experimental condi-tions, a value of γ < F. Interpretation of the nucleus size
Fig. 6(c) shows how γ varies with µ/D for different ζ values. While both µ/D and ζ are likely to dependon the physicochemical environment under which self-assembly occurs, and the underlying properties of thebiopolymer system, µ/D exhibits an additional depen-dence on the free monomer concentration m , which caneasily be varied in experimental contexts. It is useful tonote that different ζ values result in highly characteristicchanges in concentration-dependence (Fig 6(c)), so thatthe ζ value intrinsic to a particular biopolymer systemcan be determined by experimentallly varying µ/D andmeasuring the changing concentration-dependence of as-sembly of the stable polymer. Thus, while it is possibleto determine n c from the values of γ observed in the lim-its where µ/D → ±
2, our model predictions provide aneasier alternative. Instead, we anticipate that it will bepossible to determine n c and ζ from the characteristiccurve shape with which γ varies at intermediate valuesof µ/D . This may be preferable as it eliminates the needto explore extreme values of µ/D , which is challengingdue to lack of sensitivity of experimental methods at lowconcentrations, and the likely presence of additional com-peting processes at high concentrations. This approachcan also be generalized to cases similar to the NCP model[21] in which multiple intermediates are present. In thisscenario, we would typically expect the Oosawa behav-ior ( γ = n c /
2) to be recovered in the µ/D → − γ = n c / ( N + 1)) to be re-covered in the µ/D → n c as being equal to the minimum sizeof a metastable intermediate j , future models which con-sider additional complexities of the conversion dynamicsmay arrive at different interpretations for n c . For ex-ample, a previous coarse-grained simulation study [22]which investigated a scenario similar to our cooperativeswitching regime identified n c with the minimum size8at which conversion can initiate within an intermediate,rather than the absolute minimum size of an intermedi-ate. While we take these sizes to be the same in ourmodel, future work which treats these quantities as dif-ferent may provide analytical justification for this result. VI. CONCLUSIONS
We have developed a minimal model of two-stepbiopolymer nucleation which accounts for conformationalconversion of metastable intermediates by an initiation-propagation (IP) mechanism. In our framework, the sta-ble polymer phase initiates locally within the intermedi-ate, and propagates by autocatalytic nearest-neighbor in-teractions between the assembled monomers. This mech-anism is well-supported by experimental studies showingprogressive ordering of intermediates on slow timescales[7, 15, 24–27], and computational studies suggesting anautocatalytic mechanism for structural conversion of themonomers [20, 24, 28, 29]. The inclusion of conforma-tional ordering on timescales comparable to condensa-tion of the intermediates causes our model to exhibitrich dynamics and a nonequilibrium critical point, whichare not predicted by models with an implicit orderingstep [19–23]. Furthermore, by explicitly considering thedynamics of conformational conversion, our frameworkunifies existing models, and shows that their apparentlyunconnected behaviors occur as different limits of a sin-gle underlying mechanism. The IP mechanism has thepotential to explain a wide variety of experimental be-haviors which are not predicted by existing models, suchas the formation of large, conformationally mixed inter-mediates [6, 7, 14, 15], and nucleation of the stable poly-mer phase at a rate independent of the free monomerconcentration [55]. In our model, these phenomena arisenaturally as a result of the competing dynamics of self-assembly and conversion, resolving the apparent para- dox that concentration-independent nucleation rates sug-gested critical nuclei with a non-physical size [56].Biopolymer nucleation intermediates are responsi-ble for toxicity in a huge range of human diseases[16, 18, 71–74], and the nonequilibrium critical pointpredicted by our model may explain the extremesensitivity of many polymerization disorders to smallchanges in self-assembly conditions [105–107]. As aresult, our model provides a deeper understanding ofthe origin of these toxic intermediates, and suggestspromising therapeutic strategies to treat diseases causedby their formation. Moreover, future theoretical studieswhich build on the results presented here are likely toyield analytical solutions for the size distribution andcomposition of nucleation intermediates, which willprovide a powerful tool in therapeutic development.Biopolymers are also increasingly exploited by materialsscientists [108–110], and the development of a generalmodel of their nucleation will allow us to regulate theirformation and prevent unwarranted side reactions. Inthis paper, we have shown that explicit conformationalconversion dynamics, such as those addressed by ourinitiation-propagation mechanism, are essential to sucha theory. We thus believe that our work representsa crucial step towards a general theory of two-stepbiopolymer nucleation, with important implications forboth experiments and human disease.
ACKNOWLEDGMENTS
This work was funded by the BBSRC (grant numberBB/P002927/1). AIPT gratefully acknowledges financialsupport from the BBSRC and the University of Sheffield.BC would like to thank the hospitality of MPIPKS, Dres-den where part of this work was completed. [1] J. A. Hardy and G. A. Higgins, Alzheimer’s Disease:The Amyloid Cascade Hypothesis, Science , 184(1992).[2] M. G. Spillantini, A. R. Crowther, R. Jakes, N. J.Cairns, P. L. Lantos, and M. Goedert, Filamentous α -synuclein inclusions link multiple system atrophy withParkinson’s disease and dementia with Lewy bodies,Neurosci. Lett. , 205 (1998).[3] L. Pauling, H. A. Itano, S. J. Singer, and I. C. Wells,Sickle cell anemia, a molecular disease, Science , 543(1949).[4] J. D. Harper, S. S. Wong, C. M. Lieber, and P. T.Lansbury, Observation of metastable Abeta amyloidprotofibrils by atomic force microscopy, Chem. Biol. ,119 (1997).[5] T. R. Serio, A. G. Cashikar, A. S. Kowal, G. J. Saw-icki, J. J. Moslehi, L. Serpell, M. F. Arnsdorf, andS. L. Lindquist, Nucleated Conformational Conversion and the Replication of Conformational Information bya Prion Determinant, Science. , 1317 (2000).[6] O. Galkin, W. Pan, L. Filobelo, R. E. Hirsch, R. L.Nagel, and P. G. Vekilov, Two-step mechanism of homo-geneous nucleation of sickle cell hemoglobin polymers,Biophys. J. , 902 (2007).[7] S. Chimon, M. A. Shaibat, C. R. Jones, D. C. Calero,B. Aizezi, and Y. Ishii, Evidence of fibril-like β -sheetstructures in a neurotoxic amyloid intermediate ofAlzheimer’s β -amyloid., Nat. Struct. Mol. Biol. , 1157(2007).[8] S. Auer, C. M. Dobson, and M. Vendruscolo, Character-ization of the nucleation barriers for protein aggregationand amyloid formation, HFSP J. , 137 (2007).[9] J. Lee, E. K. Culyba, E. T. Powers, and J. W. Kelly,Amyloid- β forms fibrils by nucleated conformationalconversion of oligomers, Nat. Chem. Biol. , 602 (2011).[10] A. Levin, T. O. Mason, L. Adler-Abramovich, A. K. Buell, G. Meisl, C. Galvagnion, Y. Bram, S. A. Strat-ford, C. M. Dobson, T. P. J. Knowles, and E. Gazit,Ostwald’s rule of stages governs structural transitionsand morphology of dipeptide supramolecular polymers.,Nat. Commun. , 5219 (2014).[11] R. Pellarin, E. Guarnera, and A. Caflisch, Pathways andIntermediates of Amyloid Fibril Formation, J. Mol. Biol. , 917 (2007).[12] B. Urbanc, M. Betnel, L. Cruz, G. Bitan, and D. B.Teplow, Elucidation of amyloid β -protein oligomeriza-tion mechanisms: Discrete molecular dynamics study,J. Am. Chem. Soc. , 4266 (2010).[13] E. Chatani, R. Inoue, H. Imamura, M. Sugiyama,M. Kato, M. Yamamoto, K. Nishida, and T. Kanaya,Early aggregation preceding the nucleation of insulinamyloid fibrils as monitored by small angle X-ray scat-tering, Sci. Rep. , 15485 (2015).[14] M. Zhu, P. O. Souillac, C. Ionescu-Zanetti, S. A. Carter,and A. L. Fink, Surface-catalyzed amyloid fibril forma-tion, J. Biol. Chem. , 50914 (2002).[15] J. Luo, S. K. T. S. W¨arml¨ander, A. Gr¨aslund, andJ. P. Abrahams, Alzheimer peptides aggregate into tran-sient nanoglobules that nucleate fibrils, Biochemistry , 6302 (2014).[16] Y. Shin and C. P. Brangwynne, Liquid phase con-densation in cell physiology and disease, Science. ,eaaf4382 (2017).[17] F. Chiti and C. Dobson, Protein Misfolding, FunctionalAmyloid, and Human Disease, Annu. Rev. Biochem. ,333 (2006).[18] I. Benilova, E. Karran, and B. De Strooper, The toxicA β oligomer and Alzheimer’s disease: an emperor inneed of clothes., Nat. Neurosci. , 349 (2012).[19] A. Vitalis and R. V. Pappu, Assessing the contributionof heterogeneous distributions of oligomers to aggrega-tion mechanisms of polyglutamine peptides, Biophys.Chem. , 14 (2011).[20] S. Auer, P. Ricchiuto, and D. Kashchiev, Two-step nu-cleation of amyloid fibrils: Omnipresent or not?, J. Mol.Biol. , 723 (2012).[21] G. A. Garcia, S. I. A. Cohen, C. M. Dobson, andT. P. J. Knowles, Nucleation-conversion-polymerizationreactions of biological macromolecules with prenucle-ation clusters, Phys. Rev. E , 032712 (2014).[22] A. ˇSari´c, T. C. Michaels, A. Zaccone, T. P. Knowles, andD. Frenkel, Kinetics of spontaneous filament nucleationvia oligomers: Insights from theory and simulation, J.Chem. Phys. , 211926 (2016).[23] C. T. Lee and E. M. Terentjev, Mechanisms and ratesof nucleation of amyloid fibrils, J. Chem. Phys. ,105103 (2017).[24] S. Auer, C. M. Dobson, M. Vendruscolo, and A. Mar-itan, Self-templated nucleation in peptide and proteinaggregation, Phys. Rev. Lett. , 258101 (2008).[25] D. W. Li, S. Mohanty, A. Irb¨ack, and S. Huo, Formationand growth of oligomers: A monte carlo study of anamyloid tau fragment, PLoS Comput. Biol. , e1000238(2008).[26] S. Kumar and J. B. Udgaonkar, Conformational Con-version May Precede or Follow Aggregate Elongation onAlternative Pathways of Amyloid Protofibril Formation,J. Mol. Biol. , 1266 (2009).[27] S. W. Chen, S. Drakulic, E. Deas, M. Ouberai, F. A.Aprile, R. Arranz, S. Ness, C. Roodveldt, T. Guilliams, E. J. De-Genst, D. Klenerman, N. W. Wood, T. P. J.Knowles, C. Alfonso, G. Rivas, A. Y. Abramov, J. M.Valpuesta, C. M. Dobson, and N. Cremades, Structuralcharacterization of toxic oligomers that are kineticallytrapped during α -synuclein fibril formation., Proc. Natl.Acad. Sci. U. S. A. , E1994 (2015).[28] F. Ding, J. J. LaRocque, and N. V. Dokholyan, Di-rect observation of protein folding, aggregation, anda prion-like conformational conversion, J. Biol. Chem. , 40235 (2005).[29] J. Lipfert, J. Franklin, F. Wu, and S. Doniach, Proteinmisfolding and amyloid formation for the peptide GN-NQQNY from yeast prion protein sup35: Simulation byreaction path annealing, J. Mol. Biol. , 648 (2005).[30] R. J. Glauber, Time-Dependent Statistics of the IsingModel, J. Math. Phys. , 294 (1963).[31] See Supplemental Material at [*URL to be inserted] fora more detailed description of our mathematical analysisand simulations.[32] F. Oosawa and M. Kasai, A theory of linear and heli-cal aggregations of macromolecules, J. Mol. Biol. , 10(1962).[33] F. Oosawa and S. Asakura, Thermodynamics of thePolymerization of Protein (Academic Press, New York,1975).[34] R. Peierls and M. Born, On Ising’s model of ferromag-netism, Math. Proc. Camb. Philos. Soc. , 477 (1936).[35] L. Onsager, Crystal statistics. I. A two-dimensionalmodel with an order-disorder transition, Phys. Rev. ,117 (1944).[36] C. Goldsbury, J. Kistler, U. Aebi, T. Arvinte, and G. J.Cooper, Watching amyloid fibrils grow by time-lapseatomic force microscopy, J. Mol. Biol. , 33 (1999).[37] T. Ban, M. Hoshino, S. Takahashi, D. Hamada,K. Hasegawa, H. Naiki, and Y. Goto, Direct observa-tion of A β amyloid fibril growth and inhibition, J. Mol.Biol. , 757 (2004).[38] M. Sleutel, I. Van Den Broeck, N. Van Gerven, C. Feuil-lie, W. Jonckheere, C. Valotteau, Y. F. Dufrˆene, andH. Remaut, Nucleation and growth of a bacterial func-tional amyloid at single-fiber resolution, Nat. Chem.Biol. , 902 (2017).[39] Y. Xu, M. S. S. Safari, W. Ma, N. P. Schafer, P. G.Wolynes, and P. G. Vekilov, Steady, Symmetric, and Re-versible Growth and Dissolution of Individual Amyloid- β Fibrils, ACS Chem. Neurosci. , 2967 (2019).[40] J. Straub and D. Thirumalai, Toward a Molecular The-ory of Early and Late Events in Monomer to Amy-loid Fibril Formation, Annu. Rev. Phys. Chem. , 437(2011).[41] M. Iljina, G. A. Garcia, M. H. Horrocks, L. Tosatto,M. L. Choi, K. A. Ganzinger, A. Y. Abramov,S. Gandhi, N. W. Wood, N. Cremades, C. M. Dobson,T. P. J. Knowles, and D. Klenerman, Kinetic model ofthe aggregation of alpha-synuclein provides insights intoprion-like spreading, Proc. Natl. Acad. Sci. U. S. A. ,E1206 (2016).[42] F. Ferrone, Analysis of protein aggregation kinetics.,Methods Enzymol. , 256 (1999).[43] S. I. Cohen, M. Vendruscolo, M. E. Welland, C. M.Dobson, E. M. Terentjev, and T. P. Knowles, Nucleatedpolymerization with secondary pathways. I. Time evo-lution of the principal moments, J. Chem. Phys. ,065105 (2011). [44] T. P. J. Knowles, C. A. Waudby, G. L. Devlin, S. I. A.Cohen, A. Aguzzi, M. Vendruscolo, E. M. Terentjev,M. E. Welland, and C. M. Dobson, An analytical so-lution to the kinetics of breakable filament assembly,Science , 1533 (2009).[45] T. C. T. Michaels, G. A. Garcia, and T. P. J. Knowles,Asymptotic solutions of the Oosawa model for thelength distribution of biofilaments, Proc. Natl. Acad.Sci. U. S. A. , 194906 (2014).[46] A. B. Bortz, M. H. Kalos, and J. L. Lebowitz, A newalgorithm for Monte Carlo simulation of Ising spin sys-tems, J. Comput. Phys. , 10 (1975).[47] D. Gillespie, A general method for numerically simu-lating the stochastic time evolution of coupled chemicalreactions, J. Comput. Phys. , 403 (1976).[48] A. F. Voter, Introduction to the Kinetic Monte CarloMethod, in Radiation Effects in Solids , edited by K. E.Sickafus and E. A. Kotomin (Springer, NATO Publish-ing Unit, Dordrecht, The Netherlands, 2005) Chap. 1,pp. 1–23.[49] R. J. Morris, K. Eden, R. Yarwood, L. Jourdain, R. J.Allen, and C. E. MacPhee, Mechanistic and environ-mental control of the prevalence and lifetime of amyloidoligomers, Nat. Commun. , 1891 (2013).[50] K. Eden, R. Morris, J. Gillam, C. E. MacPhee, and R. J.Allen, Competition between primary nucleation and au-tocatalysis in amyloid fibril self-assembly, Biophys. J. , 632 (2015).[51] H.-J. Bartsch, Handbook of Mathematical Formulas (Academic Press, New York, 1974).[52] A. Lomakin, D. S. Chung, G. B. Benedek, D. A.Kirschner, and D. B. Teplow, On the nucleation andgrowth of amyloid beta-protein fibrils: detection of nu-clei and quantitation of rate constants., Proc. Natl.Acad. Sci. U. S. A. , 1125 (1996).[53] G. Meisl, X. Yang, E. Hellstrand, B. Frohm, J. B.Kirkegaard, S. I. A. Cohen, C. M. Dobson, S. Linse,and T. P. J. Knowles, Differences in nucleation behav-ior underlie the contrasting aggregation kinetics of theA β
40 and A β
42 peptides, Proc. Natl. Acad. Sci. U. S.A. , 9384 (2014).[54] J. Habchi, S. Chia, C. Galvagnion, T. C. Michaels,M. M. Bellaiche, F. S. Ruggeri, M. Sanguanini, I. Idini,J. R. Kumita, E. Sparr, S. Linse, C. M. Dobson, T. P.Knowles, and M. Vendruscolo, Cholesterol catalysesA β
42 aggregation through a heterogeneous nucleationpathway in the presence of lipid membranes, Nat. Chem. , 673 (2018).[55] J. E. Gillam and C. E. MacPhee, Modelling amyloidfibril formation kinetics: mechanisms of nucleation andgrowth, J. Phys.: Condens. Matter , 373101 (2013).[56] J. P. Bernacki and R. M. Murphy, Model Discriminationand Mechanistic interpretation of Kinetic Data in Pro-tein Aggregation Studies, Biophy. J. , 2871 (2009).[57] J. S. Langer, Statistical Theory of the Decay ofMetastable States, Ann. Phys. (Leipzig) , 258 (1969).[58] P. H¨anggi, P. Talkner, and M. Borkovec, Reaction-ratetheory: fifty years after Kramers, Rev. Mod. Phys. ,251 (1990).[59] H. A. Kramers, Brownian Motion in a Field of Force andthe Diffusion Model of Chemical Reactions, Physica ,284 (1940).[60] N. Ben-Tal, D. Sitkoff, I. A. Topol, A. S. Yang, S. K.Burt, and B. Honig, Free energy of amide hydrogen bond formation in vacuum, in water, and in liquidalkane solution, J. Phys. Chem. B , 450 (1997).[61] S.-Y. Sheu, D.-Y. Yang, and E. W. Schlag, Energeticsof hydrogen bonds in peptides, Proc. Natl. Acad. Sci.U. S. A. , 12684 (2003).[62] A. N. Naganathan and V. Mu˜noz, Scaling of FoldingTimes with Protein Size, J. Am. Chem. Soc. , 480(2005).[63] S. I. Cohen, R. Cukalevski, T. C. Michaels, A. ˇSari´c,M. T¨ornquist, M. Vendruscolo, C. M. Dobson, A. K.Buell, T. P. Knowles, and S. Linse, Distinct thermody-namic signatures of oligomer generation in the aggre-gation of the amyloid- β peptide, Nat. Chem. , 523(2018).[64] V. N. Uversky, Mysterious Oligomerization of the amy-loidogenic proteins, FEBS Journal , 2940 (2010).[65] J. P. Colletier, A. Laganowsky, M. Landau, M. Zhao,A. B. Soriaga, L. Goldschmidt, D. Flot, D. Cascio, M. R.Sawaya, and D. Eisenberg, Molecular Basis of amyloid- β polymorphism., Proc. Natl. Acad. Sci. U. S. A. ,16938 (2011).[66] R. Kodali and R. Wetzel, Polymorphism in the interme-diates and products of amyloid assembly, Curr. Opin.Struct. Biol. , 48 (2001).[67] V. N. Uversky, Under-folded proteins: conformationalensembles and their roles in protein folding, function,and pathogenesis, Biopolymers , 870 (2013).[68] M. M. Babu, R. van der Lee, N. Sanchez de Groot, andJ. Gsponer, Intrinsically disordered proteins: regulationand disease, Curr. Opin. Struct. Biol. , 432 (2011).[69] T. Vavouri, J. I. Semple, R. Garcia-Verdugo, andB. Lehner, Intrinsic protein disorder and interactionpromiscuity are widely associated with dosage sensitiv-ity, Cell , 198 (2009).[70] Y. Miller, B. Ma, and R. Nussinov, Polymorphism inAlzheimer A β amyloid organization reflects conforma-tional selection in a rugged energy landscape, ChemRev. , 4820 (2010).[71] M. P. Lambert, A. K. Barlow, B. A. Chromy, C. Ed-wards, R. Freed, M. Liosatos, T. E. Morgan, I. Ro-zovsky, B. Trommer, K. L. Viola, P. Wals, C. Zhang,C. E. Finch, G. A. Krafft, and W. L. Klein, Diffusible,nonfibrillar ligands derived from A β , 6448 (1998).[72] S. Lesn´e, M. T. Koh, L. Kotilinek, R. Kayed, C. G.Glabe, A. Yang, M. Gallagher, and K. H. Ashe, A spe-cific amyloid- β protein assembly in the brain impairsmemory, Nature , 352 (2006).[73] K. Ono, M. M. Condron, and D. B. Teplow,Structure-neurotoxicity relationships of amyloid β -protein oligomers, Proc. Natl. Acad. Sci. U. S. A. ,14745 (2009).[74] W. F. Xue, A. L. Hellewell, E. W. Hewitt, and S. E.Radford, Fibril fragmentation in amyloid assembly andcytotoxicity: when size matters, Prion , 20 (2010).[75] M. Ahmed, J. Davis, D. Aucoin, T. Sato, S. Ahuja,S. Aimoto, J. I. Elliott, W. E. Van Nostrand, and S. O.Smith, Structural conversion of neurotoxic amyloid-beta(1-42) oligomers to fibrils., Nat. Struct. Mol. Biol. , 561 (2010).[76] G. Bitan, M. Kirkitadze, A. Lomakin, S. Vollers,G. Benedek, and D. Teplow, Amyloid β -protein (A β )assembly: A β
40 and A β
42 oligomerize through distinct pathways, Proc. Natl. Acad. Sci. U. S. A. , 330(2003).[77] R. Sabat´e and J. Estelrich, Evidence of the existence ofmicelles in the fibrillogenesis of beta-amyloid peptide.,J. Phys. Chem. B , 11027 (2005).[78] W. Yong, A. Lomakin, M. Kirkitadze, D. Teplow, S.-H. Chen, and G. Benedek, Structure determination ofmicelle-like intermediates in amyloid β -protein fibril as-sembly by using small angle neutron scattering, Proc.Natl. Acad. Sci. U. S. A. , 150 (2002).[79] M. M. Murr and D. E. Morse, Fractal intermediatesin the self-assembly of silicatein filaments, Proc. Natl.Acad. Sci. U. S. A. , 11657 (2005).[80] D. G. Kendall, The generalized ’birth and death’ pro-cess, Ann. Math. Stat. , 1 (1948).[81] D. R. Cox and H. D. Miller, The Theory of Stochas-tic Processes (Chapman & Hall, CRC, Boca Raton,Florida, 1965).[82] J. Meinhardt, C. Sachse, P. Hortschansky, N. Grigorieff,and M. F¨andrich, A β (1-40) fibril polymorphism impliesdiverse interaction patterns in amyloid fibrils, J. Mol.Biol. , 869 (2007).[83] J. S. Pedersen, C. B. Andersen, and D. E. Otzen, Amy-loid structure - one but not the same: the many levels offibrillar polymorphism, FEBS Journal , 4591 (2010).[84] R. Tycko, Solid state NMR studies of amyloid fibrilstructure, Annu. Rev. Phys. Chem. , 279 (2011).[85] M. G. Iadanza, M. P. Jackson, E. W. Hewitt, N. A.Ranson, and S. E. Radford, A new era for understandingamyloid structures and disease, Nature Rev. Mol. CellBiol. , 755 (2018).[86] J. L. Jim´enez, E. J. Nettleton, M. Bouchard, C. V.Robinson, C. M. Dobson, and H. R. Saibil, The protofil-ament structure of insulin amyloid fibrils, Proc. Natl.Acad. Sci. U. S. A. , 9196 (2002).[87] F. Hasecke, T. Miti, C. Perez, J. Barton, D. Sch¨olzel,L. Gremer, C. S. R. Gr¨uning, G. Matthews, G. Meisl,T. P. J. Knowles, D. Willbold, P. Neudecker, H. Heise,G. Ullah, W. Hoyer, and M. Muschol, Origin ofmetastable oligomers and their effects on amyloid fib-ril self-assembly, Chem. Sci. , 5937 (2018).[88] J. A. Varela, M. Rodriguez, S. De, P. Flagmeier,S. Gandhi, C. M. Dobson, D. Klenerman, and S. F. Lee,Optical structural analysis of individual α -synucleinoligomers, Angewandte Chemie Int. Ed. , 4886(2018).[89] N. Narayan, A. Orte, R. W. Clarke, B. Bolognesi,S. Hook, K. A. Ganzinger, A. Meehan, M. R. Wilson,C. M. Dobson, and D. Klenerman, The extracellularchaperone clusterin sequesters oligomeric forms of theamyloid- β (1-40) peptide, Nature Struct. Mol. Biol. ,79 (2012).[90] P. Narayan, K. M. Holmstr´om, D. ˆa. Kim, D. J. Whit-comb, M. R. Wilson, P. St. George-Hyslop, N. W.Wood, C. M. Dobson, K. Cho, A. Y. Abramov, andD. Klenerman, Rare individual amyloid- β oligomers acton astrocytes to initiate neuronal damage, Biochemistry , 2442 (2014).[91] N. Cremades, S. I. A. Cohen, E. Deas, A. Y. Abramov,A. Y. Chen, A. Orte, M. Sandal, R. W. Clarke,P. Dunne, F. A. Aprile, C. W. Bertoncini, N. W. Wood,T. P. J. Knowles, C. M. Dobson, and D. Klenerman,Direct observation of the interconversion of normal andtoxic forms of α -synuclein, Cell , 1048 (2012). [92] V. N. Uversky and A. K. Dunker, Controlled Chaos,Science , 1340 (2008).[93] D. Otzen and P. H. Nielsen, We find them here, we findthem there: functional bacterial amyloid, Cell. Mol. LifeSci. , 910 (2008).[94] E. Papaleo, G. Saladino, M. Lambrughi, K. Lindorff-Larsen, F. L. Gervasio, and R. Nussinov, The role ofprotein loops and linkers in conformational dynamicsand allostery, Chem. Rev. , 6391 (2016).[95] F. A. Ferrone, J. Hofrichter, H. R. Sunshine, and W. A.Eaton, Kinetics studies on photolysis-induced gelationof sickle-cell hemoglobin suggest a new mechanism., Bio-phys. J. , 361 (1980).[96] A. M. Ruschak and A. D. Miranker, Fiber-dependentamyloid formation as catalysis of an existing pathway,Proc. Natl. Acad. Sci. U. S. A. , 12341 (2007).[97] S. I. A. Cohen, S. Linse, L. M. Luheshi, E. Hellstrand,D. A. White, L. Rajah, D. E. Otzen, M. Vendruscolo,C. M. Dobson, and T. P. J. Knowles, Proliferation ofamyloid- β
42 aggregates occurs through a secondary nu-cleation mechanism, Proc. Natl. Acad. Sci. U. S. A. ,9758 (2013).[98] T. C. T. Michaels, A. ˇSari´c, J. Habchi, S. Chia, G. Meisl,M. Vendruscolo, C. M. Dobson, and T. P. J. Knowles,Chemical kinetics for bridging molecular mechanismsand macroscopic measurements of amyloid fibril forma-tion, Annu. Rev. Phys. Chem. , 273 (2018).[99] M. T¨ornquist, T. C. T. Michaels, K. Sanagavarapu,X. Yang, G. Meisl, S. I. A. Cohen, T. P. J. Knowles,and S. Linse, Secondary nucleation in amyloid forma-tion, Chem. Commun. , 8667 (2018).[100] T. C. T. Michaels, A. J. Dear, and T. P. J. Knowles,Scaling and dimensionality in the chemical kinetics ofprotein filament formation, Int. Rev. Phys. Chem. ,679 (2016).[101] F. Carbonell, Y. Iturria-Medina, and A. C. Evans,Mathematical modeling of protein misfolding mecha-nisms in neurological diseases: a historical overview,Front. Neurol. , 38 (2018).[102] F. Ferrone, Assembly of A β proceeds via monomeric nu-clei, J. Mol. Biol. , 287 (2015).[103] N. S. Bieler, T. P. J. Knowles, D. Frenkel, and R. V´acha,Connecting Macroscopic Observables and MicroscopicAssembly Events in Amyloid Formation Using CoarseGrained Simulations, PLoS Comput. Biol. , e1002692(2012).[104] A. ˇSari´c, A. K. Buell, G. Meisl, T. C. T. Michaels, C. M.Dobson, S. Linse, T. P. J. Knowles, and D. Frenkel,Physical determinants of the self-replication of proteinfibrils, Nature Phys. , 874 (2016).[105] R. G. Snell, J. C. MacMillan, J. P. Cheadle, I. Fen-ton, L. P. Lazarou, P. Davies, M. E. MacDonald, J. F.Gusella, P. S. Harper, and D. J. Shaw, Relationshipbetween trinucleotide repeat expansion and phenotypicvariation in Huntington’s disease, Nat. Genet. , 393(1993).[106] C. A. Lemere, J. K. Blusztajn, H. Yamaguchi, T. Wis-niewski, T. C. Saido, and D. J. Selkoe, Sequence of De-position of Heterogeneous Amyloid β -Peptides and APOE in Down Syndrome: Implications for Initial Eventsin Amyloid Plaque Formation, Neurobiol. Dis. , 16(1996).[107] J. D. F. Wadsworth, E. A. Asante, M. Desbruslais, J. M.Linehan, S. Joiner, I. Gowland, J. Welch, L. Stone, S. E. Lloyd, A. F. Hill, S. Brandner, and J. Collinge, Humanprion protein with valine 129 prevents expression of vari-ant CJD phenotype., Science , 1793 (2004).[108] J. Wang, K. Liu, R. Xinga, and X. H. Yan, Peptide self-assembly: thermodynamics and kinetics, Chem. Soc.Rev. , 5589 (2016). [109] E. Gazit, Self-assembled peptide nanostructures: thedesign of molecular building blocks and their techno-logical utilization, Chem. Soc. Rev. , 1263 (2007).[110] G. Wei, Z. Su, N. P. Reynolds, P. Arosio, I. W. Hamley,E. Gazit, and R. Mezzenga, Self-assembling peptide andprotein amyloids: from structure to tailored function innanotechnology, Chem. Soc. Rev.46