The deposition of pathological protein aggregates in the brain plays a central role in cognitive decline and structural damage associated with neurodegenerative diseases. In Alzheimer’s disease, the formation of amyloid-β plaques and neurofibrillary tangles of the tau protein is associated with the appearance of symptoms and pathology. Detailed models for the specific mechanisms of aggregate formation, such as nucleation and elongation, exist for aggregation in vitro where the total protein mass is conserved. However, in vivo, an additional class of mechanisms that clear pathological species is present and is believed to play an essential role in limiting the formation of aggregates and preventing or delaying the emergence of disease. A key unanswered question in the field of neuro-degeneration is how these clearance mechanisms can be modeled and how alterations in the processes of clearance or aggregation affect the stability of the system toward aggregation. Here, we generalize classical models of protein aggregation to take into account both production of monomers and the clearance of protein aggregates. We show that, depending on the specifics of the clearance process, a critical clearance value emerges above which accumulation of aggregates does not take place. Our results show that a sudden switch from a healthy to a disease state can be caused by small variations in the efficiency of the clearance process and provide a mathematical framework to explore the detailed effects of different mechanisms of clearance on the accumulation of aggregates.

Alzheimer’s disease (AD) and other related neurodegenerative diseases are associated with the assembly of specific proteins into pathological fibrillar aggregates. Alzheimer’s disease, in particular, is characterized by the aggregation of amyloid-β (Aβ) into plaques and of the tau protein into neurofibrillary tangles (NFTs). The role of Aβ in Alzheimer’s is thought to be so central to the disease that it is the basis of the so-called “amyloid-β hypothesis,”1–3 stating that the accumulation and deposition of oligomeric or fibrillar amyloid-β peptide are the main cause of the disease. This hypothesis has provided a guide for most of the AD research over the last 20 years. However, recent experimental evidence, and the failure of several drug trials, has led to renewed scrutiny of this foundational assumption, and next generation therapeutic intervention strategies that target more specifically the low molecular weight oligomers are showing more promise.4 

The production of Aβ is a natural process related to neuronal activity. Indeed, Aβ is a normal metabolic waste by-product5,6 that is typically removed from intracellular and extracellular compartments by several clearance mechanisms.7,8 In healthy subjects, waste proteins are broken down by enzymes, removed by cellular uptake, crossing the blood–brain barrier, or efflux to cerebrospinal fluid compartments where they eventually reach arachnoid granulations, or lymphatic vessels. While healthy clearance mechanisms, working in harmony, avert the buildup of toxic Aβ plaques and tau NFT, their impairment or dysfunction can lead to pathology.8 The specifics of in vivo clearance mechanisms remain a topic of clinical debate; however, the kinetics enabling proteins to amass into pathological aggregates can be carefully, and systematically, studied in vitro and under varied conditions. The production of Aβ, at a high level, is mediated by a membrane protein, the amyloid precursor protein (APP). APP is typically cleaved by α-secretase, and the resulting products do not aggregate. However, APP can also be cleaved by β-secretase and γ-secretase, a process that results in soluble monomeric APP fragments of different sizes. The most common variants are Aβ38, Aβ40, and Aβ42. While monomeric Aβ38 is not prone to further aggregation, Aβ40 and Aβ42, containing additional amino acids at the C terminus, are the main isoforms of interest in the study of AD pathology.

Protein aggregation pathways are, in general, complex and involve multiple steps.9 In fact, it has recently been shown10 that the aggregation properties of Aβ40, which is more abundant, differ from those of the more aggregate-prone Aβ42 even under the same conditions. A theoretical framework of chemical kinetics and aggregation theory11–13 has been combined with careful, systematic in vitro experiments performed under differing conditions, such as varied concentration or pH. This approach based on chemical kinetics has elucidated effective pathways and mechanisms for nucleation, aggregation, and fragmentation14 and produced a deep understanding of key properties, underlying the formation of aggregates under ideal conditions, with the potential for therapeutic intervention.15,16

Here, we extend this framework and develop a mathematical approach to describe the effects of clearance and monomer production on the kinetics of aggregation. We apply the framework to the study of Aβ. To accomplish this objective, we extend the current theory describing Aβ aggregation in vitro, which has been validated against experiment, to include monomer production and oligomer clearance terms. In particular, we study two different clearance mechanisms: one where all aggregated species are cleared at the same rate (size-independent clearance) and the other where different aggregated species may be cleared at different speeds (size-dependent clearance). In the former case, we show the full system reduces to three equations amenable to a systematic analysis. We identify a critical value of clearance above which the production of aggregates does not take place. Our results offer further evidence in support of two main hypotheses that clearance mechanisms play a crucial role in neurodegenerative disease initiation and progression and that therapies enhancing clearance above a prescribed, critical value may serve as a possible intervention strategy. In particular, we will exhibit the existence of critical clearance values, which are consistent with the observation of disease onset when natural clearance mechanisms within the brain have degraded through aging.

Our model for protein aggregation-dynamics includes the key molecular steps: heterogeneous primary nucleation, homogeneous primary nucleation, secondary nucleation, linear elongation, and clearance (cf. Fig. 1). These mechanisms lead to a general class of mathematical models that can describe a wide range of aggregating systems in vitro. In particular, by including heterogeneous primary nucleation terms, a source term for new nuclei, which can become independent of monomer concentration, is present;17 we include this here in addition to the usual monomer-dependent homogeneous primary nucleation. Thus, in such a model, the importance of interfaces in the initiation of nucleation is sufficiently accounted for. In the model, each aggregate of a given size is represented by a population. In general, each population, with aggregates of size i, will be represented by an indexed concentration. We use the special notation m(t) for the monomer population i = 1, while all other aggregate concentrations are denoted by pi(t) for i2,3,. The master equations are then

dmdt=γλ1mnckhρ(m)ncknmnc2k+mPn2k2σ(m)M+2koffP,
(1)
dpidt=λipi+δi,nc(knmnc+khρ(m))+2k+m(pi1pi)+2koff(pi+1pi)+δi,n2k2σ(m)M,
(2)

where we require inc > 1 and δi,j is Kronecker’s delta (δi,j = 1 if i = j and 0 otherwise) and

σ(m)=mn2KMKM+mn2,ρ(m)=mncKPKP+mnc,P=i=ncpi,M=i=ncipi.
(3)

Here, P and M are the first two moments of the population distribution. They represent the total number and total mass concentrations of aggregates, respectively. The functions σ(m) and ρ(m) account for the ability of both secondary nucleation and heterogeneous primary nucleation to saturate. As fibrils can generally grow from both ends, the conventional definition of the rates is with respect to the concentration of fibril ends, which is twice the number concentration of aggregates, giving rise to a factor of 2, in (1) and (2). In these equations, the parameters represent the following effects, sketched in Fig. 1: γ: (constant) monomer production such as by secretase mitigated cleavage of APP, driving mass influx; λi: clearance of aggregates of size i via, e.g., cerebrospinal-fluid mediated, blood–brain barrier mediated, or cellular degradation processes; kh: heterogeneous primary nucleation rate constant; kn: homogeneous primary nucleation of aggregates; nc: size of the smallest stable fibril (the nucleus) formed by either process of primary nucleation (nc > 1); k+: linear elongation transforming aggregate from size i to i + 1; k2: secondary nucleation of aggregates of size n2nc; KM: saturation constant of the secondary nucleation; KP: saturation constant of the heterogeneous primary nucleation; and koff: depolymerization by loss of one monomer from the fibril end. When fitting experimental data, it is often the case that only one of heterogeneous primary nucleation or homogeneous primary nucleation is required to explain the data. We include both here for completeness. Finally, a general note on the choice of model is that the aim here is to provide a tractable model with a focus on the kind of behavior that emerges when aggregate removal processes are active. We have therefore chosen a level of coarse-graining that is sufficient to describe the majority of currently available experimental data. In particular, we use one rate constant per process to describe elongation, de-polymerization, and nucleation. In reality, an aggregated system is likely to be heterogeneous and there may be slight variations of the rate constants, for example, with fibril length.18 Nevertheless, our model of aggregation in vitro captures the key physics as evidenced by the its successful application to modeling numerous experimental systems well to within the measurement accuracy.14,19 However, one should keep in mind that the rate constants extracted from experimental analysis necessarily constitute ensemble-averaged quantities, into which such small heterogeneities are subsumed.

FIG. 1.

Mechanisms included in the master equations [(1) and (2)]. We consider multiple effects for the formation of aggregates into our systems with rate constants ki. The constants corresponding to the transfer of mass to and from the external system are represented by Greek letters (γ and λi).

FIG. 1.

Mechanisms included in the master equations [(1) and (2)]. We consider multiple effects for the formation of aggregates into our systems with rate constants ki. The constants corresponding to the transfer of mass to and from the external system are represented by Greek letters (γ and λi).

Close modal

The aggregation model [(1) and (2)] and its many variations have been used to describe the results of many in vitro experiments.10,15,20 The analyses of multiple experimental datasets have shown that the exponents nc = n2 = 2 for Aβ40 and for Aβ42 and that these processes show little sign of saturation at low μM concentrations of protein when aggregating in phosphate-buffered saline (PBS) solution.10,20 In contrast, for Aβ42 aggregating in HEPES buffer, primary nucleation appears to be fully saturated and its rate shows no monomer dependence.15 While n2 = 2 still holds under those conditions, primary nucleation appears to proceed purely via heterogeneous primary nucleation that is fully saturated, i.e., KPm. In this limit, khρ(m) → khKP. For the discussions and derivations in this manuscript, we take the view of Aβ experiments10,20 so that nc = n2 = 2 and that heterogeneous primary nucleation, if it is present, is fully saturated, i.e., khρ(m) → khKP = k0, where we have defined k0, the rate of heterogeneous primary nucleation in the fully saturated regime, for simplicity. Adaptation to other aggregating systems is straightforward, and all numerical results are qualitatively similar.

The primary purpose of this manuscript is to describe the qualitative impact of clearance mechanisms on the dynamics of protein aggregation. Because primary nucleation is generally not a dominant process in the aggregation of disease associated proteins, the particular choice of nucleation mechanism, such as heterogeneous vs homogeneous, does not affect the results. Examples of fitted Aβ model parameters are listed in Table I. PBS and HEPES refer to the buffers used in the corresponding experiments. In these experiments, aggregation proceeds much faster than depolymerization and koff = 0 is found to be a good fit to describe the dynamics. However, from a theoretical point of view, we note that koff = 0 violates detailed balance and implies that there is no non-vanishing stationary distribution in the absence of clearance and production terms. Here, we will first follow experimental data and take koff = 0. Then, we will show that the addition of this small term does not change our results. Therefore, we will use the fitted experimental parameters given in Table I. Clearance and production have not been investigated experimentally; thus, we leave them as free parameters. Notably, a primary contribution of the current work is in determining particular values of these parameters when a qualitative change of the dynamics occurs.

TABLE I.

Typical best-fit parameters for the model. PBS and HEPES refer to the buffer used in the experiments.

ParametersMechanismAβ40 PBS20 Aβ42 PBS10 Aβ42 HEPES21 Units
k0 Heterogeneous nucleation 1.60 × 10−11 M h−1 
kn Homogeneous nucleation 5.80 × 10−3 1.20 × 10−1 M1nch1 
nc Homogeneous nucleation Unitless 
k2 Secondary nucleation 1.10 × 107 3.60 × 107 2.10 × 1014 M−2h−1 
n2 Secondary nucleation Unitless 
KM Saturation 3.60 × 10−11 3.60 × 10−12 2.30 × 10−17 M2 
k+ Elongation 1.10 × 109 1.10 × 1010 1 × 1010 M−1h−1 
koff Depolymerization h−1 
m0 Initial monomer c. 3 × 10−6 3 × 10−6 3 × 10−6 
λcrit Critical clearance 0.72 2.45 17.0 h−1 
λ̃crit Perfect bifurcation 0.72 2.47 17.0 h−1 
α Nonlinear coefficient 312 042 647 390 2.837 26 ×106 M−1 h−1 
τ1 Exponential time scale 1.4 0.4 0.06 
τ2 Amplification time scale 12.6 2.5 0.4 
λcrit(1) Critical clearance ν = 1 7.8 × 10−5 9.2 × 10−5 4.8 × 10−3 h−1 
λcrit(0) Critical clearance ν = 0 0.72 2.47 17.0 h−1 
λcrit(1) Critical clearance ν = −1 13.2 × 103 13.2 × 104 12 × 104 h−1 
ParametersMechanismAβ40 PBS20 Aβ42 PBS10 Aβ42 HEPES21 Units
k0 Heterogeneous nucleation 1.60 × 10−11 M h−1 
kn Homogeneous nucleation 5.80 × 10−3 1.20 × 10−1 M1nch1 
nc Homogeneous nucleation Unitless 
k2 Secondary nucleation 1.10 × 107 3.60 × 107 2.10 × 1014 M−2h−1 
n2 Secondary nucleation Unitless 
KM Saturation 3.60 × 10−11 3.60 × 10−12 2.30 × 10−17 M2 
k+ Elongation 1.10 × 109 1.10 × 1010 1 × 1010 M−1h−1 
koff Depolymerization h−1 
m0 Initial monomer c. 3 × 10−6 3 × 10−6 3 × 10−6 
λcrit Critical clearance 0.72 2.45 17.0 h−1 
λ̃crit Perfect bifurcation 0.72 2.47 17.0 h−1 
α Nonlinear coefficient 312 042 647 390 2.837 26 ×106 M−1 h−1 
τ1 Exponential time scale 1.4 0.4 0.06 
τ2 Amplification time scale 12.6 2.5 0.4 
λcrit(1) Critical clearance ν = 1 7.8 × 10−5 9.2 × 10−5 4.8 × 10−3 h−1 
λcrit(0) Critical clearance ν = 0 0.72 2.47 17.0 h−1 
λcrit(1) Critical clearance ν = −1 13.2 × 103 13.2 × 104 12 × 104 h−1 

One key fact to consider before we embark on solving the dynamics is what boundary conditions, in the form of protein homeostasis, we impose on the system. In the aggregation of purified protein in vitro, there are no sources or sinks for protein mass, so the total mass of protein is conserved. In the present case, a living system, both sources and sinks exist and different possibilities arise. Generally, organisms will strive to maintain stable conditions and balance these sources and sinks. We consider two extreme cases and show that our conclusions hold in both regimes and are thus likely to apply in a range of biological systems. First, we will consider a system in which the monomer production rate is constant over time (see Secs. III and IV). Second, we will consider a system where instead the monomer concentration is constant in time, so the monomer production rate may vary (see Sec. V). Most biological systems are likely to lie somewhere between these two extremes.

In the case of size-independent clearance, we have λi = λ > 0 for all i. Our main question is to understand the role of the clearance term. In particular, we will establish that if clearance is sufficiently large, the formation of aggregates does not take place. We will then estimate this critical clearance value, denoted λcrit, under the different experimental conditions considered.

In the size-independent case, a previously established11 but remarkable feature of the system (1) and (2) is that a closed system of equations for the first two moments P and M and the monomer concentration m can be obtained as

dPdt  =λP+k0+knm2+k2σ(m)M,
(4)
dMdt=λM+2k0+2(k+mkoff)P+2knm2+2k2σ(m)M,
(5)
dmdt=γλm2k02(k+mkoff)P2knm22k2σ(m)M,
(6)

where σ(m) = m2KM/(KM + m2), and we have chosen n2 = 2. The total mass of the system Mtot = M + m satisfies, by summing (5) and (6), the evolution equation

dMtotdt=λMtot+γ.
(7)

This equation implies that the total mass in the system evolves to a stable steady state Mtot = γ/λ with a typical time scale 1/λ. To simplify the analysis, we will further assume that, initially, the system has already reached this state before the dynamics of aggregation starts. To do so, we chose the following unseeded initial conditions:

M(0)=P(0)=0,m(0)=m0=γ/λ,
(8)

and thus, the total mass of the system is conserved for all time Mtot(t) = m0. The term “unseeded” refers to the fact that, initially, there is no aggregated protein in the system (hence, no seed). This condition assumes a lack of aggregated species in a healthy in vivo state. Indeed, it is observed that soluble Aβ monomers are found in healthy individuals of all ages, while aggregates larger than monomers are generally correlated with Alzheimer’s disease progression.22 An extra advantage of this approach is that it fixes the constant γ = m0λ.

Before we study the system in full generality, it is useful to consider the overall dynamics of the system for a typical set of parameters for the aggregation of Aβ40 given in the first column of Table I. We will use this set of parameters for all our examples. The other datasets are qualitatively equivalent, and the values of various derived quantities are given in Table I. As shown in Fig. 2, the typical behavior of the system from an unseeded initial condition is for the aggregate mass to increase up to a finite value M, while the monomer concentration decreases to m, in a typical sigmoid-like behavior. We observe that, in the absence of clearance, the monomer population is completely converted to aggregates (λ = 0, dashed curves in Fig. 2). Conversely, for large clearance, almost no conversion takes place (λ = 1, dotted curves in Fig. 2). Some of the monomers are converted (solid curves for λ = 0.2 in Fig. 2) for the case of moderate clearance. Of particular interest, for our discussion, is the change of behavior at some critical value λcrit of the clearance λ where aggregation becomes negligible.

FIG. 2.

Typical dynamics of the monomer (blue) and aggregate (red) concentration (in units of micromolar) for different values of the clearance (λ in h−1), evaluated using the parameters for Aβ40 from Table I and λ = 0 (large dashed), λ = 0.2 (solid), and λ = 1 (small dashed). Asymptotic values for λ = 0.2 are shown with dotted lines.

FIG. 2.

Typical dynamics of the monomer (blue) and aggregate (red) concentration (in units of micromolar) for different values of the clearance (λ in h−1), evaluated using the parameters for Aβ40 from Table I and λ = 0 (large dashed), λ = 0.2 (solid), and λ = 1 (small dashed). Asymptotic values for λ = 0.2 are shown with dotted lines.

Close modal

To derive a value for λcrit, we determine the dependence of the asymptotic states m on λ. Considering steady state m = m, P = P, and M = M in (4) and (5), one expresses the latter two states as a function of parameters, λ and m. These relations are substituted into (6) to produce the implicit equation H(m, λ) = 0 with

H(m,λ)=2k+knm5m42k+k2KM2λkn+2knkoffm3λ22k+k2m0KM+2k2λKM2k+knKM2k2koffKM2k+k0m22k0λ2k2λm0KM+2k2m0koffKM2λknKM+2knkoffKM+2k0koff+λ2m0+mλ2KM+2k+k0KM+2λk0KM2k0koffKMλ2m0KM.
(9)

For instance, for the same parameter values as in Fig. 2, we show in Fig. 3 the values of m as a function of λ. We observe a sharp transition for a critical value of the clearance parameter λ. There are three necessary conditions for λcrit: first that λcrit is non-negative, second that m is maximal, and third that the value of m coincides with m0. The last two conditions can be realized by computing the derivative of the expression H(λ, m) = 0 evaluated at m = m. Therefore, λcrit is given by the positive root of L(λ) = 0, where

L(λ)=Hmm=m0=4koffKMmk2m0+kn+6k+knKMm2+6k2koffKMm28k+k2KMm3+6k+k2m0KMm2+2k+k0KM8knkoffm3+10k+knm44k0koffm+6k+k0m2+λ4KMmk2m0+kn6k2KMm2+8knm3+4k0m+λ2KM+3m22m0m.
(10)

For Aβ-40, the critical clearance, as shown in Fig. 3, is λcrit = 0.72 h−1. Critical clearance rates for the other experimental datasets are given in Table I for comparison.

FIG. 3.

Perfect (red) and imperfect (blue) transcritical bifurcation, shown for the parameters for Aβ40. Unstable (dashed) and stable (solid) equilibrium solutions. For clearance faster than the critical value, λ > λcrit, essentially all protein is present in monomeric form, whereas for lower clearance rates, λ < λcrit, a significant fraction of the total protein is aggregated. In this case, we have for the perfect bifurcation λ̃crit0.72 and m̃0.7+3.2λ̃. The dashed curves indicate unstable equilibrium solutions, and the solid curves denote stable equilibrium,

FIG. 3.

Perfect (red) and imperfect (blue) transcritical bifurcation, shown for the parameters for Aβ40. Unstable (dashed) and stable (solid) equilibrium solutions. For clearance faster than the critical value, λ > λcrit, essentially all protein is present in monomeric form, whereas for lower clearance rates, λ < λcrit, a significant fraction of the total protein is aggregated. In this case, we have for the perfect bifurcation λ̃crit0.72 and m̃0.7+3.2λ̃. The dashed curves indicate unstable equilibrium solutions, and the solid curves denote stable equilibrium,

Close modal

In a neighborhood of λcrit, m, as a function of λ, undergoes a sharp transition. This transition is not a bifurcation in the strict sense, but, in the parlance of dynamical systems, it can be described as an imperfect transcritical bifurcation when heterogeneous nucleation and homogeneous nucleation terms can be understood as an imperfection and are sufficiently small with respect to the elongation. More specifically, when k0/(k+m02)1 and kn/k+ ≪ 1, the system is well approximated by k0 = 0 and kn = 0. In this limiting case, the fixed point (P, M, m) = (0, 0, m0) for the system (4)(6) undergoes a (perfect) transcritical bifurcation at λ̃critλcrit that can be obtained by locally expanding m in λ to find

m̃=m0+1α(λλ̃crit)+O(λλ̃crit)2,
(11)

where λ̃crit is specified by the formula

λ̃crit=m0k2KMm0k2m0+2k+KM2koffKM+2m02k+m0koff+k2m0KMKM+m02
(12)

and α is defined by the expression

α=m0k2KM2λ̃crit+3k+m02koffλ̃crit2λ̃critKM+m02k2m02KM.
(13)

When the clearance is close to the critical value, the linear approximation to the perfect bifurcation is a reasonable approximation for the imperfect bifurcation, as can be appreciated in Fig. 3 where λ̃crit0.72 and m̃0.7+3.2λ̃. By analogy with epidemiology, we define a dimensionless neurodegenerative reproduction number,

R0=λ̃critλ,
(14)

such that for R0 < 1 the aggregate level is negligible and grows to a finite value for R0 > 1.

The existence of a critical clearance rate shows that in the healthy regime, i.e., for sufficiently large values of clearance, the system (1) and (2) with size-independent clearance can support a small, endemic, population of aggregates. The formation of a significant aggregate population, in this case, occurs only when the system’s clearance rate, λ, drops sufficiently below the critical clearance rate λcrit. We can explore the dynamics close to the bifurcation by considering the normal form of the system for the perfect system [(4)(6) with k0 = kn = 0] near λ=λ̃crit. The general method to obtain the normal form of a transcritical bifurcation for an arbitrary smooth vector field is given in  Appendix A. Applying these ideas, we can approximate the full system by

Ṗ=(λλ̃crit)P+αvPP2,
(15)
Ṁ=(λλ̃crit)MαM2,
(16)
ṁ=(λλcrit)(mm0)+α(mm0)2,
(17)

where α is given by (13) and

vP=λ̃crit2k+m0koff+λ̃crit.
(18)

Figure 4 shows a comparison of the total aggregate mass evolution, vs time, obtained for the imperfect unseeded system, the perfect seeded system, and the normal form. As expected, the agreement is excellent as long as the system is close enough to the bifurcation point.

FIG. 4.

Aggregate mass concentration M(t) as a function of time for the unseeded (black dashed) system (4)(6), the perfect seeded system (neglecting homogeneous and heterogeneous primary nucleation) (red), and the normal form approximation of M(t) (dotted blue). The initial conditions were selected so that the initial growth rates matched the initial growth rate of the unseeded system. The red curve was generated with unseeded initial conditions, the black dashed curve was computed using seeded compatible initial conditions given by (P(0), M(0), m(0)) = (S, S/2, m0S) with S = 2.4 × 10−13, and the blue dashed curve was generated by solving (16) with M(0) = S/2. The parameters are for the Aβ40 values of Table I and λ = 1/2.

FIG. 4.

Aggregate mass concentration M(t) as a function of time for the unseeded (black dashed) system (4)(6), the perfect seeded system (neglecting homogeneous and heterogeneous primary nucleation) (red), and the normal form approximation of M(t) (dotted blue). The initial conditions were selected so that the initial growth rates matched the initial growth rate of the unseeded system. The red curve was generated with unseeded initial conditions, the black dashed curve was computed using seeded compatible initial conditions given by (P(0), M(0), m(0)) = (S, S/2, m0S) with S = 2.4 × 10−13, and the blue dashed curve was generated by solving (16) with M(0) = S/2. The parameters are for the Aβ40 values of Table I and λ = 1/2.

Close modal

Next, we consider the effect of clearance on size distribution. First, we take koff = 0, which the data suggest is a good approximation to the system unless m approaches 0. Since we are interested in the asymptotic size distribution, we can assume that m = m in Eqs. (1) and (2), in which case, we have simply that

pi=2k+mλ+2k+mpi1=δ0pi1,pi=δ0i2p2,  i>2.
(19)

Using the definition of M = i>1ipi, we obtain

p2=M(1δ0)22δ,pi=Mδ0i2(1δ0)22δ0  i>2.
(20)

This analysis is not valid for λ → 0. In that case, the total mass of the system is systematically transferred to larger and larger particles, and in the long-time limit, all finite aggregate concentrations tend to vanish and the trivial distribution is pi = pi−1 = p2 = 0. However, in that limit, the assumption koff = 0 is not justified anymore as even a small value of koff allows for a non-trivial size distribution. Indeed, with koff ≠ 0, we have the following recurrence relation for ci:

0=λipi+2k+m(pi1pi)+2koff(pi+1pi),i>2,
(21)

with a single bounded solution of the form

pi=δi2p2,  i>2
(22)

with

δ=k+m2koff+12+λ4koff122koff+2k+m+λ  24koff24k+mkoff.
(23)

An asymptotic expression of δ for small and large values of koff gives

δ=δ012λkoff2mk++λ  2+Okoff2for   koff<mk+δ02mk++λ2koff+Okoff2for   koff>mk+.
(24)

We see that unless λ = 0, the role of koff, when sufficiently small, is negligible. We conclude that clearance (or depolymerization) is sufficient to obtain a non-degenerate size distribution (Fig. 5).

FIG. 5.

The effect of the parameter koff on the size distribution can be appreciated by computing δ/δ0 as a function of koff. We see that for koff < mk+, the role of koff is negligible. The dashed curves are given by the asymptotic approximation (24). The parameters are for the Aβ40 values of Table I and λ = 1/2.

FIG. 5.

The effect of the parameter koff on the size distribution can be appreciated by computing δ/δ0 as a function of koff. We see that for koff < mk+, the role of koff is negligible. The dashed curves are given by the asymptotic approximation (24). The parameters are for the Aβ40 values of Table I and λ = 1/2.

Close modal

Next, we generalize the above treatment and allow the clearance rate of an aggregate to depend on its size. In this case, there is no simple, closed equation for the moments, as in Sec. III, and we must study the full system. Here, we make a key assumption about the dependence of the clearance on the aggregate size. We assume that there exists a critical aggregate size, N, such that all aggregates of size N, or greater, are too large to be cleared. Explicitly, this assumption implies that λi = 0, ∀iN. We also assume that koff = 0 and nc = n2 = 2, and then (1) and (2) can be written as

dM̃dt=i=2N1λiipi+2k0+2knm2+2k+mP+2k2σ(m)M̃,
(25)
dmdt=λ1(m0m)2k02knm22k+mP2k2σ(m)M̃,
(26)
dp2dt=λ2p2+k0+knm22k+mp2+k2σ(m)M̃,
(27)
dpidt=λipi+2k+m(pi1pi),  i>2,
(28)

where P=i=2pi and M̃=i=2ipi. The unseeded initial conditions for this system are

m(0)=m0,pi(0)=0 for2i,  M̃(0)=0.
(29)

In general, under the assumption of a constant monomer production rate, there is no guarantee of mass conservation. For instance, if λiλ1∀i > 2 and there is at least one i ≥ 2 such that λi < λ1, then the overall mass of proteins will increase in time as shown in  Appendix B.

To study the dynamics of (25)(28), we introduce a finite system with equivalent dynamics. Here, we follow Ref. 23 (see also Ref. 24) and introduce a super-particle, denoted qN, which represents the concentration of all aggregates of size greater than or equal to N,

qN=i=Npi.
(30)

Since λi = 0 for all iN, we can take the limit of the partial sums of (28) to obtain

dqNdt=i=N2k+m(pi1pi)=limji=Nj2k+m(pi1pi)=2k+mpN12k+limjmpj.
(31)

Since the monomer concentration m remains bounded, for any fixed time, the last term of (31) tends to zero as j and the super particle concentration satisfies the equation

dqNdt=2k+mpN1.
(32)

We will distinguish the finite system with a super-particle from the infinite system (25)(28) by introducing the notation qi = pi for i < N. Defining Q=i=2Nqi, and using (32), the corresponding super-particle system is defined by

dMdt=i=2N1λiiqi+2k0+2knm2+2k+mQ+2k2σ(m)M,
(33)
dmdt=λ1(m0m)2k02knm22k+mQ2k2σ(m)M,
(34)
dq2dt=λ2q2+k0+knm22k+mq2+k2σ(m)M,
(35)
dqidt=λiqi+2k+m(qi1qi),i=2,,N1,
(36)
dqNdt=2k+mqN1.
(37)

The unseeded conditions for (33)(37) are

m(0)=m0,qi(0)=0,  for2iN,M(0)=0.
(38)

For unseeded initial conditions, the dynamics of the finite system is equivalent to the infinite one in the following sense: First, note that Q̇=Ṗ; this follows directly from the definition of Q, P, and (30). Thus, Q and P will agree for all time. In turn, (25) and (33) coincide when the initial data (29) and (38), respectively, are used; thus, M̃(t)=M(t) in this case. Finally, by definition, pi = qi for 2 ≤ i < N and (30) and (31) have already established that solving (37) produces qN(t)=i=Npi(t), provided that the initial conditions agree. The above establishes an important fact that we rely on for the rest of the section; solving (25)(28) with initial conditions (29) and (33)(37) with initial conditions (38) yields

Q(t)=P(t),M(t)=M̃(t),pi(t)=qi(t) for 2i<N and qN(t)=i=1pi(t).
(39)

We remark, however, that M(t), defined as the solution of (33), is the total aggregate mass of both (25)(28) and (33)(37), due to (39), for the unseeded initial conditions (38); however, M(t) cannot be constructed a posteriori from the knowledge of qi(t) where i = 2, 3, , N in the same manner that M̃(t) can be retrieved from the knowledge of pi(t). That is, we have M(t)i=2Niqi(t). Indeed, in the closure process of reducing the full system to a finite one, we lost information regarding the mass of individual particles making up the superparticle. Nevertheless, both the evolution of aggregate mass of the full system and the size distribution (up to size N) can be obtained by studying the finite system (33)(37).

Given a constant monomer production rate, systems such as (25)(28) or (33)(37), with size-dependent clearances, do not conserve mass in general (see  Appendix B) and the aggregate mass may increase with time. We study in more detail the particular choice

λi=λ/i fori=1,2,,N1,
(40)

which expresses the modeling assumption that aggregates become increasingly difficult to clear as their size increases. An example of the dynamics of the system (33)(37) is shown in Fig. 6. We observe two different behaviors. Initially, up to a time τ2, the system mostly behaves like the conservative no-clearance model (λ = 0) even for large values of clearance. This behavior is markedly different than the one observed in Fig. 2. Second, for larger times, t > τ2, the monomer mass always decreases and the aggregate mass always increases as predicted from our general analysis. We observe that larger clearance leads to faster aggregate mass creation. This is due to the fact that we balance monomer clearance by monomer production, so an increase in clearance automatically implies an increase in production. The question is then to understand the transition between the two regimes as well as the small and large time behaviors of all species. The situation where the clearance rate is higher relative to the production rate is shown in Fig. 7. Biologically, this corresponds, for example, to two individuals in which production proceeds at the same rate, but clearance is more efficient in one than in the other. While there is no qualitative change in behavior, i.e., the total concentration is still unbounded, both the speed at which aggregates accumulate and the steady state monomer concentration are lower in the system with a higher clearance rate. Thus, in a system with constant (i.e., time-invariant) monomer production rate and a clearance rate that decreases with size, whether increasing clearance will be beneficial or not will depend on whether it causes a simultaneous increase in the rate of monomer production or not. In either case, the aggregate mass is not bounded and only its rate of increase can be affected by changes in the clearance rate.

FIG. 6.

Aggregate mass dynamics for the size-dependent clearance λi = λ/i. The monomer concentration [m(t), blue lines] and total aggregate mass [M(t), red lines] are shown in (a) for clearance rates (in h−1): λ = 0 (dashed), λ = 0.2 (solid), and λ = 1 (dotted) and monomer production rates (γ = m0λ) and in (b) for clearance rates (in h−1): λ = 0.2 (solid) and λ = 1 (dotted) and a monomer production rate that is the same for both curves (γ = 0.2m0). The parameters are for the Aβ40 values of Table I, λ = 1/2 and N = 20.

FIG. 6.

Aggregate mass dynamics for the size-dependent clearance λi = λ/i. The monomer concentration [m(t), blue lines] and total aggregate mass [M(t), red lines] are shown in (a) for clearance rates (in h−1): λ = 0 (dashed), λ = 0.2 (solid), and λ = 1 (dotted) and monomer production rates (γ = m0λ) and in (b) for clearance rates (in h−1): λ = 0.2 (solid) and λ = 1 (dotted) and a monomer production rate that is the same for both curves (γ = 0.2m0). The parameters are for the Aβ40 values of Table I, λ = 1/2 and N = 20.

Close modal
FIG. 7.

Aggregate mass dynamics for the size-dependent clearance λi = λ/i for various clearance efficiencies, λ. The monomer concentration [m(t), blue dashed lines] and total aggregate mass [M(t), red dashed lines] are shown for various values of λ under a constant production, γ = 0.6, assumption [whereby γ/λ1 = γ/λ = m(0)].

FIG. 7.

Aggregate mass dynamics for the size-dependent clearance λi = λ/i for various clearance efficiencies, λ. The monomer concentration [m(t), blue dashed lines] and total aggregate mass [M(t), red dashed lines] are shown for various values of λ under a constant production, γ = 0.6, assumption [whereby γ/λ1 = γ/λ = m(0)].

Close modal

On long time scales, i.e., long enough so that the monomer concentration begins to decrease, the monomer production, aggregation, and nucleation processes result in an increase in subsequent aggregated species and, therefore, in the overall aggregate mass M. The asymptotic behavior of the system aggregate mass M is observed to depend entirely on the production rate, γ = λm0, as

M(t)tγt.
(41)

This behavior is illustrated in a log-plot in Fig. 8; the characteristic time scale, τ2, indicates the time at which the monomer mass begins to decay. Once the asymptotic behavior of M has been established, the equations can be balanced asymptotically by the following dynamics:

mαmt2/3,  qNαNt2/3,qiαit1/3,i=1,,N1,
(42)

where the symbol “∼” is understood as the long-time asymptotic behavior and αi are constants. This asymptotic behavior shows that the super-particle dominates the long-term dynamics; thus, PqN for large times. Physically, in the long-time limit, the monomer population, renewed by the continuous production, is quickly promoted to the super-particle through linear aggregation.

FIG. 8.

Long time (in h) concentration (in mol) dynamics of (33)(37) (λi = λ/i, in h−1) with N = 20 and Aβ40 parameters (Table I, third column); curves for λ = 0.2 (solid) with asymptotic slopes (dotted). Time runs from 0 h to 5000 h. The parameters are for the Aβ40 values of Table I and N = 20.

FIG. 8.

Long time (in h) concentration (in mol) dynamics of (33)(37) (λi = λ/i, in h−1) with N = 20 and Aβ40 parameters (Table I, third column); curves for λ = 0.2 (solid) with asymptotic slopes (dotted). Time runs from 0 h to 5000 h. The parameters are for the Aβ40 values of Table I and N = 20.

Close modal

We observe in Fig. 6 that the early time behavior is not greatly perturbed by altering the clearance rate. Hence, we can obtain characteristic time scales for the amplification of the aggregate mass by considering the limit λ → 0+. In this case, the early evolution of the aggregate mass is governed by the dynamics of (4)(6) with λ = 0. There are two characteristic time scales of importance. First, the time scale τ1 associated with the exponential growth of the aggregate mass in early time via the inverse of the positive linear eigenvalue, μ = 1/τ1, corresponding to the linearization of (4)(6) around the healthy state m = m0, M = P = 0. The linear eigenvalue is given by the positive root of

μ2+μ4m0kn2k2m02KMKM+m022k+k2m03KMKM+m02+4k+m02kn=0.
(43)

Second, there is a time scale τ2 where both nucleation and amplification are balanced. It is given by the time for the linearized solution for M(t) to reach m0. Hence, τ2 is the solution of

m02=m02kn+k0KM+m022knKM+2m02knk2m0KM1eτ2/τ12.
(44)

For example, for the first parameter set (Aβ40) used for Figs. 68, these times are τ1 ≈ 1.4 h and τ2 ≈ 12.6 h. The value of τ2 is a rudimentary estimate for the time of amplification; it is a lower bound for the typical time scale of growth (see Fig. 6). Nevertheless, in Fig. 8, we see that τ2 can indeed act as an indicator for the onset of decay for the monomer mass. A more refined estimate can be obtained by using the approximate solution for the full dynamics given in Ref. 14.

In Secs. III and IV, we investigated a system in which the monomer production rate was constant, i.e., invariant with time. The other biologically reasonable assumption is that the monomer concentration instead remains at a constant level m0. In this case, the system adjusts the monomer production rate to compensate for any loss of monomer, through degradation but also through aggregation. Experimental observations for example of the concentration of PrP, the main protein that aggregates in prion disease, in mice with prion disease support such a constant monomer.25 We assume that, regardless of other parameters, (1) is instead specified by

dmdt=0
(45)

so that, with unseeded initial conditions, we have m(t) = m0 for all time. Assuming again no depolymerization, no fragmentation, and dimer nucleation, the master equations now read

dp2dt=λ2p2+k0+knm022k+m0p2+k2σ0M,
(46)
dpidt=λipi+2k+m0(pi1pi),i>2,
(47)

where σ0 = σ(m0) and M=i=2ipi is the total aggregate mass. This is an infinite system of linear ordinary differential equations. For this system, we consider three types of clearance: the size-independent case in addition to two different size-dependent paradigms. All three clearance relations can be summarily presented by a power-law of the form

λi=λiν.
(48)

When ν = 0, we recover the size-independent case; when ν = −1, we recover the size-dependent diminishing clearance formulation used in Sec. IV; and finally, the case of ν = 1 corresponds to improved clearance, with increasing size, which could arise due to, for instance, antibody binding. Depending on the two parameters λ and ν, the solution to this system may have a steady state or increase indefinitely. The question is then to identify the critical values at which this transition happens.

We start with the simple case of size-independent clearance ν = 0; this is the analog to Sec. III for a constant free monomer assumption [cf. (45)] The moments (cf. Sec. III) are specified by a simple pair of linear equations given by

dPdt  =λP+k0+knm02+k2σ0M,
(49)
dMdt=λM+2k0+2k+m0P+2knm02+2k2σ0M,
(50)

which can be written as

q̇=Aq+b,
(51)

where q = (P, M)T, b=(k0+knm02,2k0+2knm02)T, and

A=λab2aλ=λk2σ02k+m02k2σ0λ.
(52)

The constant solution sole steady state for this system is q = −A−1b; q is positive and finite if

λ>a+a2+ab=k2σ0+k22σ02+2k2σ0k+m0=λcrit(0).
(53)

This condition naturally provides a value for the critical clearance. Specifically, the largest linear eigenvalue for the system is κ=λcrit(0)λ; solutions converge to q exponentially in time (as eκt) for λ>λcrit(0) and grow unbounded for λλcrit(0) (see Fig. 9).

FIG. 9.

The accumulation of aggregate mass in the case of a constant monomer concentration and a constant clearance rate. To illustrate the stability of this system to the introduction of aggregates, the curves are shown for a reaction with 10% seed, i.e., M(0) = 0.1m0 and P(0) = M(0)/1000. The same qualitative behavior emerges also in the absence of seeds, but differences appear over longer time scales. The curves are shown for the system without clearance (large dashed) and for values of the clearance rate just above (small dashed) and just below (solid) its critical value. The parameters are for the Aβ40 values of Table I.

FIG. 9.

The accumulation of aggregate mass in the case of a constant monomer concentration and a constant clearance rate. To illustrate the stability of this system to the introduction of aggregates, the curves are shown for a reaction with 10% seed, i.e., M(0) = 0.1m0 and P(0) = M(0)/1000. The same qualitative behavior emerges also in the absence of seeds, but differences appear over longer time scales. The curves are shown for the system without clearance (large dashed) and for values of the clearance rate just above (small dashed) and just below (solid) its critical value. The parameters are for the Aβ40 values of Table I.

Close modal

The values given in Table I for the different parameters show that this estimate is indistinguishable from the case studied in Sec. III, which is explained by the fact that at the bifurcation point, the monomer population is constant in both cases.

We now turn our attention to the general case where the clearance terms are no longer size-independent. Then, the master equations do not yield a closed system for the moments. Nevertheless, due to the simplicity introduced by m(t) = m0 being constant, we can find conditions for the existence of a fixed-point solution, (p2*,p3*,) to (46) and (47). If such a steady state pi* for i > 2 exists, it must satisfy the recurrence relation

pi*=δipi1*,  δi=bb+λi=2k+m02k+m0+λi.
(54)

We note that each of the recursion coefficients, δi, is now dependent on i via λi. Define a sequence of real numbers, indexed by i, as

Δi=j=3iδj,i>2.
(55)

We define Δ2 = 1, and the ith steady state is expressible, for all i ≥ 2, through its recurrence relation as

pi*=Δip2*,i2.
(56)

Defining

Δ=j=3jΔj,
(57)

the steady state for the total aggregate mass solution M* is then given by

M*=i=2iΔip2*=Δp2*,
(58)

and an application of (46), at steady state, gives the value of p2* as

p2*=k0+knm02λ2+2k+m0k2σ0Δ.
(59)

Therefore, for a fixed point to exist, we need the three following conditions to be satisfied:

C1: limiΔi=0,
(60)
C2: Δ=i=2iΔi converges,
(61)
C3: k2σ0Δλ22k+m0<0.
(62)

An analysis of the case ν = 0 recovers the previous condition, and it can then be verified directly that conditions C1–C3 are satisfied, as expected, for λ>λcrit(0).

1. Enhanced clearance: ν = 1

For ν = 1, we have (see  Appendix C) Δ(1) = 2 + b/λ, and the steady population of dimers, whenever it exists, is given by

p2*=λ(k0knm02)2(k+m0+λ)(λk2σ0).
(63)

Hence, condition C3 leads to λ>λcrit(1) with

λcrit(1)=k2σ0.
(64)

We note that the above implies that the critical clearance depends only on the secondary nucleation process and, in particular, not the process of elongation [cf. λcrit(0) in (53)].

2. Reduced clearance: ν = 1

For ν = −1, the situation is not as simple. Condition C1 is verified, but C2 leads to λ > 2b for which

Δ(1)=λ+2bΓλb2Γλb+2Γλb22bΓλb2,
(65)

where Γ(·) is the usual Gamma function. Condition C3 is satisfied if λ>λcrit(1) where λcrit(1) is the positive solution of

fλb=1+2k+m0k2σ0, with fλb=Γλb2Γλb+2Γλb2.
(66)

This equation always has a solution as f: z ∈ [2, ] → f(z) is such that f′(z) < 0, f(z)z2, and f(z) →z1. For the parameters listed in Table I, 2k+m0/k2σ0 ≫ 1, in which case, we can approximate the function f(z) close to z = 2 by f(z) ≈ 6/(z − 2), which leads to the critical value

λcrit(1)=8k+m0(2k2σ0+k+m0)k2σ0+2k+m0.
(67)

This last relation can be further simplified by realizing that k+m0k2σ0, which leads to

λcrit(1)=4k+m0.
(68)

For the parameters given in Table I, this last approximation of the critical clearance gives the correct value [compared to (66)] to six digits. Note that, in contrast to the critical clearance rate for enhanced clearance [cf. (64)], (67) depends only on the elongation rate k+. In particular, in a reduced clearance regime, a change in the rate of secondary nucleation has no effect on the clearance rate required to keep the system stable. The general trend that can be observed from Table I is that λcrit(1)>λcrit(0)>λcrit(1), as expected.

3. Further reduced clearance: ν = 2

Finally, for ν = −2, skipping computational details, we find that

limnΔn(2)=14πλbλb2+5λb+4cschπλb,
(69)

which is positive for all finite positive values of λ. Hence, condition C1 is not satisfied and there is no constant solution or critical value of the clearance that would limit unbounded growth of the aggregate mass. We note that we have neglected the effect of fragmentation. For ν < 0, the effect of fragmentation is the creation of smaller aggregates that increase the overall expansion of the protein population but also boosts clearance. Indeed, since smaller aggregates are more likely to be cleared and we expect a reduction in the critical value of clearance as well as the possibility of a finite value of clearance for ν = −2 or smaller, as shown in the work of Meisl.26 Comparing the different critical clearance values given in Table I for the three values of ν, it is clear that that the choice of clearance law has a significant impact on the clearance values as they differ, from the smallest to the largest, by nine orders of magnitude. Hence, enhancing or inhibiting the clearance mechanism may be extremely important to the overall increase in aggregate mass.

We have assessed the impacts of monomer production and aggregate clearance on the kinetics of protein aggregation using a theoretical model, cf. (1) and (2), based on previous work that has been very successful in describing experimental data.10,20,21 Our findings suggest that clearance, regardless of the specific assumptions of the system, may fundamentally modulate aggregation kinetics (see Table II). We consider the two extreme regimes of constant monomer production and constant monomer concentration, which are both biologically reasonable and mathematically tractable. In most real biological systems, the situation is likely to be somewhere between these two extremes, with some feedback mechanisms adjusting the protein expression rate, but monomer concentrations not staying completely constant either. Which regime is more representative of the real system will also depend on the protein under consideration: For proteins that fulfill a functional role, such as the tau protein, the situation of a constant monomer concentration may be a more biologically relevant one, whereas for peptides that are essentially waste products, such as Aβ, the assumption of a constant production rate may be closer to reality.

TABLE II.

Summary of the behavior and dependence of the critical rate on the system parameters for the different regimes.

MonomerSize-dependenceCritical rateLong time
Constant production λi = λ Equation (12) Always bounded 
 λi = λ/i Does not exist Not bounded 
Constant concentration λi = λ Equation (53) Bounded if λ < λcrit 
 λi = λ/i Equation (67) Bounded if λ < λcrit 
 λi =  Equation (64) Bounded if λ < λcrit 
 λi = λ/i2 Only for fragmenting system Not bounded 
MonomerSize-dependenceCritical rateLong time
Constant production λi = λ Equation (12) Always bounded 
 λi = λ/i Does not exist Not bounded 
Constant concentration λi = λ Equation (53) Bounded if λ < λcrit 
 λi = λ/i Equation (67) Bounded if λ < λcrit 
 λi =  Equation (64) Bounded if λ < λcrit 
 λi = λ/i2 Only for fragmenting system Not bounded 

In the case of a size-independent clearance rate, we showed that aggregation is controlled, directly, by a critical clearance threshold, regardless of the assumptions about monomer production. Clearance above this level provides for a robust environment, which is, essentially, free of protein aggregates; clearance below this level triggers an instability and a propensity toward aggregate mass accumulation. When monomer production is constant, once aggregation is triggered, the soluble monomer population is diminished as further aggregates form and the maximal level of aggregate formation is mediated by the clearance level.

We then generalized this treatment by considering the reasonable in vivo hypothesis that the clearance can depend on the aggregate size i. First, we explored this clearance paradigm using a simple inverse proportionality law λi = λ/i, whose biological interpretation is that aggregates become more difficult to clear as they increase in size. The resulting set of equations for this type of clearance does not yield a finite system for the moments. Thus, a super-particle system, with identical trajectories in the presence of unseeded initial conditions, has been advanced as a means of study. In the presence of any aggregation effects, the system immediately begins accumulating aggregates, even from unseeded initial conditions. Moreover, when the monomer production rate is constant, the mass is not conserved and the aggregate mass grows in time, without bound. The clearance, however, determines the asymptotic rate of increase in the aggregate mass as a function of time with M(t) ∼ λm0t. Thus, the implications of such a size-dependent clearance are quite different than the size-independent case as aggregate accumulation can only be delayed, not entirely prevented, even in an unseeded system. By contrast, when the monomer concentration is constant, a critical rate also exists for clearance of the form λi = λ/i, although the dependence of its value on the other parameters of the system differs from that of size-independent clearance.

Finally, we also considered the case when the clearance rate increases with aggregate size, λi = , as one might encounter when the process is mediated by antibodies and the clearance rate of each aggregate size is dependent on the probability of at least one antibody recognizing the aggregate. A critical clearance rate can be determined, although again the dependence of its value on the other parameters of the system differs from that of other clearance mechanisms. Remarkably, these results therefore suggest that, depending on the specific size-dependence of clearance, the processes of elongation and secondary nucleation contribute to the value of the critical clearance to different degrees. In turn, this implies that a more detailed understanding of what clearance processes dominate is crucial in order to predict which aggregation process should be targeted for inhibition in order to reduce the critical clearance rate.

Overall, the role of clearance in aggregation kinetics is highly non-trivial. However, we show that clearance may play an important role in the aggregation kinetics of amyloid-β, determining whether the system is stable or enters a run-away aggregation state. Our work provides the theoretical and mathematical foundations for the study of clearance of protein aggregates in living systems. We envision that the extension of this approach and its application to the analysis of experimental data will help shed light on the complexities of protein aggregation in vivo and the factors that determine the resilience of a system toward aggregation.

The work of A.G. was supported by the Engineering and Physical Sciences Research Council (Grant No. EP/R020205/1). The work of T.B.T. was supported partially by the John Fell Oxford University Press Research Fund (Grant No. 000872) (Project No. BKD00160), and that of A.G. was supported partially by the Engineering and Physical Sciences Research Council (Grant No. EP/R020205/1).The work of G.M. and T.P.J.K. was supported by the Cambridge Home and EU Scholarship Scheme, the Cambridge Centre for Misfolding Diseases, and the Frances and Augustus Newman Foundation and the European Research Council under the European Union’s Seventh Framework Programme (Grant No. FP7/2007-2013) through the ERC grant PhysProt (Grant Agreement No. 337969).

The data that support the findings of this study are available within the article.

Here, we derive the normal form of a transcritical bifurcation for a general dynamical system. We consider an autonomous n-dimensional C2 vector field of the form

ẋ=f(x,λ),xRn,
(A1)

and assume that there exists a constant solution x0 such that f(x0, λ) = 0 and a different equilibrium solution in a neighborhood of the critical value λ0. The conditions for the existence of a transcritical bifurcation at the critical value λ0 are given by Sotomayor’s theorem,27 and the reduced form of the system takes close to that value can be captured by normal form theory.28–31 Here, we use multiple scale analysis to obtain a convenient form of the reduced equations. The result in itself is not original, but it may not be obvious to find a direct reference for either the statement or the proof. Therefore, its inclusion may be helpful to the reader.

Using multiple-scale expansion, we expand the solution as

x=x0+εx1+ε2x2+,λ=λ0+ελ1,
(A2)

where x0 is constant and xi = xi(t, τ) with i > 1 and τ = εt being the slow time.32 The expansion of the vector field close to second order is

f=f0+εDf0x1+fλ,0λ1+ε2Df0x2+12Hf0(x1,x1)+λ1Dfλ,0x1+,
(A3)

where f0 = f(x0, λ0) indicates that f is evaluated at the point (x0, λ0) and

(Df)ij=fixj,  Df0=Df(x0,λ0),
(A4)
fλ=fλ,  fλ,0=fλ(x0,λ0),
(A5)
(Hf)ijk=2fixjxk,  Hf0=Hf(x0,λ0),
(A6)
(Dfλ)ij=2fixjλ,  Dfλ,0=Dfλ(x0,λ0).
(A7)

If the system has a bifurcation of co-dimension one at λ0, then Df0 has rank n − 1 and the vectors w and v,

wDf0=0,Df0v=0,
(A8)

define the left and right null spaces of Df0. The generic condition for a transcritical bifurcation to occur is

wfλ,0=0.
(A9)

To order O(ε), the differential equation reads

ẋ1=Df0x1+λ1fλ,0,
(A10)

and we are interested in the solution

x1=c(τ)v,
(A11)

whose existence is guaranteed by the condition w · fλ,0 = 0. To second order O(ε2), we have

ẋ2+c(τ)v=Df0x2+c212Hf0(v,v)+cλ1Dfλ,0v.
(A12)

The Fredholm alternative gives a condition for the existence of a solution of this inhomogeneous system,

w(c(τ)v)=wc212Hf0(v,v)+cλ1Dfλ,0v,
(A13)

which gives the equation

c(τ)=βλ1c+αc2,
(A14)

where

α=121vwwHf0(v,v),
(A15)
β=1vwwDfλ,0v.
(A16)

Taking into account that ελ1 = λλ0 and defining y = εc, the local solution is x = x0 + yv, where

ẏ=β(λλ0)y+αy2
(A17)

is the normal form of a transcritical bifurcation at λ = λ0. The local evolution of the variables for which vi ≠ 0 is given by

ẋi=β(λλ0)(xix0,i)+αvi(xix0,i)2.
(A18)

For unseeded initial conditions, we can show that the total mass of the system is not conserved. Assume that, for all 2 ≤ i, we have λiλ1 and assume that there exists some index j, with 2 ≤ j, such that the inequality is strict (i.e., λj < λ1). In this case, we have

dM̃dt>λ1i=2N1ipi+2k0+2knm2+2k+mP+2k2σ(m)M̃>λ1i=2ipi+2k0+2knm2+2k+mP+2k2σ(m)M̃=λ1M̃+2k0+2knm2+2k+mP+2k2σ(m)M̃.
(B1)

Likewise, for i = 2, we have a similar inequality

dp2dt=λ2p2+k0+knm22k+mp2+k2σ(m)M̃>λ1p2+k0+knm22k+mp2+k2σ(m)M̃
(B2)

and likewise for i > 2. The above observation shows that the system (25)(28) grows faster than the constant-clearance case system where λi = λ1 for every i1,2,. We note that, as in Sec. III, the total system mass for (25)(28) is M̃tot=M̃+m; this follows from the common definition of M̃, here, and M [see (3)]. Adding (25) to (26) and using (B1) give

dM̃totdt>λ1m0M̃tot.
(B3)

In the presence of the unseeded initial conditions (29), we have that M̃tot(0)=m0 so that the left-hand side of (B3) is strictly positive and mass conservation is violated at the outset. Now, let Mtotλ1 denote the total mass of the constant clearance case λi = λ1 for all i. We know that, in the presence of unseeded initial conditions, a system with constant clearance conserves mass so that

dMtotλ1dt=0.

From (B1) and (B2), which holds analogously for i > 2 and for i = 1, we have equality, we can conclude that

dM̃totdtdMtotλ1dt=0
(B4)

for unseeded initial conditions. Taking together, (B3) implies that the system (25)(28), with unseeded initial conditions, initially gains mass, while (B4) shows that it can never lose mass. Therefore, not only does (25)(28) not conserve mass, but it can never return to the state of initial unseeded mass.

For ν = 1, case (55) takes the form

Δi=m=3ibb+mλ=(b+λ)(b+2λ)b2bλib+λλi1,i3,
(C1)

where the subscript (x)i = x(x + 1) (x + 2)⋯(x + i − 1) denotes the ascending factorial (i.e., the Pochhammer symbol). Defining ξ = −1, then (C1) is satisfied provided

limiξi(ξ+1)i=0.
(C2)

The function ξi((ξ+1)i)1 is monotonically decreasing in both ξ and i, and condition C1 is satisfied for any ξ > 0. Using ξ = −1, expression (C1) implies

Δ(1)=2+(λ+λξ)(2λ+λξ)λ2ξ2i=3iξi(ξ+1)i=2+(λ+λξ)(2λ+λξ)λ2ξ2ξ32+3ξ+ξ2=2+ξ.

Thus, we have

Δ(1)=2+bλ,
(C3)

and it follows that p2*, for ν = 1, is determined by the formula

p2*=λ(k0knm02)2(k+m0+λ)(λk2σ0).
(C4)
1.
J.
Hardy
and
G.
Higgins
, “
Alzheimer’s disease: The amyloid cascade hypothesis
,”
Science
256
,
184
186
(
1992
).
2.
J.
Hardy
and
D.
Allsop
, “
Amyloid deposition as the central event in the aetiology of Alzheimer’s disease
,”
Trends Pharmacol. Sci.
12
,
383
388
(
1991
).
3.
D. J.
Selkoe
and
J.
Hardy
, “
The amyloid hypothesis of Alzheimer’s disease at 25 years
,”
EMBO Mol. Med.
8
,
595
608
(
2016
).
4.
S.
Linse
,
T.
Scheidt
,
K.
Bernfur
,
M.
Vendruscolo
,
T. P. J.
Knowles
, and
O.
Hansson
, “
Kinetic fingerprints differentiate the mechanisms of action of anti-aβ antibodies
,”
Nat. Struct. Mol. Biol.
27
,
1125
1133
(
2020
).
5.
A.
Bacyinski
,
M.
Xu
,
W.
Wang
, and
J.
Hu
, “
The paravascular pathway for brain waste clearance: Current understanding, significance and controversy
,”
Front. Neuroanat.
11
,
101
(
2017
).
6.
H.
Benveniste
,
X.
Liu
,
S.
Koundal
,
S.
Sanggaard
,
H.
Lee
, and
J.
Wardlaw
, “
The glymphatic system and waste clearance with brain aging: A review
,”
Gerontology
65
,
106
119
(
2019
).
7.
J. M.
Tarasoff-Conway
,
R. O.
Carare
,
M. J.
de Leon
 et al., “
Clearance systems in the brain–implications for Alzheimer disease
,”
Nat. Rev. Neurol.
11
,
457
470
(
2015
).
8.
S.-H.
Xin
,
L.
Tan
,
X.
Cao
,
J.-T.
Yu
, and
L.
Tan
, “
Clearance of amyloid beta and tau in Alzheimer’s disease: From mechanisms to therapy
,”
Neurotoxic. Res.
34
,
733
748
(
2018
).
9.
G.
Meisl
,
L.
Rajah
,
S. A. I.
Cohen
,
M.
Pfammatter
,
A.
Šarić
,
T. P. J.
Knowles
et al., “
Scaling behaviour and rate-determining steps in filamentous self-assembly
,”
Chem. Sci.
8
,
7087
7097
(
2017
).
10.
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 behavior underlie the contrasting aggregation kinetics of the Aβ40 and Aβ42 peptides
,”
Proc. Natl. Acad. Sci. U. S. A.
111
,
9384
9389
(
2014
).
11.
S. I.
Cohen
,
M.
Vendruscolo
,
M. E.
Welland
,
C. M.
Dobson
,
E. M.
Terentjev
, and
T. P. J.
Knowles
, “
Nucleated polymerization with secondary pathways. I. Time evolution of the principal moments
,”
J. Chem. Phys.
135
,
065105
(
2011
).
12.
S. I.
Cohen
,
M.
Vendruscolo
,
C. M.
Dobson
, and
T. P.
Knowles
, “
Nucleated polymerization with secondary pathways. II. Determination of self-consistent solutions to growth processes described by non-linear master equations
,”
J. Chem. Phys.
135
,
065106
(
2011
).
13.
S. I.
Cohen
,
M.
Vendruscolo
,
C. M.
Dobson
, and
T. P.
Knowles
, “
Nucleated polymerization with secondary pathways. III. Equilibrium behavior and oligomer populations
,”
J. Chem. Phys.
135
,
065107
(
2011
).
14.
G.
Meisl
,
J. B.
Kirkegaard
,
P.
Arosio
,
T. C. T.
Michaels
,
M.
Vendruscolo
,
C. M.
Dobson
,
S.
Linse
, and
T. P. J.
Knowles
, “
Molecular mechanisms of protein aggregation from global fitting of kinetic models
,”
Nat. Protocols
11
,
252
(
2016
).
15.
R.
Frankel
,
M.
Törnquist
,
G.
Meisl
,
O.
Hansson
,
U.
Andreasson
,
H.
Zetterberg
,
K.
Blennow
,
B.
Frohm
,
T.
Cedervall
,
T. P.
Knowles
 et al, “
Autocatalytic amplification of Alzheimer-associated Aβ42 peptide aggregation in human cerebrospinal fluid
,”
Commun. Biol.
2
,
365
(
2019
).
16.
F.
Kundel
,
L.
Hong
,
B.
Falcon
,
W. A.
McEwan
,
T. C. T.
Michaels
,
G.
Meisl
,
N.
Esteras
,
A. Y.
Abramov
,
T. J. P.
Knowles
,
M.
Goedert
 et al, “
Measurement of tau filament fragmentation provides insights into prion-like spreading
,”
ACS Chem. Neurosci.
9
,
1276
1282
(
2018
).
17.
A. J.
Dear
,
G.
Meisl
,
T. C. T.
Michaels
,
M. R.
Zimmermann
,
S.
Linse
, and
T. P. J.
Knowles
, “
The catalytic nature of protein aggregation
,”
J. Chem. Phys.
152
,
045101
(
2020
).
18.
S. S.
Rogers
,
P.
Venema
,
L. M. C.
Sagis
,
E.
van der Linden
, and
A. M.
Donald
, “
Measuring the length distribution of a fibril system: A flow birefringence technique applied to amyloid fibrils
,”
Macromolecules
38
,
2948
2958
(
2005
).
19.
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 solution to the kinetics of breakable filament assembly
,”
Science
326
,
1533
1537
(
2009
).
20.
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 of amyloid-β42 aggregates occurs through a secondary nucleation mechanism
,”
Proc. Natl. Acad. Sci. U. S. A.
110
,
9758
9763
(
2013
).
21.
S.
Linse
,
T.
Scheidt
,
K.
Bernfur
,
M.
Vendruscolo
,
C.
Dobson
,
S.
Cohen
,
E.
Sileikis
,
M.
Lundquist
,
F.
Qian
,
T.
O’Malley
 et al, “
Kinetic fingerprint of antibody therapies predicts outcomes of Alzheimer clinical trials
,” bioRxiv:815308 (
2019
).
22.
D. L.
Brody
,
H.
Jiang
, and
N.
Wildburger
, “
Non-canonical soluble amyloid-beta aggregates and plaque buffering: Controversies and future directions for target discovery in Alzheimer’s disease
,”
Alz. Res. Therapy
9
,
62
(
2017
).
23.
M.
Bertsch
,
B.
Franchi
,
N.
Marcello
,
M. C.
Tesi
, and
A.
Tosin
, “
Alzheimer’s disease: A mathematical model for onset and progression
,”
Math. Med. Biol.
34
,
193
(
2016
).
24.
S.
Fornari
,
A.
Schäfer
,
E.
Kuhl
, and
A.
Goriely
, “
Spatially-extended nucleation-aggregation-fragmentation models for the dynamics of prion-like neurodegenerative protein-spreading in the brain and its connectome
,”
J. Theor. Biol.
486
,
110102
(
2020
).
25.
C. E.
Mays
,
J.
van der Merwe
,
C.
Kim
,
T.
Haldiman
,
D.
McKenzie
,
J. G.
Safar
, and
D.
Westaway
, “
Prion infectivity plateaus and conversion to symptomatic disease originate from falling precursor levels and increased levels of oligomeric PrPSc species
,”
J. Virol.
89
,
12418
12426
(
2015
).
26.
G.
Meisl
, “
Modelling protein aggregation in vitro and in vivo
,” D.Phil. dissertation (
University of Cambridge
,
2016
).
27.
J.
Sotomayor
, “
Generic bifurcations of dynamical systems
,” in
Dynamical Systems
(
Elsevier
,
1973
), pp.
561
582
.
28.
J.
Guckenheimer
and
P.
Holmes
,
Nonlinear Oscillations, Dynamical Systems and Bifurcations of Vector Fields
(
Springer-Verlag
,
New York
,
1983
).
29.
S.
Wiggins
,
Global Bifurcations and Chaos
(
Springer-Verlag
,
New York, Berlin
,
1988
).
30.
A.
Goriely
,
Integrability and Nonintegrability of Dynamical Systems
(
World Scientific Publishing Company
,
2001
).
31.
A.
Goriely
, “
Painlevé analysis and normal forms theory
,”
Physica D
152-153
,
124
144
(
2001
).
32.
A.
Newell
, “
Envelope equations
,” in
Lectures in Applied Mathematics, Vol. 15, Nonlinear Wave Motion
(American Mathematical Society,
Providence, RI
,
1974
), pp.
157
163
.