Numerical and Analytical Simulation of the Growth of Amyloid-β Plaques

Alzheimer's disease (AD) is a devastating neurodegenerative disorder, impacting nearly 50 × 10⁶ individuals globally. While recent therapies approved by the U.S. Food and Drug Administration may decelerate its progression, there are currently no definitive cures for AD [1–4]. It is crucial to develop mechanistic models of processes involved in AD to pave the way for potential future treatments.

Senile plaques, primarily composed of accumulated amyloid-β (Aβ) peptides, represent a key characteristic of AD [5]. A significant surge of interest in Aβ aggregation was generated by the amyloid cascade hypothesis, which posited that the formation of Aβ aggregates initiates all other pathological processes in AD [6]. While this hypothesis has undergone extensive investigation, it has also faced recent criticism from some authors [6–9] due to the failure of many Aβ-targeted clinical trials for AD [10]. However, lifelong Aβ reduction may mitigate dementia risk [11]. Multitarget drugs, addressing Aβ among other factors, may be promising for AD treatment. Furthermore, Aβ serves as a crucial biomarker in AD diagnosis [12–15].

Previous mechanistic models of AD [16–18] primarily focused on elucidating intraneuronal mechanisms contributing to the development of AD, namely, the formation of neurofibrillary tangles composed of aggregated tau proteins and the generation of Aβ monomers at the cell membrane. However, these investigations did not simulate the gradual extracellular buildup of Aβ plaques. In contrast, the present study aims to replicate the gradual growth of an Aβ plaque, a process that can span decades in humans. This process involves the generation of amyloid precursor protein within neurons and its subsequent cleavage by β- and γ-secretases at lipid membranes [14,19–21]. The majority of these newly formed Aβ monomers are released into the extracellular environment, where they hold the potential to aggregate [1]. The purpose of this paper is to develop a mechanistic model predicting the radial increase of an Aβ plaque with time.

where A refers to a monomeric protein and B represents a protein that has undergone amyloid conversion. The kinetic constants, denoted as $k1$ and $k2$⁠, correspond to the rates of nucleation and autocatalytic growth, respectively [22]. The primary nucleation process, as defined by Eq. (1), exclusively involves monomers. Conversely, secondary nucleation, described by Eq. (2), involves both monomers and pre-existing aggregates of the same peptide [24].

In Ref. [22], the Finke–Watzky (F–W) model was utilized to fit previously published data on Aβ aggregation [25,26]. Typically, these data are presented for two isoforms, Aβ40 and Aβ42, with the latter being more prone to aggregation. In this context, the F–W model is applied to simulate the conversion of monomers, whose concentration is denoted as $CA$⁠, into aggregates, whose concentration is denoted as $CB$⁠. These aggregates encompass various forms of Aβ oligomers, protofibrils, and fibrils [5]. Aβ aggregates formed through this process assemble in amyloid plaques. It is worth noting that the simplified nature of the F–W model does not allow for differentiation between various aggregate types and sizes. The current model does not simulate the process of adhesive Aβ fibrils assembling into an Aβ plaque, which is assumed to occur faster than the formation of Aβ fibrils.

The first term on the right-hand side of Eq. (3) simulates the rate of conversion of Aβ monomers into aggregates while the second term simulates the rate of this conversion through autocatalytic growth. These expressions for the first and second terms arise from applying the law of mass action to the reactions described by Eqs. (1) and (2), respectively. The third term on the right-hand side of Eq. (3) accounts for the rate at which Aβ monomers degrade. The fourth term represents the rate of Aβ monomer production, as Aβ monomers are continuously generated through the cleavage of amyloid precursor protein [28,29]. The production of Aβ monomers occurs as a surface phenomenon at lipid membranes [5,30]. As the precise area of the lipid membrane contributing to the production of Aβ monomers per CV is unknown, the production was normalized by the area of a single CV face, $L2$⁠. The flux $qA,0$ should be interpreted as the average flux of Aβ monomers into the CV per unit area of a single face. The total rate at which monomers enter the CV is $qA,0L2$⁠. Utilizing the lumped capacitance method, the rate of Aβ monomers entering the CV is converted into volumetric monomer generation, providing the same number of monomers per unit time, $qA,0L2/L3$⁠. This conversion explains why the fourth term on the left-hand side of Eq. (3) is inversely related to L.

The first term on the right-hand side of Eq. (4) represents the rate at which Aβ aggregates are produced through nucleation while the second term simulates the rate of their production via autocatalytic growth. These terms have equal magnitudes, but opposite signs when compared to the first and second terms in Eq. (3). This is because, in the F-W model, the rate of aggregate production equals the rate of monomer disappearance. The third term in Eq. (4) represents the rate of Aβ aggregate degradation, assuming this degradation is significantly slower than that of monomers [31].

The sole independent variable in the model is time, t. The dependent variables are listed in Table 1 while the parameters employed in the model are listed in Table 2.

Table 3 presents the dimensionless independent variable while Table 4 summarizes the dimensionless dependent variables. Table 5 provides the dimensionless parameters used in the model.

Equations (6) and (7) constitute a pair of ordinary differential equations. Matlab's ODE45 solver (Matlab R2020b, MathWorks, Natick, MA) was employed to solve these equations subject to initial conditions given by Eq. (8), as it is thoroughly validated. To ensure accuracy of the solution, the error tolerance parameters, RelTol and AbsTol, were set to 1 × 10⁻¹⁰.

The increase in $CA∗+CB∗$ over time is attributed to the continuous generation of Aβ monomers.

where $erfi(x)$ is the imaginary error function.

As $t∗→∞$⁠, Eq. (16) predicts that $CA∗→0$⁠. It is important to note that Eqs. (16) and (17) are useful for larger times. However, they are not applicable when $t∗→0$⁠, for instance, Eq. (16) yields $CA∗(0)=qA,0∗$⁠, which contradicts the initial condition specified in Eq. (8a).

where $MW$ denotes the average molecular weight of an Aβ monomer.

where $a21=1021μM μm3mol$ denotes the conversion factor from $mol μm3$ to $μM$⁠.

Equation (26) explains the deceleration of Aβ plaque growth over time implying smaller (younger) plaques tend to undergo faster growth compared to larger plaques, a phenomenon observed in previous research [42].

suggesting that at small times, the radius of an Aβ plaque depends on the kinetic constant that describes the nucleation of Aβ aggregates.

Interestingly, Eq. (25) predicts that the radius of an Aβ plaque is independent of the kinetic constants describing the rate of Aβ plaque agglomeration, $k1$ and $k2$⁠. However, Eq. (25) is valid only at large times. At small times, the plaque growth is described by Eq. (28), which suggests that the initiation of the plaque formation process is driven by kinetic factors. The kinetic-driven initial stage of plaque growth agrees with the fast appearance of plaques in a mouse model observed in Ref. [43]. However, the continued growth of the plaques is restricted by the rate of Aβ monomer production and, if the half-lives of the monomers and aggregates are finite, by the rate at which the monomers and aggregates degrade. Plaque growth restriction by Aβ peptide production at extended times is consistent with observations reported in Ref. [44].

The positive values of $SL rABP$ and $SqA,0 rABP$ indicate that an increase in $L$ and $qA,0$ will result in an increase in $rABP$⁠.

where $Δ T1/2,A=10−3T1/2,A$ represents the step size. Computations were conducted with various step sizes to assess whether the sensitivity coefficients remained unaffected by changes in the step size.

Figures 2–7 were generated using the parameter values listed in Table 2, unless otherwise specified. Dysfunctional autophagic/lysosomal pathways [50,51] and impaired ubiquitin-proteasome systems [52] may potentially contribute to AD, as these defects can lead to the accumulation of misfolded proteins. Therefore, one of the scenarios investigated in the figures involves dysfunctional protein degradation machinery, potentially resulting in infinitely large half-lives of Aβ monomers and/or aggregates.

In the scenario where Aβ monomers have an infinite half-life, the concentration of Aβ aggregates, $CB$⁠, increases linearly with time, consistent with the approximate solution provided by Eq. (17) (Fig. 2(a)). However, it is worth noting that this approximate solution is not valid at small times within the range of [0–0.04 year] (Fig. 2(b)). When the machinery responsible for degrading Aβ monomers works properly (⁠$T1/2,A=4.61×103$ s, $T1/2,B→∞$⁠), the computational results show an approximately 2-year delay before the linear increase in Aβ aggregate concentration begins (Fig. 2(a)). This suggests that functional proteolytic machinery, which degrades Aβ monomers, can postpone the formation of Aβ plaques by several years [53]. If Aβ aggregates degrade, in addition to monomer degradation (⁠$T1/2,A=4.61×103$ s, $T1/2,B=2.31×104$ s), no aggregate growth is observed, as depicted in Fig. 2.

Notably, the concentrations of Aβ aggregates, denoted as $CB$⁠, plotted in Fig. 2, are computed without considering the removal of aggregates from the cytosol due to their deposition into the amyloid plaques. The model assumes that all Aβ aggregates eventually become part of the plaques, see Eqs. (22) and (25).

In the scenario where Aβ monomers have an infinite half-life, the concentration of Aβ monomers, $CA$⁠, decreases as time progresses (Figs. 3(a) and 3(b)). This decline is a consequence of the increasing concentration of Aβ aggregates (Fig. 2(a)), leading to faster autocatalytic conversion of monomers into aggregates as described by Eq. (4). The approximate solution given by Eq. (16) closely matches the numerical solution (Fig. 3(a)), but only for times greater than 0.04 years (Fig. 3(b)). In the scenario where the machinery responsible for degrading Aβ monomers functions normally (⁠$T1/2,A=4.61×103$ s, $T1/2,B→∞$⁠), there is approximately a 2-year delay in the decrease of Aβ monomer concentration. This delay is due to the postponed production of Aβ aggregates (as shown in Fig. 2(a)). These aggregates play a role in the autocatalytic conversion of monomers into aggregates (refer to Eq. (3)), which accelerates the monomer decay. If both monomers and Aβ aggregates degrade (⁠$T1/2,A=4.61×103$ s, $T1/2,B=2.31×104$ s), the concentration of monomers stays at a constant value independent of time (Fig. 2). The increase in monomer concentration due to aggregate decay is explained by the catalytic role of aggregates in converting monomers to aggregates.

In the scenario of infinite half-life of Aβ monomers, the radius of a growing Aβ plaque, $rABP$⁠, increases in proportion to the cube root of time, $t1/3$⁠, in agreement with the approximate solution provided in Eq. (25) (Fig. 4(a)). The agreement between the approximate and numerical solutions is excellent, except for when t is within the smaller time range of [0–0.04 year] (Fig. 4(b)). This implies that during early stages, the growth of the Aβ plaque is determined by the kinetic conversion of Aβ monomers into aggregates, resulting in a sigmoidal increase in the plaque's radius, consistent with observations in Ref. [44]. Conversely, during later stages, the growth is determined by Aβ production, as assumed in the approximate solution.

If the machinery responsible for Aβ monomer degradation functions normally (⁠$T1/2,A=4.61×103$ s, $T1/2,B→∞$⁠), the computational result predicts a slower increase rate for the radius of an Aβ plaque over time (Fig. 4(a)). This suggests that the cube root hypothesis is applicable only when at least 0.04 years have passed since the beginning of Aβ aggregation, and the Aβ monomer degradation mechanism is not functioning correctly. If Aβ aggregates also degrade (⁠$T1/2,A=4.61×103$ s, $T1/2,B=2.31×104$ s), the growth of Aβ plaque does not occur (Fig. 4).

Figure S1(a) available in the Supplemental Materials on the ASME Digital Collection investigates how the growth of the plaque radius is affected by slower Aβ monomer degradation than physiologically relevant. All lines are computed for $T1/2,B→∞$⁠. In addition to the physiologically relevant half-life of monomers (⁠$T1/2,A=4.61×103$ s), situations were studied where the half-life is twice and ten times longer than physiologically relevant. The situation corresponding to an infinite half-life of monomers (⁠$T1/2,A→∞$⁠) was also plotted. The increase of $T1/2,A$ by a factor of two from physiologically relevant (⁠$T1/2,A=2×4.61×103$ s) brings the plaque radius halfway closer to the scenario characterized by $T1/2,A→∞$⁠. When $T1/2,A$ is increased by a factor of ten from physiologically relevant (⁠$T1/2,A=10×4.61×103$ s), the growth of the plaque radius is nearly identical to the situation with an infinite half-life of monomers (see Fig. S1(a) available in the Supplemental Materials).

Figure S1(b) available in the Supplemental Materials investigates how the growth of the plaque radius is influenced by slower Aβ aggregate degradation than physiologically relevant. All curves are calculated for $T1/2,A→∞$⁠. In addition to the physiologically relevant half-life of aggregates (⁠$T1/2,B=2.31×104$ s), situations were examined where the half-life was one hundred and ten thousand times longer than physiologically relevant. The scenario corresponding to an infinite half-life of aggregates (⁠$T1/2,B→∞$⁠) was also depicted. For the assumed physiologically relevant value of aggregate half-life (⁠$T1/2,B=2.31×104$ s), the plaque radius rapidly increases to approximately 0.1 μm and then remains constant. If $T1/2,B$ is increased by a factor of one hundred from physiologically relevant (⁠$T1/2,B=102×2.31×104$ s), the plaque radius rapidly increases to approximately 0.9 μm and then remains constant. When $T1/2,B$ is increased by a factor of ten thousand from physiologically relevant (⁠$T1/2,B=104×2.31×104$ s), the plaque radius increases similarly to the situation with an infinite half-life of aggregates but remains approximately 20% below the curve corresponding to $T1/2,B→∞$ (see Fig. S1(b) available in the Supplemental Materials).

In Fig. 5, a hypothetical scenario is explored where a therapeutic intervention halts the production of Aβ monomers after 10 years. In the absence of Aβ aggregate degradation, plaque growth ceases after 10 years due to the lack of new monomer production (Fig. 5(a)). However, if Aβ aggregates are degraded (e.g., by autophagy [54]), the cessation of monomer supply results in the reduction (⁠$T1/2,B=104×2.31×104$ s) or complete disappearance (⁠$T1/2,B=2.31×104$ s and $T1/2,B=102×2.31×104$ s) of the Aβ plaques (Fig. 5(b)).

Figure 6(a) illustrates the relationship between the radius of Aβ plaques after 20 years of growth and the half-distance between these plaques, L. When $T1/2,A→∞$⁠, $rABP$ increases with L almost linearly. This occurs because, for a uniform volumetric production rate of Aβ monomers, a smaller number of plaques (larger distances between them) allows each individual plaque to access a greater supply of Aβ, thus enabling it to grow to a larger size. If $L=500$μm, the plaque's diameter after 20 years of growth is approximately 50 μm, which agrees with the representative diameter of an Aβ plaque observed in Ref. [33]. In the case of infinite half-life of Aβ monomers, the approximate solution provided by Eq. (25) precisely matches the numerical solution (Fig. 6(a)). In the scenario where Aβ monomer degradation functions normally (⁠$T1/2,A=4.61×103$ s, $T1/2,B→∞$⁠), an increase in L first results in an increase in $rABP$⁠, and subsequently (after reaching approximately $L=250$μm) leads to a decrease in $rABP$ (as shown in Fig. 6(a)). This occurs due to an increase in L that leads to a proportional increase in the supply of Aβ monomers into the CV through a single face, scaling with L². The number of Aβ monomers degraded within the CV scales with L³, signifying that Aβ monomer degradation becomes the prevailing process within the CV at a certain point. If Aβ aggregates also degrade (⁠$T1/2,A=4.61×103$ s, $T1/2,B=2.31×104$ s), there is no growth of Aβ plaque, regardless of L (Fig. 6(a)).

Figure 6(b) illustrates how the radius of Aβ plaques after 20 years of growth is influenced by the rate of Aβ monomer production. When $T1/2,A→∞$ and $T1/2,B→∞$⁠, $rABP$ is directly proportional to $qA,01/3$ because a higher production rate of Aβ monomers results in a greater supply of Aβ available for plaque growth (Fig. 6(b)). In the scenario where Aβ monomers have an infinite half-life, the approximate solution given by Eq. (25) precisely agrees with the numerical solution (Fig. 6(b)). However, if the machinery responsible for Aβ monomer degradation is functioning normally (⁠$T1/2,A=4.61×103$ s, $T1/2,B→∞$⁠), the plaque grows nearly at the same rate as it does for $T1/2,A→∞$⁠. This is because Aβ monomer production is continuous, and the monomers convert into Aβ plaques before they can be broken down. If Aβ aggregates also undergo degradation (⁠$T1/2,A=4.61×103$ s, $T1/2,B=2.31×104$ s), there is no growth of the Aβ plaque, irrespective of $qA,0$ (Fig. 6(b)).

Figure 6(c) depicts the correlation between the Aβ plaque radius after 20 years of growth and the half-life of Aβ monomers. If $T1/2,B→∞$ (the line marked by squares), the numerical solution produced a sigmoidal curve. For the smallest value of $T1/2,A$ in Fig. 6(c) (10 s), the plaque radius is zero. It gradually increases with the half-life of monomers, reaching a value corresponding to the infinite half-life of monomers at $T1/2,A=105$ s. Beyond this point, it aligns with the plaque radius predicted by the approximate solution, given by Eq. (25). If $T1/2,B=2.31×104$ s (the line marked by triangles), the plaque radius shows no growth when $T1/2,A<106$ s and very slow growth for larger values of $T1/2,A$⁠, approaching approximately 1 μm for $T1/2,A=1012$ s (Fig. 6(c)).

Figure 6(d) depicts the correlation between the Aβ plaque radius after 20 years of growth and the half-life of Aβ aggregates. If $T1/2,A→∞$ (the line marked by squares), for $T1/2,B<102$ s, the plaque radius is zero. It gradually increases with the half-life of aggregates, reaching a value corresponding to the infinite half-life of monomers at $T1/2,B=1010$ s. Beyond this point, it aligns with the plaque radius predicted by the approximate solution given by Eq. (25). If $T1/2,A=4.61×103$ s (the line marked by triangles), the plaque radius shows no growth when $T1/2,B<106$ s and S-shaped growth for larger values of $T1/2,B$⁠, approaching approximately 4.7 μm for $T1/2,B=1012$ s (Fig. 6(d)).